Inverse Problem of Parameter Identification for Extended Compositional Gradient Model

In this study we address simulation of initial distribution of fluid composition in hydrocarbon reservoirs, particularly in scenarios where slow vertical segregation in a low-permeability near-critical reservoir results in an incomplete relaxation to the quasi-equilibrium state. To effectively analyze these complex scenarios, we propose an extended compositional gradient model that takes into account incomplete system relaxation. Practical application of the method involves solution of the forward problem of compositional profiling, with the model parameters adjusted according to real fluid samples. We formulate the inverse problem of parameter identification for the extended compositional gradient model and solve it using gradient-free optimization methods. The adjusted parameters are specific enthalpies of fluid components and the relaxation parameter. The overarching goal of parameter identification is to closely match the predictions of the model with real data, thereby enhancing the reliability of forecasts regarding distribution of fluid composition. We believe that the results of this study will significantly influence both understanding and simulation techniques for distribution of hydrocarbon components in natural reservoirs, particularly in situations characterized by incomplete system relaxation.


Introduction
In any detailed study of hydrocarbon reservoirs, it is essential first to assess the spatial distribution of the initial fluid composition, both horizontally and vertically.This distribution often reflects the quasiequilibrium state of the multicomponent hydrocarbon system under the influence of gravity and geothermal fields.This phenomenon, referred to as the compositional gradient, usually leads to an increase in fluid density and concentration of heavier components with depth.
However, for near-critical hydrocarbon deposits located in deep low-permeability reservoirs the vertical fluid segregation could be so slow that the system hasn't yet reached the quasi-equilibrium state.For such scenarios, an extended compositional gradient is necessary.
In this paper we discuss a method for solving the forward problem of compositional gradient in cases of incomplete system relaxation to the quasi-equilibrium state.The computational algorithm is similar to the classical compositional gradient calculations with an extra parameter to control relaxation.Practical application of the method involves adjustment to the real fluid samples data by addressing the inverse problem of parameter identification solved using gradient-free optimization methods.Parameters of interest are the specific enthalpies of fluid components and the relaxation parameter, which are uncertain but control the initial distribution of pressure and fluid composition in the extended model.Therefore, the inverse problem provides the parameters for the forward problem, enabling more reliable calculation of the initial fluid composition versus depth.The objective is to ensure that the model predictions match real data as closely as possible, which is illustrated with actual fluid samples from a gas-condensate reservoir.Different forms of the misfit function and minimization algorithms are examined for the parameter identification process.The choice of weights in the misfit function is shown to be an effective tool for tuning the inverse problem solution to meet engineering requirements for a specific hydrocarbon reservoir.
We believe that this study will make a significant contribution to understanding and modeling the hydrocarbon reservoirs, particularly those where the initial fluid system has not reached the quasiequilibrium state.We aim to improve the comprehension of the compositional gradient in hydrocarbon reservoirs and enhance the accuracy of predicting the initial distribution of fluid components with depth.

Methods
The problem of compositional gradient calculation (the compositional profiling) under incomplete system relaxation is relevant to deep near-critical hydrocarbon reservoirs in low-permeability formations.In such conditions, gravitational segregation of components may not have been fully achieved at the start of the reservoir development.The method used in this study is based on the concept of system relaxation from the initially homogeneous state towards the quasi-equilibrium state in the gravity and geothermal fields.[1] It is assumed that the reservoir was originally filled with a fluid of homogeneous composition, after which the system relaxed over a period of time towards the quasi-equilibrium state.This process is modeled by the relaxation-type multiplier in the right-hand side of the conventional system of equations for the difference of chemical potentials of individual fluid components between two vertical positions ℎ 1 and ℎ 2 [2].The relaxation parameter  characterizes the ratio of the time elapsed since the start of reservoir formation  to the characteristic relaxation time  [1]: =   ,  = 1, … ,  , Here  is the pressure,  ̅ is the vector of molar component concentrations (molar composition) of the fluid;  1 and  2 , ℎ 1 and ℎ 2 are the temperature and depth values for the first and second vertical positions considered;   ,   and   are the fugacity, the molecular weight and the absolute partial molar enthalpy of component ;  and  are the molecular weight and the absolute molar enthalpy of the mixture (fluid) [1,2].
When  is zero, the right-hand side of equation ( 1) vanishes, reflecting the situation of equal fugacities at any two different vertical positions.This corresponds to the constancy of pressure and fluid composition with depth, i.e., to the initial homogeneous state.At large values of  the system represented by equation (1) converges to the classical model of compositional gradient with gravity and thermodiffusion [2].That is, as time increases, a relaxation is assumed towards the classical quasiequilibrium compositional gradient.
The component fugacities for given pressure and mixture composition are determined by solving a cubic equation of state (EoS), for example, the Peng-Robinson EoS widely accepted in phase behavior simulations for petroleum fluids.The system of equations (1) at each depth step is solved with respect to the pressure and fluid composition using the Newton's method.Additional procedures are used to locate and treat the gas-liquid (gas-oil) contact.The algorithm is illustrated below (Figure 1).More details can be found in [1].

Figure 1. Algorithm of the compositional gradient calculations
The inverse problem of parameter identification is required to calibrate the model to actual experimental data from fluid samples.In this study, the identified parameters of the compositional gradient model were:   -absolute enthalpy of the components in the ideal gas state at T=273.15K, and  -the relaxation parameter controlling incomplete quasi-equilibrium segregation.These parameters were selected due to their high level of uncertainty and strong influence on the simulation results.
The misfit criterion to be minimized is defined as a weighted sum of squares of the normalized residuals for measured data, including the molar fractions (concentrations) of fluid components, the pressure, and the  5+ : where  ′ () is the calculated concentration of the i-th component, () is the corresponding experimental (determined from fluid sample) concentration, and the -th component is excluded due to normalization of fluid composition;  ′ and  are the calculated and experimental pressure values,  5+ ′ and  5+ are the calculated and experimental values of the potential content of  5+ hydrocarbons, and   ,   and   5+ are the weights of individual terms in the overall misfit criterion.Potential content of  5+ is a widely adopted indicator of 'richness' of the gas-condensate fluid in liquid-type hydrocarbons.It is defined as the mass content of  5+ components per m 3 of gas at standard conditions.
Because of the complex model structure and implicit dependencies between the measured data and the model parameters, gradient-free optimization methods are used to solve the inverse problem.
It is typical for real-world situations that the number of measured data from samples is insufficient, and values included in the misfit function can be correlated (due to normalization of fluid composition and specifics of measurements).This leads to the typical ill-posedness of inverse problems, manifested in the non-uniqueness and non-robustness of the obtained solution.This is particularly important in the context of optimization, where non-uniqueness can lead to various and possibly conflicting solutions minimizing the overall misfit criterion.Therefore, additional constraints or strategies may need to be applied to select the optimal set of parameters for model matching.
Specifically, for the problem at hand, we will in certain cases limit the number of parameters subject to identification to ensure uniqueness and the robustness of the solution.The list of components for which enthalpy parameters are adjusted is determined taking into account the specifics of the dataset and depending on the available additional constraints or a priori information about the system.This helps to control the degrees of freedom in the optimization problem and reduces the risk of convergence to unphysical solutions.

Results and Discussion
To study the inverse problem solution in realistic situations, we use a dataset based on three fluid samples of a real oil-gas-condensate reservoir.The simulated system is similar to the conditions of the Achimov deposits of the Urengoy oil-gas-condensate field.The fluid is represented by a mixture of 24 components characterized in Table 1:  1 ,  2 and  3 are the compositions of the three samples,  is the component molecular weight,   and   are the critical temperature and pressure, ω is the acentric factor,   , , ,  and  are the parameters used for enthalpy calculations (see [1] for details).The reference point for the forward problem calculations is set at a depth of 3754.6 m, with the composition  1 , pressure 633 bars, and temperature 108.9°C.Geothermal (temperature) gradient is 0.029°C/m.The compositional profiling is performed for the depth range of 3614 to 4050 m, with the step of 1 m.For the inverse problem, the compositions  2 and  3 are known at the depths ℎ 2 = 3829.12m and ℎ 3 = 4001.25 m, with the values of  5+ 2 = 932.69g/cm 3 ,  5+ 3 = 1269.14g/cm³.For the point ℎ 3 , the pressure of 650 bars is also known.
With these data, a series of inverse problem calculations were carried out.During the study, computational experiments with multiple optimization parameters were conducted to find a set of possible solutions.Each calculation was performed with unique parameter settings to explore various scenarios and identify the most effective configurations to enhance the optimization results.In particular, we examined various sets of weights in the misfit function to investigate how weights can impact the inverse problem solution.
The optimization process was carried out using different optimizers such as Bayes Optimizer [3], Particle Swarm Optimizer [4], Differential Evolution [5], CMA-ES [6] from the Nevergrad library [7], and BOBYQA [8] with Py-BOBYQA implementation [9].A limit of 100 iterations was set as the stopping criterion.As experiments show, the optimizers often find different solutions, which underline incorrectness of the inverse problem and possible lack of efficiency of some gradient-free methods for the problems considered.
The enthalpy values were adjusted through multipliers:   =  0  , where  0  is the original value and  is the multiplier to be identified.Thus, the set of optimization variables is: Figure 2 presents the graphs of  5+ and pressure versus depth for the set of solutions to the inverse problem of parameter identification found by different optimization algorithms.The plots of the total misfit function and its individual terms are shown in Figure 3 for the best current optimization solution versus the iteration number.The parameters of the best solutions are given in Table 2.The values of the misfit function and its individual terms for the best solutions are provided in Table 3.
From Figure 2 it can be seen that the shape of the calculated  5+ curves does not match the experimental data from fluid samples.Therefore, it is not enough to control the misfit for composition, and misfit for  5+ should be explicitly included with nonzero weight into the optimization function.
Figures 2 and 3 show that the behavior of the optimization algorithms and the solutions found are quite different.This indicates the presence of multiple local minima in the inverse problem with close values of the total misfit, but different values of the misfit for pressure, composition, and  5+ .

Experiment 2
This experiment is similar to Experiment 1, but with the normalized residual for  5+ introduced into the misfit function: The results obtained for Experiment 2 are given in Figures 4 and 5, Tables 4 and 5.As seen from Figure 4, the identified models from Experiment 2 better match the experimental data of the fluid samples.However, since the total misfit value is being optimized, it's important to pay attention to the order of magnitude of individual terms and adjust the weights to equalize the relative contribution of different measurements.

Experiment 3
Experiment 3 is similar to Experiment 2, but with adjusted weight on pressure misfit, so that it was same order of magnitude as the misfit for composition:   = 1,   = 10000,   5+ = 1.The results are presented in Figures 6 and 7, Tables 6 and 7.

Experiment 4
In this experiment we use the same settings as in Experiment 3, but increase the weight on  5+ :   = 1,   = 10000,   5+ = 1000.Figures 8 and 9, and Tables 8 and 9 summarize the results.It's worth noting again that the optimization methods converge to significantly different points, while they are roughly similar in terms of the final value of the misfit function.Moreover, even similar solutions in terms of  5+ and pressure can be quite different when composition for individual components is considered (see blue and violet plots in Figures 8 and 9 as an example).

Discussion
The presented computational experiments demonstrate some negative features of the inverse problem: non-unique solutions corresponding to local minima of the total misfit and high sensitivity to the choice of weights in the misfit function.However, by adjusting the weights in the misfit function one can control the model's behavior according to the most certain data and properties of interest.The obtained solutions reflect the range of possible tuning: more accurate matching of the data on  5+ and composition, or on pressure.It's also recommended to include more points with data in the inverse problem than in the experiments considered.The absence of solutions that accurately match the fluid samples data on  5+ and pressure at all depths indicates the need for additional verification of the experimental data, as well as their compliance with the assumptions embedded in the composition gradient models.

Conclusions
In this paper we discussed the extended compositional gradient model for simulation of initial distribution of hydrocarbon composition in reservoirs with incomplete vertical fluid segregation.The forward and inverse problems were formulated, and their solution was illustrated by computational experiments on a realistic dataset based on fluid samples of an oil-gas-condensate reservoir.
Though the study underlined the intrinsic negative features of the inverse problem, it proved that the choice of weights in the misfit function can help to effectively tune the model according to the practical requirements of a certain hydrocarbon reservoir.

Figure 2 .Figure 3 .
Figure 2. The results of calculations for models obtained by different optimization algorithms (Experiment 1)

Figure 4 .Figure 5 .
Figure 4.The results of calculations for models obtained by different optimization algorithms (Experiment 2)

Figure 6 .Figure 7 .
Figure 6.The results of calculations for models obtained by different optimization algorithms (Experiment 3)

Figure 8 .Figure 9 .
Figure 8.The results of calculations for models obtained by different optimization algorithms (Experiment 4)

Table 1 .
Input data for fluid compositions and component parameters

Table 2 .
Parameters of the best solutions for Experiment 1

Table 3 .
Misfit function values at the best solutions for Experiment 1

Table 4 .
Parameters of the best solutions for Experiment 2

Table 5 .
Misfit function values at the best solutions for Experiment 2

Table 6 .
Parameters of the best solutions for Experiment 3

Table 7 .
Misfit function values at the best solutions for Experiment 3

Table 8 .
Parameters of the best solutions for Experiment 4

Table 9 .
Misfit function values at the best solutions for Experiment 4