1 Introduction

Tight/shale oil resource extensively distributes around the world. It is regarded as a promising resource to provide fossil energy in the future. Therefore, how to efficiently exploit these resources comes into our focus. However, the permeability is on the order of 10–3 to 10–1 millidarcy in this kind of reservoir. Industrial capacity cannot be yielded even using conventional horizontal wells. Extensive field experience shows that hydraulic fracturing can effectively improve the development effect of unconventional resources (Holditch and Tschirhart 2005; Liang et al. 2017). Nevertheless, the oil recovery is still lower than 10% by depletion-drive after hydraulic fracturing (Manrique et al. 2010; Kathel and Mohanty 2013; Teklu et al. 2016). Therefore, numerous attempts have been made to enhance oil recovery. In some oilfields, water flooding has been conducted. However, if the well spacing is larger, the producers are not responding under a very high injection pressure (Song and Yang 2013; Kong et al. 2016; Mansour et al. 2017). Conversely, if the well spacing is smaller, water channels faster to the producers through the hydraulic fractures (Thomas et al. 2014; Li et al. 2015; Pourabdollah 2018). Consequently, huff-n-puff by fracturing fluids or produced water has been tested in several tight oil pilots, and positive effects have been achieved (Li et al. 2015; Li 2015; Tuero et al. 2017).

Generally, the oil production mechanism of water soaking is attributed to imbibition (Wang et al. 2018). Imbibition is significantly important in tight oil reservoir since capillary force is more dominant in such Formation (Yang et al. 2018). Most studies suggest that oil can be driven out by imbibition even in mixed-wet samples (Cai et al. 2014; Dutta et al. 2012; Habibi et al. 2015). Dutta et al. (2012) reported that although the permeability of tight rock is very low, the small characteristic radius suppresses a stronger effect of capillarity; the impacts of permeability and porosity on imbibition should be taken into account. Habibi et al. (2015) proposed that pore surface usually contains both oil-wet and water-wet minerals in tight rock, so the remaining oil might be trapped in small oil-wet pores. Wang et al. (2012) found that surfactant can alter the oil-wet or mixed-wet cores toward water-wet to displace out more oil than brine.

Spontaneous imbibition plays a significant role in water uptake into porous media during long shut-in which could last for a few days in some cases. After a long time, when the water-phase pressure that pushes the water into the formation has already dissipated, capillary forces dominate (Rangel-German and Kovscek 2002). However, forced imbibition is different from spontaneous imbibition, whose driving force is capillary force alone. Forced imbibition is achieved by injecting fluid into a sample at a constant pressure higher than that of the sample displacement pressure (Roychaudhuri et al. 2014). The wetting phase firstly enters the capillary for a certain length to enhance the elastic energy. After pressure equilibrium, spontaneous imbibition happens to expel the oil. Roychaudhuri et al. (2014) conducted forced imbibition to study the efficacy of surfactant under high-pressure scenarios and found that both spontaneous and forced imbibition experiments should be applied to evaluate the effect of surfactant on liquid dynamics in tight shales. Actually, the water-phase pressure in macro-fractures is much higher than that in matrix after hydraulic fracturing; the pressure propagates during soaking. Moreover, higher pressure is favorable to recover the conductivity of both macro-fractures and natural fractures by water stimulation. Therefore, when production restarts after soaking, elastic energy of both rock and fluid will contribute to oil production as that in fractured reservoirs (Finkbeiner et al. 2010; Omosebi and Igbokoyi 2012; Valiveti et al. 2015), which increases the complexity of oil production mechanisms.

For the reasons outlined above, there are several issues to be clarified for the strategy as follows: (1) what is/are the driving force(s) to expel oil out from matrix? (2) If the oil is produced by the mechanisms of elastic energy and imbibition, how much does each one contribute? (3) How long should the cycle soaking take? and (4) what are the production performances of huff-n-puff? The answers are the keys to optimizing soaking strategy and improving the development effect. In order to figure out the answers, the whole process of pressurizing, high-pressure soaking, and depressurizing was firstly discussed, and a mechanistic model including the essential mechanisms was established and verified. The confusing problems above are resolved by the model with an excellent matching result. This work could help us understand the mechanisms and capacity of oil production by high-pressure soaking more clearly and better design EOR strategy in tight/shale oil reservoirs.

2 Forced imbibition during high-pressure soaking

Spontaneous imbibition is usually classified as co-current and counter-current. The driving force of co-current imbibition is capillary force and gravity (or called buoyancy), yet that of counter-current imbibition is capillary force. However, imbibition under high-pressure scenario is different. Since the inverse Bond number NB−1 is very large in tight oil reservoirs, the impact of gravity can be neglected. The process can be manifested by the sketch maps of forced co-current and counter-current imbibition in a variable capillary, as shown in Figs. 1 and 2, respectively. Figure 1a shows the initial state before fracturing. When fracturing fluid or water is injected, the pressure increases and water enters the capillary from both ends, which is called high-pressure transmission, as shown in Fig. 1b. After high-pressure transmission, the equilibrium water-phase pressure is Prh; thus, the oil-phase pressure is Prh + Pc1 at the narrow end, and that is Prh + Pc2 at the wide end. Because the radius of narrow end is smaller, Pc1 is larger than Pc2. Consequently, water enters the capillary from the narrow end, while oil and antecedent water flow out from the wide end, which is the well-known mechanism of co-current imbibition, as shown in Fig. 1c. After soaking for a while, the well opens with a lower bottom hole pressure (BHP). The pressure difference gradually increases between the internal and external fluids as the pressure decreases in fractures. In addition, as the fluid expands and pores shrink, a large amount of fluid is driven out of the capillary, which is called elastic energy release as shown in Fig. 1d. There might be questions why the water proportion increases. The reason is that the permeability at the wide end is larger than that of the narrow end, so oil flows faster from the wide end than water from the narrow end. As a result, more oil gets expelled, and the water saturation increases in the capillary. According to the above analysis, the oil expelled due to pressure decrease is regarded as the contribution of elasticity. The process and mechanisms are similar for counter-current imbibition as in Fig. 2. We can see that it is more favorable for water entering and driving oil out from matrix.

Fig. 1
figure 1

Sketch map of pressurizing-soaking-depressurizing process for co-current imbibition scenario (KH, high permeability; KL, low permeability; Pr, reservoir pressure; Prl, low reservoir pressure; Prh, high reservoir pressure; Pc1, capillary pressure at the narrow end; Pc2, capillary pressure at the wide end)

Fig. 2
figure 2

Sketch map of pressurizing-soaking-depressurizing process for counter-current imbibition scenario

3 Simulation of forced imbibition with high-pressure soaking

In order to understand the process and mechanisms of oil production by high-pressure soaking in tight formations more clearly, a simulation method is developed. As can be seen from the above analysis, both imbibition and elastic energy have non-negligible contributions to oil production, so the capillary force, specific gravity, and compressibility of rock and fluids should be included. In tight/shale reservoirs, the rock fractured into discrete blocks around the wellbore is surrounded by the fracturing fluid or water, and the matrix block with its adjacent fractures are regarded as a unit. The oil can be displaced out from the matrix by imbibition, gravity, and pressure difference between matrix and fractures. Because the fractures are filled with proppants, they are regarded as highly permeable porous media. Therefore, the diffusivity equation of water can be written as

$$\nabla \cdot \left[ {\frac{{KK_{{{\text{rw}}}} \rho_{{\text{w}}} }}{{\mu_{{\text{w}}} }}\left( {\nabla p_{{\text{w}}} + \gamma_{{\text{w}}} {\kern 1pt} \nabla D} \right)} \right]{ + }\delta_{{\text{w}}} \rho_{{\text{w}}} q_{{\text{w}}} { = }\frac{\partial }{\partial t}\left( {\rho_{{\text{w}}} S_{{\text{w}}} \phi } \right).$$
(1)

The diffusivity equation of oil under both imbibition and pressure difference is

$$\nabla \cdot \left\{ {\frac{{KK_{{{\text{ro}}}} \rho_{{\text{o}}} }}{{\mu_{{\text{o}}} }}\left[ {\nabla p_{{\text{w}}} { + }\nabla p_{{{\text{cow}}}} (S_{{\text{w}}} ) + \gamma_{{\text{o}}} \nabla D} \right]} \right\}{ + }\delta_{{\text{o}}} \rho_{{\text{o}}} q_{{\text{o}}} { = }\frac{\partial }{\partial t}\left( {\rho_{{\text{o}}} S_{{\text{o}}} \phi } \right),$$
(2)

where K is the absolute permeability; Krw and Kro are the relative permeability of water and oil, respectively; ρo and ρw are the density of oil and water, respectively; μo and μw are the viscosity of oil and water, respectively. pw is the water phase pressure; pcow(Sw) is the capillary pressure between oil and water at a certain water saturation; γo and γw are the gravity of oil and water, respectively; D is the buried depth; qo and qw are the volume flow of source/sink term; So and Sw are the saturation of oil and water, respectively; ϕ is the porosity; t is time; δo and δw are Haviside function of oil and water phase, respectively.

The auxiliary equations are as follows:

$$p_{{{\text{cow}}}} (S_{{\text{w}}} ){ = }p_{{\text{o}}} - p_{{\text{w}}} ,$$
(3)
$$S_{{\text{o}}} { + }S_{{\text{w}}} { = }1.$$
(4)

The compressibility of rock and fluid is the major source of elastic energy. The state equations of rock and fluid are as follows:

$$\rho_{{\text{l}}} = \rho_{{{\text{l0}}}} {\text{e}}^{{C_{{\text{l}}} \left( {p_{{\text{l}}} - p_{0} } \right)}} ,$$
(5)
$$\phi = \phi_{{0}} { + }C_{{\text{r}}} \left( {p_{{\text{l}}} - p_{0} } \right),$$
(6)

where l represents the liquid phase (i.e., water or oil); ρl0 is the density of phase l at reference pressure; ϕ0 is the porosity at reference pressure; Cl is the compressibility of oil or water; Cr is the compressibility of rock; pl is the pressure of phase l; p0 is the reference pressure.

As reported by Corey (1954), the relation between capillary force and water saturation can be approximately described as

$$p_{{{\text{cow}}}} \left( {S_{{\text{w}}} } \right) = C_{{{\text{ow}}}} \sqrt {\frac{1}{{S_{{\text{w}}} }}} ,$$
(7)

where Cow is named as capillary index. To simplify the process, Eq. (7) can be written as

$$p_{{{\text{cow}}}} \left( {S_{{\text{w}}} } \right) = C_{{{\text{ow}}}} \sqrt {\frac{1}{{a\overline{S}_{{\text{w}}} + b}}} ,$$
(8)

where \(\overline{S}_{{\text{w}}}\) is the normalized water saturation

$$\overline{S}_{{\text{w}}} = \frac{{S_{{\text{w}}} - S_{{{\text{wc}}}} }}{{1 - S_{{{\text{o}}{\text{r}}}} - S_{{{\text{wc}}}} }},$$
(9)

where Swc is the connate water saturation; Sor is the residual oil saturation; \(a = 1 - S_{{\text{or}}} - S_{{{\text{wc}}}}\) and \(b = S_{{{\text{wc}}}}\).

Generally, some of tight oil reservoirs are water-wet, so imbibition can happen even if there is no surfactant in the fracturing fluid. However, in some mixed-wet or oil-wet reservoirs, surfactant is usually added into water or the fracturing fluid to promote imbibition by changing the wettability (Begum et al. 2017; Alvarez et al. 2017; Meng et al. 2018; Huang et al. 2020). In order to reflect this effect, the capillary force ratio of water to surfactant can be approximated using the Young–Laplace equation assuming constant pore diameter for water and surfactant imbibition,

$$p_{\cos } = \frac{{p_{{{\text{cow}}}} \sigma_{{{\text{os}}}} \cos \theta_{{{\text{os}}}} }}{{\sigma_{{{\text{ow}}}} \cos \theta_{{{\text{ow}}}} }},$$
(10)

where pcos is the capillary force between oil and surfactant; σow is the interficial tension (IFT) between oil and water; σos is the IFT between oil and surfactant; θow is the contact angle between oil and water; θos is the contact angle between oil and surfactant.

Besides, the capillary force is significantly related to the properties of porous media (Habibi et al. 2015). Hence, the capillary force ratio in different porous media for the same fluid is

$$p_{{{\text{cow}}}} \left( {K,\phi } \right) = p_{{{\text{co}} {\text{wref}}}} \sqrt {\frac{{\phi K_{{{\text{ref}}}} }}{{\phi_{{{\text{ref}}}} K}}} ,$$
(11)

where pcowref is the capillary force between oil and water in reference media; Kref and ϕref are the permeability and porosity of referenced media, respectively.

Combining Eqs. (10) and (11), the ratio of the capillary index for different systems is

$$C_{{{\text{os}}}} = \frac{{C_{{{\text{ow}}}} \sigma_{{{\text{os}}}} \cos \theta_{{{\text{os}}}} }}{{\sigma_{{{\text{ow}}}} \cos \theta_{{{\text{ow}}}} }},$$
(12)

where Cow is the capillary index for the oil–water system; Cos is the capillary index for the oil-surfactant system.

Analogously, the ratio of the capillary index for water imbibition in different media is

$$C_{{{\text{ow}}}} \left( {K,\phi } \right) = C_{{{\text{owref}}}} \sqrt {\frac{{\phi K_{{{\text{ref}}}} }}{{\phi_{{{\text{ref}}}} K}}} ,$$
(13)

where Cowref is the capillary index for the oil–water system in referenced media.

The impacts of surfactant in any porous media can be included by combining Eqs. (12) and (13):

$$C_{{{\text{os}}}} { = }C_{{{\text{owref}}}} \frac{{\sigma_{{{\text{os}}}} \cos \theta_{{{\text{os}}}} }}{{\sigma_{{{\text{ow}}}} \cos \theta_{{{\text{ow}}}} }}\sqrt {\frac{{\phi K_{{{\text{ref}}}} }}{{\phi_{{{\text{ref}}}} K}}} .$$
(14)

The added surfactant also affects the residual oil saturation by changing IFT. Because the flow rate is very low for imbibition, the relation between the residual oil saturation and IFT can be simplified to the following form

$$S_{{{\text{or}}}} = \frac{{S_{{{\text{orw}}}} }}{{1 + {\psi \mathord{\left/ {\vphantom {\psi {\sigma_{{{\text{ow}}}} }}} \right. \kern-\nulldelimiterspace} {\sigma_{{{\text{ow}}}} }}}}.$$
(15)

Here, Sorw and ψ can be attained by matching the experimental data.

Another important parameter for imbibition is the relative permeability. The added surfactant will affect the relative permeability by changing the residual oil saturation (Babadagli 2003; Lu et al. 2014; Goudarzi et al. 2015). Referring to the relative permeability formula of wetting and non-wetting phases achieved by Burdine (1953) and combining Eq. (8), the relative permeability equations for the oil–water system can be obtained by integration,

$$K_{{{\text{rw}}}} = K_{{{\text{rwor}}}} \frac{{\overline{S}_{{\text{w}}}^2 (a\overline{S}_{{\text{w}}}^2 + 2b\overline{S}_{{\text{w}}} )}}{a + 2b},$$
(16)
$$K_{{{\text{ro}}}} = K_{{{\text{rowc}}}} \frac{{(1 - \overline{S}_{{\text{w}}} )^{2} (a + 2b - a\overline{S}_{{\text{w}}}^2 - 2b\overline{S}_{{\text{w}}} )}}{a + 2b},$$
(17)

where Krwor is the water-phase relative permeability at residual oil saturation; Krowc is the oil-phase relative permeability at connate water saturation. Because a and b are the functions of residual oil and connate water saturation, the impact of surfactant on relative permeability is included.

4 Solution and validation

4.1 Solution method

In this work, the IMPES method is used to solve the pressure and saturation equations. Combining Eqs. (1) and (2), the following equation is achieved:

$$\frac{1}{{\rho_{{\text{w}}} }}\nabla \cdot \left\{ {\frac{{KK_{{{\text{rw}}}} \rho_{{\text{w}}} }}{{\mu_{{\text{w}}} }}\left[ {\nabla p_{{\text{o}}} - \nabla p_{{{\text{cow}}}} (S_{{\text{w}}} ) + \gamma_{{\text{w}}} {\kern 1pt} \nabla D} \right]} \right\}{ + }\delta_{{\text{w}}} q_{{\text{w}}} { + }\frac{1}{{\rho_{{\text{o}}} }}\nabla \cdot \left[ {\frac{{KK_{{{\text{ro}}}} \rho_{{\text{o}}} }}{{\mu_{{\text{o}}} }}\left( {\nabla p_{{\text{o}}} + \gamma_{{\text{o}}} \nabla D} \right)} \right]{ + }\delta_{{\text{o}}} q_{{\text{o}}} = \phi c_{{\text{t}}} \frac{{\partial p_{{\text{o}}} }}{\partial t}$$
(18)

with

$$c_{{\text{t}}} = c_{{\text{p}}} + c_{{\text{o}}} S_{{\text{o}}} + c_{{\text{w}}} S_{{\text{w}}} ,$$
$$c_{{\text{p}}} = \frac{1}{\phi }\frac{\partial \phi }{{\partial p_{{\text{o}}} }},$$
$$c_{{\text{o}}} = \frac{1}{{\rho_{{\text{o}}} }}\frac{{\partial \rho_{{\text{o}}} }}{{\partial p_{{\text{o}}} }},$$
$$c_{{\text{w}}} = \frac{1}{{\rho_{{\text{w}}} }}\frac{{\partial \rho_{{\text{w}}} }}{{\partial p_{{\text{o}}} }}.$$

Euler differential method is applied to establish the difference scheme,

$$\frac{1}{{\rho_{{\text{w}}} }}\Delta T_{{\text{w}}} \Delta \left( {p_{{\text{o}}} - p_{{{\text{cow}}}} + \gamma_{{\text{w}}} {\kern 1pt} D} \right){ + }\delta_{{\text{w}}} Q_{{\text{w}}} + \frac{1}{{\rho_{{\text{o}}} }}\Delta T_{{\text{o}}} \Delta \left( {p_{{\text{o}}} + \gamma_{{\text{o}}} D} \right){ + }\delta_{{\text{o}}} Q_{{\text{o}}} = \frac{{\left( {V_{{\text{p}}} c_{{\text{t}}} } \right)^{n} }}{\Delta t}\left( {p_{{\text{o}}}^{n + 1} - p_{{\text{o}}}^{n} } \right)$$
(19)

with

$$V_{{\text{p}}} = V_{{\text{b}}} \phi ,$$
$$Q_{{\text{w}}} = V_{{\text{b}}} q_{{\text{w}}} ,$$
$$Q_{{\text{o}}} = V_{{\text{b}}} q_{{\text{o}}} .$$

Thus, the pressure equation in terms of the oil phase pressure po is:

$$a_{i,j,k} p_{{{\text{o}}i - 1,j,k}}^{n + 1} + c_{i,j,k} p_{{{\text{o}}i,j - 1,k}}^{n + 1} + e_{i,j,k} p_{{{\text{o}}i,j,k - 1}}^{n + 1} + g_{i,j,k} p_{{{\text{o}}i,j,k}}^{n + 1} + b_{i,j,k} p_{{{\text{o}}i + 1,j,k}}^{n + 1} + d_{i,j,k} p_{{{\text{o}}i,j + 1,k}}^{n + 1} + f_{i,j,k} p_{{{\text{o}}i,j,k + 1}}^{n + 1} = h_{i,j,k}$$
(20)

with

$$a_{i,j,k} = \frac{{T_{{{\text{w}}i{ - }\frac{1}{2},j,k}} }}{{\rho_{{{\text{w}}i,j,k}} }} + \frac{{T_{{{\text{o}}i{ - }\frac{1}{2},j,k}} }}{{\rho_{{{\text{o}}i,j,k}} }},$$
$$b_{i,j,k} = \frac{{T_{{{\text{w}}i{ + }\frac{1}{2},j,k}} }}{{\rho_{{{\text{w}}i,j,k}} }} + \frac{{T_{{{\text{o}}i{ + }\frac{1}{2},j,k}} }}{{\rho_{{{\text{o}}i,j,k}} }},$$
$$c_{i,j,k} = \frac{{T_{{{\text{w}}i,j{ - }\frac{1}{2},k}} }}{{\rho_{{{\text{w}}i,j,k}} }} + \frac{{T_{{{\text{o}}i,j{ - }\frac{1}{2},k}} }}{{\rho_{{{\text{o}}i,j,k}} }},$$
$$d_{i,j,k} = \frac{{T_{{{\text{w}}i,j{ + }\frac{1}{2},k}} }}{{\rho_{{{\text{w}}i,j,k}} }} + \frac{{T_{{{\text{o}}i,j{ + }\frac{1}{2},k}} }}{{\rho_{{{\text{o}}i,j,k}} }},$$
$$e_{i,j,k} = \frac{{T_{{{\text{w}}i,j,k{ - }\frac{1}{2}}} }}{{\rho_{{{\text{w}}i,j,k}} }} + \frac{{T_{{{\text{o}}i,j,k{ - }\frac{1}{2}}} }}{{\rho_{{{\text{o}}i,j,k}} }},$$
$$f_{i,j,k} = \frac{{T_{{{\text{w}}i,j,k\user2{ + }\frac{{1}}{{2}}}} }}{{\rho_{{{\text{w}}i,j,k}} }} + \frac{{T_{{{\text{o}}i,j,k\user2{ + }\frac{1}{2}}} }}{{\rho_{{{\text{o}}i,j,k}} }},$$
$$g_{i,j,k} = - \left[ {a_{i,j,k} + b_{i,j,k} + c_{i,j,k} + d_{i,j,k} + e_{i,j,k} + f_{i,j,k} + {{\left( {V_{{\text{p}}} c_{{\text{t}}} } \right)_{i,j,k} } \mathord{\left/ {\vphantom {{\left( {V_{{\text{p}}} c_{{\text{t}}} } \right)_{i,j,k} } {\Delta t}}} \right. \kern-\nulldelimiterspace} {\Delta t}}} \right],$$
$$h_{i,j,k} = - \frac{{\Delta \left( {T_{{\text{w}}} \Delta p_{{{\text{cow}}}} } \right)_{i,j,k} + \Delta \left( {T_{{\text{w}}} \gamma_{{\text{w}}} \Delta D} \right)_{i,j,k} + \delta_{{{\text{w}}i,j,k}} Q_{{{\text{w}}i,j,k}} }}{{\rho_{{{\text{w}}i,j,k}} }} - \frac{{\delta_{{{\text{o}}i,j,k}} Q_{{{\text{o}}i,j,k}} + \Delta \left( {T_{{\text{o}}} \gamma_{{\text{o}}} \Delta D} \right)_{i,j,k} }}{{\rho_{{{\text{o}}i,j,k}} }} - \frac{{\left( {V_{{\text{p}}} c_{{\text{t}}} } \right)_{i,j,k} }}{\Delta t}p_{{{\text{o}}i,j,k}}^{n} ,$$
$$T_{{{\text{l}}i \pm \frac{{1}}{2},j,k}} = \frac{{\left( {AK} \right)_{{i \pm \frac{{{1}}}{{{2}}},j,k}} }}{{{\Delta x}_{{i \pm \frac{{{1}}}{{{2}}},j,k}} }}\left( {\frac{{K_{{{\text{rl}}}} \rho_{{\text{l}}} }}{{\mu_{{\text{l}}} }}} \right)_{{i \pm \frac{{1}}{{2}},j,k}} ,$$
$$T_{{{\text{l}}i,j \pm \frac{{{1}}}{{{2}}},k}} = \frac{{\left( {AK} \right)_{{i,j \pm \frac{{{1}}}{{{2}}},k}} }}{{{{\Delta y}}_{{i,j \pm \frac{{{1}}}{{{2}}},k}} }}\left( {\frac{{K_{{{\text{rl}}}} \rho_{{\text{l}}} }}{{\mu_{{\text{l}}} }}} \right)_{{i,j \pm \frac{{{1}}}{{{2}}},k}} ,$$
$$T_{{{\text{l}}i,j,k \pm \frac{{{1}}}{{{2}}}}} = \frac{{\left( {AK} \right)_{{i,j,k \pm \frac{{{1}}}{{{2}}}}} }}{{{\Delta z}}_{{i,j,k \pm \frac{{{1}}}{{{2}}}}} }\left( {\frac{{K_{{{\text{rl}}}} \rho_{{\text{l}}} }}{{\mu_{{\text{l}}} }}} \right)_{{i,j,k \pm \frac{{1}}{{2}}}} .$$

The difference equation of the oil phase is as follows:

$${\Delta}\left( {T_{{\text{o}}}^{n}{\Delta p}}_{{\text{o}}}^{n + 1} \right)_{i,j,k} + \delta_{{\text{o}}}^{n} Q_{{{\text{o}}i,j,k}}^{n} = \frac{1}{{{\Delta t}}}\left[ {\left( {V_{{\text{p}}} S_{{\text{o}}} \rho_{{\text{o}}} } \right)^{n + 1} - \left( {V_{{\text{p}}} S_{{\text{o}}} \rho_{{\text{o}}} } \right)^{n} } \right]_{i,j,k}$$
(21)

Based on Eq. (21), the oil saturation is achieved, and then, the water saturation is obtained by Eq. (4).

4.2 Validation of the numerical model

In order to verify the reliability of the simulation model, we choose a well of J Oilfield in China to simulate high-pressure soaking and depressurizing production. J Oilfield is a typical terrestrial reservoir. The reservoir depth is about 2000–4000 m with an effective thickness of 10–30 m. Natural fractures are not well developed, and the dissolved gas–oil ratio is very low. The formation pressure gradient is around 1.0–1.3 MPa/100 m. The permeability is 0.001–0.1 mD, and the oil viscosity is several to tens centipoises. The wettability is neutral-wet to slightly oil-wet, and water sensitivity is weak. Multistage fracturing was performed in horizontal well with cement plug. The development effect of depletion drive is not ideal in the pilot. Therefore, high-pressure soaking with fracturing fluid after hydraulic fracturing has been conducted in this oilfield.

JHW023 is the pilot well for high-pressure soaking in J Oilfield, and the oil production is much higher than the previous wells. The length of the horizontal well is 1,246 m, and it is fractured with 37,408 m3 of fracturing fluid. The fracturing rate reaches up to 14–15 m3/min. The fracturing pressure is around 70–75 MPa. After fracturing, it was soaked for about 60 days under high-pressure conditions, then different chokes were employed for depletion drive. Oil and water production rates are recorded; thus, constant liquid production rates can be used to match the oil production rates and water-cut, as shown in Fig. 3. The parameters of JHW023 and oil reservoir are shown in Table 1. The reservoir and fluid parameters come from the oilfield and laboratory. The parameters of hydraulic fractures, relative permeability, and compressibility are obtained by matching the practical data. During the simulation process, about 37,400 m3 water were injected into the reservoir by a horizontal well with 81 cluster of hydraulic fractures at first. Then, high-pressure soaking were followed for about 60 days. After that, depletion drive was carried out. The BHP and production data of high-pressure soaking and depletion drive stages were monitored and used to match the simulation results. Figures 4 and 5 show comparisons of oil production rate, water-cut, and BHP between field data and simulation results. We can see that the simulation results have a good agreement with the field data. Figure 6 shows the flowback percentage of the injected water. It is clear that only 20% of the fracturing fluid flows back after producing for six months, which agrees well with the typical data found in literature (King 2010). We also used the fitted model to study the contributions of imbibition and elastic energy. Figure 7 demonstrates the cumulative oil production of the scenarios with/without imbibition. Based on the data, the contributions of elasticity and imbibition are about 20% and 80%, respectively.

Fig. 3
figure 3

The specific liquid production rates and choke sizes at different times

Table 1 Parameters used for simulation of Well JHW023 in J Oilfield
Fig. 4
figure 4

Comparisons of oil production rate and water-cut between field data and simulation

Fig. 5
figure 5

Comparisons of bottom hole pressure between field data and simulation

Fig. 6
figure 6

Backflow percentage of the injected fracturing fluid

Fig. 7
figure 7

Cumulative oil productions of different scenarios

5 Results and discussion

5.1 Model applications

On the basis of oilfield production data, the injected water was retained in the reservoir. Smaller flowback percentage of the fracturing fluid indicates a better production (Abbasi et al. 2012; Williams-Kovacs et al. 2015; Cao et al. 2017). As mentioned, high-pressure soaking plays an important role in the fracturing fluid backflow. In order to clearly understand the significances of high-pressure soaking, the contributions of different mechanisms, the optimal soaking times, and the huff-n-puff performances, an ideal model based on the average values in J Oilfield is employed to conduct simulations.

The formation thickness is 15 m. The fracture spacing and fracture half-length used in the simulation are 40 m and 100 m, respectively. Therefore, a typical unit was taken from the fractured tight reservoir. Because the fracture conductivity is very high, the distribution of the injected or produced fluid in fractures is almost uniform. Therefore, the gridding system and well placement are as shown in Fig. 8. The initial reservoir pressure is 30 MPa, and the fracturing pressure is 70 MPa, which means the initial pressure in fractures is 70 MPa after fracturing. The initial oil saturation in matrix and fractures are 0.70 and 0.25, respectively. We assume that the matrix block is cut into a cuboid by the fractures, which are full of the fracturing fluid or surfactant solution. The maximum BHP is set as 70 MPa for injecting water, and the BHP is set as 15 MPa for production. The other parameters used in the simulation cases are from the fitted model in Sect. 4 and listed in Table 2. The parameters of sensitivity analysis are shown in Table 3.

Fig. 8
figure 8

Schematic diagram of the gridding system and well placement in the typical unit

Table 2 Parameters used in simulation cases
Table 3 Parameters for operational sensitivity

5.1.1 Significances of high-pressure soaking

In this section, we aimed to observe the significances of high-pressure soaking by comparing the production performances, pressure propagation, and fluid movement in different scenarios. Because the wettability of most tight rocks is mixed-wet or slightly oil-wet, the performances of high-pressure soaking for both water-wet and mixed-wet scenarios are compared as well. Figure 9 shows the comparisons of production performances of water-wet rock without soaking, water-wet rock with soaking, and mixed-wet with soaking. From Fig. 9a, the oil recovery of the rock without soaking is 1.82%, and it is the lowest. The oil recovery for mixed-wet rock with high-pressure soaking for 60 days is 2.88% and that for water-wet rock with high-pressure soaking for 60 days is 4.02%. As the results indicate, high-pressure soaking turns out very significant. For the water-wet rock without soaking, there is no enough time for imbibition and pressure propagation, and thus, less oil is expelled by imbibition and elastic energy. For the mixed-wet rock with high-pressure soaking, despite no imbibition occurring due to the adverse wettability, there is enough time for pressure propagation and more oil is expelled at depressurizing stage. For the water-wet rock with high-pressure soaking, both imbibition and pressure propagation play a part in increasing elastic energy, so the most oil is expelled. Figure 9b shows the water-cut and water backflow in different scenarios. We can see that the water-cut is as high as 100% initially for the case without soaking, and the water backflow percentage reaches about 33%. However, if it is soaked for 60 days, the water-cut is lower and sharply reduces. The backflow percentages of mixed-wet and water-wet rocks are 30.5% and 22.4%, respectively. Thus again, more water can enter the matrix to expel the oil out by high-pressure soaking.

Fig. 9
figure 9

Comparisons of production performances with and without soaking for 60 days

Despite the mechanisms can be speculated from laboratory and field experiments, we also demonstrated the visual comparisons of pressure and saturation in the whole process. Figure 10 indicates the variations of BHP during high-pressure soaking and depressurizing production. It is clear that BHP decreases rapidly at the initial stage of soaking. The balanced pressure, about 58 MPa, is reached after soaking for about 10 days. The variations of pressure distribution in the matrix from the initial state to depressurizing production stage are shown in Fig. 11. The pressure propagation mainly occurs in the first 8 days. After high-pressure soaking, the pore pressure in the matrix is much higher than the initial reservoir pressure, so the elastic energy significantly increases.

Fig. 10
figure 10

Variations of BHP at the stages of soaking and depressurizing production

Fig. 11
figure 11

Variations of the pressure distributions in the whole process (unit: MPa)

Figure 12 shows the variations of oil saturation in different regions during high-pressure soaking. The oil saturation always increases in fractures, but decreases rapidly in the matrix surface firstly and then changes slightly. The water saturation in the matrix center mainly changes on the side. Figure 13 demonstrates the variations of the mean water saturation in the fracture and matrix surface. From Fig. 13a, we can see that the variations of the mean water saturation show different features in fractures at different stages. At the high-pressure soaking stage, the water saturation gradually decreases because water enters the matrix to expel out the oil by imbibition. At the depressurizing production stage, the mean water saturation decreases sharply and then slowly because the water in fractures is easier to flow into the well. Because water enters the matrix under the high-pressure gradient, the mean water saturation increases fast in the matrix surface at first, which is corresponding to the pressurizing stage in Figs. 1b, 10, and 11. Then, the mean water saturation slowly increases during soaking, which is corresponding to high-pressure soaking in Fig. 1c. Afterwards, the mean water saturation shows a small jump at the beginning of depressurizing production and then does not change anymore. The “jump” indicates more oil flows out than water, which is corresponding to the analysis of depressurizing stage in Fig. 1d. These visible maps and statistic data well support the previous inferences.

Fig. 12
figure 12

Variations of the distribution of oil saturation in different locations during high-pressure soaking

Fig. 13
figure 13

Variations of the mean water saturation in fractures and the matrix surface

5.1.2 Contributions of different mechanisms

As the previous analysis, elastic energy and imbibition are the main mechanisms of oil production for high-pressure soaking. In the above simulations, we preliminarily analyzed the contributions of elastic energy and imbibition. In this section, we particularly investigated the impacts of the acting forces on oil production. The acting forces during elastic drive include gravity and elastic force. The acting forces during imbibition include gravity and capillary force. Therefore, four cases were conducted. The contribution of elastic force caused by pressure change is observed in Case 1. The contribution of both gravity and elastic force are included in Case 2, so the contribution of gravity could be obtained from the difference between Cases 1 and 2. The contribution of both capillary force and elastic force is included in Case 3, so the contribution of capillary force could be obtained from the difference between Cases 1 and 3. The contributions of both elastic drive and imbibition are included in Case 4 by considering the gravity, capillary, and elastic force, so the contribution of gravity to both elastic drive and imbibition could be obtained from the difference between Cases 3 and 4. Correspondingly, the contribution of gravity to imbibition can be observed from the differences among Cases 1–4. Figure 14 shows the cumulative oil production in different scenarios. We can see that the cumulative oil production of Case 2 is lower than that of Case 1, which indicates that gravity plays a slightly negative role in elastic drive. The reason is that the gravitational differentiation makes the water separate from oil in the fractures, and water is much easier to flow into the wellbore. The cumulative oil production of Case 3 is much larger than that of Case 1, which means the capillary force has a significant contribution through imbibition. Case 4 is very close to Case 3, which seemingly means gravity has no impact on high-pressure soaking. In fact, the positive role of gravity in oil production by imbibition counteracts its negative role in oil production by elastic drive. In summary, the capillary force plays a dominant role in imbibition, pressure dominates elastic drive, and gravity has slight impacts on both imbibition and elastic drive.

Fig. 14
figure 14

Cumulative oil productions of considering different mechanisms

Henceforth, we studied the impacts of oil viscosity, BHP, IFT, matrix permeability, and specific surface area on the contributions of different mechanisms, which are shown in Figs. 15, 16, 17, 18 and 19. From Fig. 15, oil viscosity has a slight impact on the contribution, since the increase in oil viscosity is unfavorable to both imbibition and elastic drive. Generally, the contribution of elastic energy increases as the oil viscosity increases. From Fig. 16, the contribution of elastic energy significantly decreases as the BHP increases. Elastic energy is dominated by the pressure difference, which has less impact on imbibition. As observed in Fig. 17, the contribution of elastic energy decreases as the IFT increases. The reason is that when the IFT is low, the capillary force is extremely low, so imbibition is very weak. As observed in Fig. 18, the contribution of elastic energy first increases as the matrix permeability increases and then decreases. Generally, the matrix permeability has a slight impact on the contribution. When the permeability is low, both the capillary force and the flow resistance are high, so imbibition contributes more than elastic energy. However, as the permeability increases, the capillary force decreases, so imbibition contributes less than elastic energy. From Fig. 19, the contribution of imbibition increases as the specific surface area increases. That is because larger specific surface area provides relatively larger contact area to conduct imbibition.

Fig. 15
figure 15

The impacts of oil viscosity on the contribution of different mechanisms

Fig. 16
figure 16

The impacts of BHP on the contribution of different mechanisms

Fig. 17
figure 17

The impacts of IFT on the contribution of different mechanisms

Fig. 18
figure 18

The impacts of matrix permeability on the contribution of different mechanisms

Fig. 19
figure 19

The impacts of specific surface area on the contribution of different mechanisms

5.1.3 Optimal soaking time

According to the above simulations, the pressure propagation and imbibition occur during high-pressure soaking. They decide the oil production at the depressurizing stage. Figure 20 shows the oil production for different soaking times. We can see that the oil production increases to a plateau after producing for several days. Both the rising rate and the oil production increase as the soaking time increases. Despite longer soaking time works better, the oil production increases slowly when the soaking time gets longer. Therefore, there must exist an optimal value of soaking time for highly efficient development. We take 10 days as a step and compare the oil production. If the soaking time is increased by 10 days, while the yield is increased by less than 5%, the best soaking time is obtained.

Fig. 20
figure 20

Oil production for different soaking times

Figure 21 shows the impact of oil viscosity on the optimal soaking time. The relation between the optimal soaking time and oil viscosity is as follows:

$$t_{{{\text{soak}}}} = 3.5\mu_{{\text{o}}} + 37.5.$$
(18)
Fig. 21
figure 21

The impacts of oil viscosity on the optimal soaking time

Figure 22 shows the impact of BHP on the optimal soaking time. The optimal soaking time linearly changes as the BHP increases, and the relation is as follows:

$$t_{{{\text{soak}}}} = 2p_{{{\text{wf}}}} + 30.$$
(19)
Fig. 22
figure 22

The impacts of bottom hole pressure on the optimal soaking time

Figure 23 shows the impact of IFT on the optimal soaking time. When the IFT is lower, imbibition is very weak. It just takes little time to reach the balanced pressure, which is about 10 days. As the IFT increases, imbibition gradually strengthens, so it needs a longer time to conduct imbibition. A logistic function can be used to approximately describe their relationship:

$$t_{{{\text{soak}}}} = \frac{65}{{1{ + }\exp \left( {1.875 - {{\sigma_{{{\text{ow}}}} } \mathord{\left/ {\vphantom {{\sigma_{{{\text{ow}}}} } 4}} \right. \kern-\nulldelimiterspace} 4}} \right)}}.$$
(20)
Fig. 23
figure 23

The impacts of interfacial tension on optimal soaking time

Figure 24 shows the impact of matrix permeability on the optimal soaking time. We can see that the optimal soaking time logarithmically declines as the permeability increases:

$$t_{{{\text{soak}}}} = 40.4 - 6.82\ln K.$$
(21)

According to Eq. (11), when the permeability is lower, the capillary force is larger. Imbibition plays an important role on oil production, so it takes a longer soaking time. As shown in Fig. 25, the optimal soaking time also logarithmically declines as the specific surface area (As) increases:

$$t_{{{\text{soak}}}} = 20.4 - 24.9\ln A_{{\text{s}}} .$$
(22)
Fig. 24
figure 24

The impacts of matrix permeability on optimal soaking time

Fig. 25
figure 25

The impacts of specific surface area on optimal soaking time

Therefore, the tight formations with conventional fracturing need a longer soaking time than that of the ones with stimulated reservoir volume (SRV).

According to the impacts of different factors on the soaking time, a general model was achieved to calculate the optimal soaking time for high-pressure soaking including the above factors:

$$t_{{\text{o}}} = 3.18\mu_{{\text{o}}} - 39.56k + 2.23\sigma_{{{\text{ow}}}} + 2p_{{{\text{wf}}}} - 34.9A_{{\text{s}}} - 28.$$
(23)

The agreement between the predicted results and experimental values is shown in Fig. 26.

Fig. 26
figure 26

The agreement between experimental and predicted values of optimal soaking times

5.1.4 Performances of huff-n-puff

After soaking by the fracturing fluid and depressurizing production, huff-n-puff would be carried out. Some simulations were performed to observe the huff-n-puff performances. The first cycle is hydraulic fracturing—high pressure soaking—depressurizing production. After that, water is injected until the bottom hole pressure reaches 70 MPa, and the soaking is performed for a period of time. After soaking, the well produces at a constant bottom hole pressure of 15 MPa. The cyclic soaking and production times were optimized one cycle after another, and the flowchart is shown in Fig. 27.

Fig. 27
figure 27

The flowchart of optimizing the optimal soaking and production times for each cycle

Figure 28 shows the oil production performances of huff-n-puff. We can see that both the oil production rate and the cycle oil production decrease from the first cycle to the last. The reason is that the fluid redistributes under the mechanisms of imbibition and elastic energy in each cycle. It can be explained by Fig. 1. After the first cycle, the distribution of fluid is shown in Fig. 1d. In the high-pressure transmission stage of the second cycle, the injected water enters from both ends of the capillary, and the length of water column increases. The difference between Pc1 and Pc2 gradually decreases to zero, so less fluid or oil could be driven out at forced imbibition stage. Correspondingly, less oil could be expelled from the wide end at the depressurizing stage. In addition, after the forward cycles, more water enters the ends of the capillary, so water will be produced firstly at the stage of the elastic energy release. For the above reasons, cycle oil production decreases from the first cycle to the last.

Fig. 28
figure 28

Oil production performances of huff-n-puff

Figure 29 shows the oil productions and the optimal soaking times in different cycles. The optimal soaking times gradually decrease from 60 days in the first cycle to 0 day after eight cycles, but the optimal production times keep for 30 days in each cycle. The reason is that, on one hand, imbibition is proceeding in the early production process, so the soaking time gradually decreases; on the other hand, since the difference between Pc1 and Pc2 decreases, the imbibition equilibrium is much easier to be reached. Therefore, a longer soaking time should be chosen at the first several cycles, and more oil can be recovered. Combining with Fig. 28, the economic benefit should be evaluated for the lowly periodic oil production after several cycles.

Fig. 29
figure 29

Oil productions and optimal soaking times in different cycles

Figure 30 shows the water production/injection performances of huff-n-puff. All of the water-cut, cycle water production, and cycle water injection are the largest in the first cycle (hydraulic fracturing and backflow stage). Both the water-cut and cycle water production gradually increase from the second cycle, but the cycle water injection is very close. The reason is that the fractures are full of water after hydraulic fracturing and part of the oil is flooded far away the wellbore, so the water-cut is high at first. After the injected water flows back, an amount of oil in the fractures and less water are produced, so the water-cut decreases. From the second cycle, the oil in fractures gradually decreases and cycle water injection is close, so the water-cut increases.

Fig. 30
figure 30

Water production/injection performances of huff-n-puff

The optimal soaking and production times in different scenarios were also studied. Figure 31 shows the impacts of matrix permeability and specific area on the optimal soaking and production times. It is clear that the optimal soaking time gradually decreases and the optimal production time does not change in different scenarios. The optimal soaking and production times decrease as the specific area increases, but increase as the matrix permeability decreases. Therefore, we can choose a shorter soaking and production time for the reservoirs with high permeability or stimulated reservoir volume (SRV).

Fig. 31
figure 31

The impacts of matrix permeability and specific area on the optimal soaking and production times

5.2 Limitations of our study

Although we have demonstrated a significant study of high-pressuring soaking, the real process is more complicated. In our study, phase entrapment during leak-off and flow-back is neglected. The injected water may have a negative influence on tight formations due to permeability reduction, which is attributed to water sensitivity. Therefore, the optimal time will slightly decrease in such circumstance. In addition, soaking is also important for recovering the reservoir temperature near the wellbore after hydraulic fracturing, especially for the crude oil with high viscosity or wax content.

6 Conclusions

A whole process of pressurizing, high-pressure soaking, and depressurizing was analyzed to understand the capacity and mechanisms of high-pressure soaking after hydraulic fracturing. A mechanistic model was established and verified. Finally, a simulation method was employed to address the concerns about the mechanisms and capacity of high-pressure soaking. The key findings are summarized as follows,

  • High-pressure soaking has a significant impact on oil production through the mechanisms of imbibition and elasticity, both of which have considerable contributions. Imbibition works at the soaking stage, and elasticity works at the depressurizing stage.

  • High-pressure soaking could make more water enter the matrix even in mixed-wet rocks. About 20% of the fracturing fluid flows back after adequate soaking in water-wet rocks. For water-wet or mixed-wet rocks, the flowback rate after soaking is less than that without soaking.

  • The contribution of imbibition increases as the increase in BHP, IFT, and specific surface area. However, it slightly decreases as the oil viscosity increases, and first decreases and then slightly increases as the increase in matrix permeability.

  • The optimal soaking time increases linearly as the oil viscosity and BHP increase. It logarithmically declines with the increase in matrix permeability and specific surface area. Moreover, it shows a rising trend as the IFT increases.

  • During huff-n-puff, cyclic oil production decreases from the first cycle to the last one. The optimal soaking time of each cycle gradually decreases, but the optimal production time is unchanged. The optimal soaking and production times decrease as the specific area increases, but increase as the matrix permeability decreases.

  • Both water production and injection in the first cycle are the largest. The cyclic water cut and water production gradually increase from the second cycle, but the cyclic water injections are very close.