Exact Solution for Non-Self-Similar Wave-Interaction Problem during Two-Phase Four-Component Flow in Porous Media

and Applied Analysis 3 the fractional flow of the aqueous phase is affected by concentration of the dissolved chemical, the fractional flow is independent of salt concentration, chemical and salt do not dissolve in oil, linear sorption for the polymer a = Γc, Henry’s sorption coefficient Γ is salinity-dependent, salt does not adsorb on the rock, and temperature is constant. The system of governing equations consists of mass balance equations for aqueous phase, for dissolved and adsorbed chemical, and for dissolved salt [8, 9]: ∂s ∂t + ∂f (s, c) ∂x = 0 (4)


Introduction
Exact self-similar solutions of Riemann problems for hyperbolic systems of conservation laws and non-self-similar solutions of hyperbolic wave interactions have been derived for various flows in gas dynamics, shallow waters, and chromatography (see monographs [1][2][3][4][5][6][7][8]). For flow in porous media, hyperbolic systems of conservation laws describe twophase multicomponent displacement [9,10] where s is the saturation (volumetric fraction) of aqueous phase and f is the water flux.Equation ( 1) is the mass balance for water and (2) is the mass balance for each component in the aqueous solution.Under the conditions of thermodynamic equilibrium, the concentrations of the components adsorbed on the solid phase (a i ) and dissolved in the aqueous phase (c i ) are governed by adsorption isotherms:  =  () , =( 1 ,  2 , . . .,   ) , =( 1 ,  2 , . . .,   ) .
Exact and semianalytical solutions of one-dimensional flow problems are widely used in stream-line simulation for flow prediction in three-dimensional natural reservoirs [10].The sequence of concentration shocks in the one-dimensional analytical solution is important for interpretation of laboratory tests in two-phase multicomponent flow in natural reservoir cores.The scalar hyperbolic equations ( 1) and ( 2),  = 0, correspond to displacement of oil by water [9,10].The ( + 1) × ( + 1) system (1) and (2) describes two-phase flow of oleic and aqueous phases with  components (such as polymer and different salts) that may adsorb and be dissolved in both phases.These flows are typical for so-called chemical enhanced oil recovery displacements, like injections of polymers or surfactants, and for numerous environmental flows [9,10].For polymer injection in oil reservoirs,  = 1 corresponds to polymer and  = 2, 3, . . .,  to different ions.
Therefore the system (1) and ( 2) is called the multicomponent polymer-flooding model [11,12].Besides, (−1)×(−1) hyperbolic system (1) and (2) describes two-phase -component displacement, which is typical for so-called gas methods of enhanced oil recovery [9,10,13,14].The processes of hot water injection with phase transitions, secondary migration of hydrocarbons with consequent formation of petroleum accumulations, enhanced geothermal energy projects, and injections into aquifers are described by the above systems.The Riemann problems correspond to continuous injection of chemical solutions or gases into oil reservoirs; the solutions are self-similar [3,9,14].The wave-interaction problems correspond to piece-wise-constant initial-boundary conditions, for which the solutions are non-self-similar [1,10,[15][16][17].The wave-interaction solutions describe injection of limited slugs (banks) of chemical solutions or gaseous solvents driven in the reservoirs by water or gas [9,10].
The system (1) and ( 2) describes two-phase multicomponent displacements in large scale approximation, where the dissipative effects of capillary pressure, diffusion, and thermodynamic nonequilibrium are negligible if they are compared with advective fluxes under the large length scale of the natural subterranean reservoirs.Travelling waves near to shock discontinuities in dissipative systems have been presented in [10,32].A semianalytical global solutions have been obtained by Geiger et al. [33] and Schmid et al. [34]; see also [16].
The particular case of so-called multicomponent polymer flooding is the dependency of the component sorption concentration of its own concentration only   ( 1 ,  2 , . . .,   ) =   (  ).Exact solutions of the Riemann problem for this case show that the concentration of each component performs the jump without shocks of other components (see the corresponding solution in the books [10,21]).Therefore, in concentration profiles, the shocks are located in order of decrease of derivatives of the sorption functions.In the case of Henry isotherms   (  ) = Γ    , the shocks are located in order of increase of Henry's sorption coefficients Γ  .
The distinguished invariant feature of ( + 1) × ( + 1) conservation law systems for two-phase multicomponent flows in porous media with sorption and phase transitions equations (1) and ( 2) is its splitting into an  ×  auxiliary system for concentrations   (, ) and a scalar hyperbolic equation for saturation (, ) [35,36].This splitting explains the simple form of Riemann problem solutions for system (1) and (2) as compared with gas dynamics or chromatography [1,2,37].
The non-self-similar solution of system (1) and ( 2),  = 2, for slug injections has been considered by Fayers [17], where the qualitative behaviour of characteristic lines and shocks has been described.The exact solutions of system (1) and ( 2) for  = 2 and 3 have been obtained in [15] (see book [10] for detailed derivations, in which the sorption of component depends on its own concentration only   =   (  ),  = 1, 2, . ..).Numerous interactions of different saturation-concentration shocks occur after the injection, resulting in appearance of moving zones with different combinations of components.However, after all interactions, different component slugs are separated from each other.As in the case of continuous injection, the slugs are finally positioned in the order of decreasing sorption isotherm derivatives (  /  ).It seems that this simplified case draws the line where the analytical solutions can be found from the analysis of system (1) and ( 2) directly.Consideration of crosseffects   =   ( 1 ,  2 , . . .,   ) in sorption functions equation (3) introduces significant difficulties into wave analysis, and even the Riemann problem cannot be solved for any arbitrary case  = 2 (see [38], where the Riemann solutions have been obtained for several particular cases).
The splitting technique reduces number of equations in (1) and ( 2) by one, allowing for exact solutions in more complex multicomponent cases [35][36][37][38][39][40].The Riemann problem with cross-effects for adsorption   =   ( 1 ,  2 ) has been solved in [39,41,42] for continuous polymer injection with varying salinity using the splitting method.In the current paper, the exact solution for non-self-similar problem of injection of polymer slug with varying salinity followed by water drive is obtained.
The structure of the text is as follows.The particular case of the general system (1) and ( 2) that is discussed in the current work is introduced in Section 2 along with physics assumptions and initial-boundary conditions for slug injection problem.The detailed description of the splitting procedure for the system is discussed along with formulation of initial and boundary conditions for the auxiliary system which is presented in Section 3. Section 4 contains derivation of the Riemann solution that corresponds to the first stage of the slug injection.The wave-interaction slug injection problem is solved in Section 5. Section 6 contains a simplified solution for the particular case where the initial chemical concentration is zero, which corresponds to the case of polymer slug injection.The paper is concluded by physical interpretation of the solution obtained for chemical slug injection with different water salinity into oilfield (Sections 7 and 8).

Formulation of the Problem
Let us discuss the displacement of oil by aqueous chemical solution with water drive accounting for different salinities of formation and injected waters.In the following text, the component  = 1 is called the polymer, and that  = 2 is called the salt.The assumptions of the mathematical model are as follows: both phases are incompressible, dispersion and capillary forces are neglected, there are two phases (oleic and aqueous phases) and two components dissolved in water (polymer and salt), water and oil phases are immiscible, chemical and salt concentrations in water are negligibly small and do not affect the volume of the aqueous phase, the fractional flow of the aqueous phase is affected by concentration of the dissolved chemical, the fractional flow is independent of salt concentration, chemical and salt do not dissolve in oil, linear sorption for the polymer  = Γ, Henry's sorption coefficient Γ is salinity-dependent, salt does not adsorb on the rock, and temperature is constant.
The system of governing equations consists of mass balance equations for aqueous phase, for dissolved and adsorbed chemical, and for dissolved salt [8,9]: where  is the water saturation,  is the fractional flow function,  is the polymer sorption isotherm, and  and  are chemical and salt concentrations, respectively.The fractional flow function (water flux) depends on the water saturation  and on the chemical concentration .The typical S-shapes of fractional flow functions  under  = const are shown in Figure 1.The fractional flow is a monotonically decreasing function of .Sorption isotherms are linear for fixed salinity (, ) = Γ().The functions and  are assumed to be bounded and smooth.
The displacement of oil by chemical slug corresponds to the following initial-boundary problem: For  < 1, during continuous injection of chemical solution with different salinity, the solution of system ( 4)-( 6) subject to initial-boundary conditions equations ( 7) and ( 8) coincides with the solution of the Riemann problem: The initial condition is denoted by  in Figure 1 and the boundary condition corresponding to injection of the slug is denoted as .
Generally (, 0) =  2 > 0 is positive.Further in the text, the component with concentration  is called "chemical, " while for the case of the absence of this component initially in the reservoir (, 0) =  2 = 0 we use the term "polymer." The solution of the Riemann problem is self-similar: (, ) = (), (, ) = (), (, ) = (),  = / and it can be found in [37,39,40].The solution of the problem (7) and (8) in the neighbourhood of the point (0, 1) in (, )plane is also self-similar.The global solution of the system (4)-( 6) subject to the initial-boundary conditions equations Figure 1: Fractional flow curves and Riemann problem solution, where  1 is the slug concentration,  2 is the initial concentration, and  * is the intermediate concentration.(7) and ( 8) is non-self-similar; it expresses the interactions between hyperbolic waves occurring from decays of Riemann discontinuities in points (0, 0) and (0, 1) in (, )-plane.
System of ( 4)-( 6) subject to the initial and boundary conditions equations ( 9) and ( 10) is solved in Section 4 using the method so-called splitting procedure [35].This procedure is explained in the next section.

Splitting Procedure
In the present section we briefly explain the splitting method for the solution of hyperbolic system of conservation laws equations ( 4)-(6).

Streamline/Potential Function and Auxiliary System.
As it follows from divergent (conservation law) form of equation for mass balance for water (1) or (4), there does exist such a potential function (, ) that that is, Equation ( 4) is merely the condition of equality of second derivatives of the potential  as taken in different orders.
It also expresses that the differential of the first order form equation ( 12) equals zero.The splitting procedure consists of changing the independent variables from (, ) to (, ) in system ( 4)-( 6).Figures 2 and 3 show the corresponding mapping [43,44].
From fluid mechanics point of view, (, ) is a potential function, which equals the volume of fluid flowing through f L (x = 0, t 1 ), c L (x = 0, t 1 ) f L (x = 0, t 1 ), c L (x = 0, t 1 )   a trajectory connecting points (0, 0) and (, ).As it follows from (12), two streamlines in Figure 3 correspond to constant values of potential, that is, there is no flux through streamlines: Equation ( 4) shows that the integral of ( 13) along the closed contour is equal to zero; that is, the volume of fluid flowing through a trajectory connecting points (0, 0) and (, ) is independent of trajectory and depends on end points only.The potential function equation ( 13) is determined in the way that one end of trajectory is fixed at point (0, 0).Let us derive the relationship between the elementary wave speeds of the system in (, ) coordinates and those of the large system in (, ).Consider the trajectory  =  0 () and its image  =  0 () by the mapping equation ( 13): Define the trajectory speeds as Let us use  as a parameter for both curves  =  0 () and  =  0 ().Taking derivation of both parts of (13) in  along trajectories and using speed definitions in (16), we obtain from which follows the relationship between elementary wave speed in planes (, ) and (, ): For example, the eigenvalues of the system of equation in (, ) plane   and in (, ), Λ  , are related by (Figure 4 [43,44]) From now on, the independent variables (, ) are used in ( 4)-( 6) instead of (, ).Expressing the differential form  from (12) as and accounting for zero differential of the form we obtain the expression for (4) in coordinates (, ) So, ( 22) is the mass balance for water; that is, it is (4) rewritten in coordinates (, ).
Let us derive ( 5) in (, ) coordinates.The conservation laws for (5) in the integral form are where Ω is a closed domain Ω ⊂  2 , so the integral of ( 23) is taken over the closed contour.
Applying the definition of the potential function equation ( 13) into (23) yields Tending the domain radius to zero and applying Green's theorem, Now let us perform change of independent variables in ( 6) in (, ) coordinates as follows: Finally, the ( + 1) × ( + 1) system of conservation laws for two-phase n component chemical flooding in porous media with adsorption can be split into an  ×  auxiliary system equations ( 25) and ( 26) and one independent lifting equation (22).The splitting is obtained from the change of independent variables (, ) to (, ).This change of coordinates also transforms the water conservation law into the lifting equation.
The auxiliary system contains only equilibrium thermodynamic variables, while the initial system contains both hydrodynamic functions (phase's relative permeabilities and viscosities) and equilibrium thermodynamic variables.
The above splitting procedure is applied to the solution of displacement of oil by polymer slug with alternated salinity in the next section.

Formulation of the Splitting Problem for Two-Phase Flow with Polymers and Salt. Introducing new variables "density"
and "flux"  and applying the splitting technique, the 3 × 3 system (4)-( 6) is transformed to the following form: ( (, ))  +  ()  = 0 (28) The auxiliary system equations ( 29) and (30) are independent of (28).The auxiliary system has thermodynamic nature since it contains only sorption function (, ) and the unknowns are the component concentrations c and .Equation ( 28) is the volume conservation for two immiscible phases.For the known auxiliary solution of ( 29) and (30), equation ( 28) is a scalar hyperbolic equation.Figure 5 shows the projection of the space of the large system into that of auxiliary system and the lifting procedure [43,44].Figure 5: Projection of the space of the large system into that of auxiliary system and the lifting procedure using the solution of auxiliary system.
The boundary conditions for slug problem equation ( 7) are reformulated for coordinates (, ) as Figure 2 shows how the initial and boundary conditions for the large system (4) and ( 6) are mapped into those for auxiliary system and the lifting equations ( 28)- (30).
The initial conditions for slug problem equation ( 8) are reformulated for coordinates (, ) as The solution of the Riemann problem for  < 1 corresponds to the following initial and boundary conditions:

Solution for the Riemann Problem
Let us discuss the solution of the problem equations ( 7) and (8) for  < 1, which is self-similar; that is, the boundary and initial conditions become ( 9) and (10).
The mass balance conditions on shocks which follow from the conservation law (Hugoniot-Rankine condition) form of the system (28)- (30) are where  is reciprocal to the shock velocity of ( 28)- (30).
As salt and polymer concentration are connected by the thermodynamic equilibrium relationship (, ), function  is discontinuous if  is discontinuous, so is .Since  is a function of  and , discontinuity of  and  yields discontinuity of .

Solution for the Lifting Equation
Points 2 and 3 are determined by the condition of equality of  and  shock speeds: Point 4 is determined by condition of equality of the shock velocities , , and : Point 4 is connected to point  by -shock:  The solution of the Riemann problem equations ( 28)- (30) with free variables (, ) is given by the following formulae: The expression  = −(, ) allows calculating saturation (, ): Figure 8 shows the solution of the system ( 28)-( 30) in (, )plane.For  < 1, the solution is self-similar; the wave interaction occurs at  > 1.
Figure 8: Solution of the auxiliary and lifting system for slug problem in (, )-plane.

Inverse Mapping: Change of Variables from (𝜑, 𝑥) to (𝑡, 𝑥).
Time  = (, ) for solution is calculated from (12) along any path from point (, ) to point (0, 0).The expression for time t in zone II is The expression for time  in zone III is Figure 9: Non-self-similar solution of the problem for wave interactions in (, )-plane.
In zone IV, integral for calculation time,  = ∫ , 0,0 (/ + /) is calculated along the characteristic in centred wave: Figure 9 shows the solution for the Riemann problem at  < 1; see Figure 10 for detailed description of the Riemann solution and profiles of unknown functions.Finally, the solution of the Riemann problem for the system (4)-( 6) is

Solution of the Slug Problem
Now let us solve the slug problem equations ( 31) and (32) for auxiliary system ( 29) and (30).The solution of Riemann problem at the point (0, 1) is given by -shock with  − =  2 and  + =  1 under constant : The solution of the auxiliary system is given by (49).So, zone I in Figure 8 corresponds to initial conditions, the solution is given by point 4 in zone II, and point 3 holds  in zone III.Centred waves equation (39) fills in zone IV.In zone V,  =  2 and  = 0. Now let us solve the lifting equation (28) with given (, ) and (, ).
The Hugoniot-Rankine condition for the rear slug front is is constant along the characteristic lines behind the rear front where point (  ,   ) is located on the rear front and is located on the same characteristic line with point (, ): The solution of lifting equation (, ) is given by different formulae in zones I-V: where the equation for rear front of the chemical slug in the auxiliary plane is From ( 54), (55) follow the expression for   (  ) in a parametric form: ( + ) =  ( + ,  1 ) −  1 +   ( + ,  1 )  ( + ,  1 ) −   ( + ,  1 ) ( + + Γ) =  ( + ,  1 ) −  1 +   ( + ,  1 ) Δ . (57) Integration of the form (41) along the rear front gives Finally, the solution of the slug problem for the system (4)-( 6) subject to initial and boundary conditions equations ( 7) and ( 8) is (Figure 9) (59) Figure 11 presents trajectories of shock fronts in (, )-plane along with profiles of unknowns , , and  at typical moments.
Here the trajectory of the rear slug front   =   () is given in a parametric form (Figure 12) where   is the abscissa of point  and   is the ordinate of point  (Figure 12).Equations (60) can be solved graphically.Straight line AB is a tangent to the fractional flow curve  =  1 , the tangent point in +.The rear front position   is determined by the interval BC at the moment determined by AC.

Particular Case for the Polymer Absence in the Reservoir before the Injection
In reality, there is no chemical initially in the reservoir during the majority of chemical enhanced oil recovery applications; that is, (, 0) = 0.For zero initial polymer concentration, the intermediate polymer concentration is equal to zero, so the points ahead and behind the -shock coincide in planes (, ) and (, ).The particular simplified solution is (Figures 13 and 14)

Fluid Mechanics Interpretation of the Solution
Following exact solution equations ( 4)-( 9), let us describe structure of two-phase flow with chemical and salt additives during chemical slug injection.
During continuous injection  < 1, the solution of chemical slug injection coincides with that of continuous chemical injection.Initial conditions equation ( 10) is shown by point  that corresponds to low initial saturation and initial concentrations of chemical  2 and of salt  = 1.
The trajectories of shocks 2->3, 3->4, and 4->  are shown in Figure 7. Shock velocities are constant, so the trajectories are straight lines.Let us fix the position  = 1 of the raw of production wells.Before arrival of the front 4->  at the moment  = 1/ 3 , oil with fraction of water   and initial concentrations of chemical and salt is produced.After arrival of the front, water-oil mixture with water fraction  4 and initial concentrations of chemical and salt is produced until the arrival of the 3->4 front at the moment  = 1/ 2 .
The corresponding profiles of saturation and concentrations are shown in Figure 10.The moment  1 for profiles is fixed in Figure 10

Figure 2 :
Figure 2: Introduction of potential function (Lagrangian coordinate) and mapping between independent variables.

Figure 3 :
Figure 3: Derivation of mass balance equation in Eulerian and Lagrangian coordinate systems.

Figure 4 :
Figure 4: Speeds of a particle in Eulerian and Lagrangian coordinates.

4. 1 .
Solution for the Auxiliary System.The solution of auxiliary system is presented in Figure 6 by sequence of c-shock from point  into intermediate point and (, )-shock into point .The corresponding formulae are as follows:

Figure 7
exhibits initial and boundary conditions for hydrodynamics lifting equation(28).Curves  = (, ) are shown for constants  =  1 ,  =  2 , and  =  * ; they are obtained from fractional flow curves  = (, ) for the same constant values of concentration .Point  corresponds to  tending to infinity and  tending to minus infinity, where the fractional flow  tends to zero.The tangent of the segment (0, 0)-(, ) tends to −  .

Figure 6 :
Figure 6: Solution of the auxiliary problem.(a) Adsorption isotherm for chemical for different water salinities and the Riemann problem solution; (b) Riemann problem solution on the plane of chemical concentration c and salinity .

Figure 7 :
Figure 7: The image of the solution in (-) plane.

Figure 11 :
Figure 11: The solution of the slug injection problem: (a) trajectories of shock fronts and characteristic lines in (, )-plane; (b) saturation profile; (c) chemical concentration profile; (d) salinity profile.
(a), allowing defining positions of all fronts in this moment.Corresponding profiles at that moment for saturation, chemical concentration, and salinity are shown in Figures 10(b), 10(c), and 10(d), respectively.The saturation profile consists of declining interval   - 2 in -wave, two oilwater banks  3 and  4 , and the initial undisturbed zone   .The chemical concentration profile consists of injected value  =  1 in -wave, intermediate value  * in  3 -bank, and initial