Next Article in Journal
Thermochemical Analysis of a Packed-Bed Reactor Using Finite Elements with FlexPDE and COMSOL Multiphysics
Next Article in Special Issue
Experimental and Numerical Studies of Fine Quartz Single-Particle Sedimentation Based on Particle Morphology
Previous Article in Journal
Coal Mine Personnel Safety Monitoring Technology Based on Uncooled Infrared Focal Plane Technology
Previous Article in Special Issue
Study of the Fluid Passing through the Screen in the Three Products Hydrocyclone Screen (TPHS): A Theoretical Analysis and Numerical Simulation
 
 
Font Type:
Arial Georgia Verdana
Font Size:
Aa Aa Aa
Line Spacing:
Column Width:
Background:
Article

A Homotopy Method for the Constrained Inverse Problem in the Multiphase Porous Media Flow

1
School of Mathematics and Statistics, Northeastern University at Qinhuangdao, Qinhuangdao 066004, China
2
School of Electrical and Electronic Engineering, Nanyang Technological University, Singapore 639798, Singapore
3
Eighth Geological Brigade, Hebei Bureau of Geology and Mineral Resources Exploration, Qinhuangdao 066000, China
*
Author to whom correspondence should be addressed.
Processes 2022, 10(6), 1143; https://doi.org/10.3390/pr10061143
Submission received: 10 May 2022 / Revised: 1 June 2022 / Accepted: 6 June 2022 / Published: 7 June 2022
(This article belongs to the Special Issue Computational Modeling of Multiphase Flow (II))

Abstract

:
This paper considers the constrained inverse problem based on the nonlinear convection-diffusion equation in the multiphase porous media flow. To solve this problem, a widely convergent homotopy method is introduced and proposed. To evaluate the performance of the mentioned method, two numerical examples are presented. This method turns out to have wide convergence region and strong anti-noise ability.

1. Introduction

The multiphase flow in porous media is of great importance in biology [1], chemistry [2], civil engineering [3], mechanics [4], and geosciences [5]. The porous media flow in many engineering problems, including reservoir engineering and reservoir simulation engineering, is often used to simulate the gas and oil flow within the reservoir. The authors in [6] proposed the fractional flow formulation of the multiphase flow equations in porous media, which leads to a nonlinear convection-diffusion saturation equation. Recently, Saeed et al. exactly analyzed the second grade fluid [7] and viscoelastic liquid with single slip assumption [8], and studied the natural convection flow [9] from in mathematical perspective. Abdeljawad et al. [10] and Firdous et al. [11], respectively, considered the MHD Maxwell fluid and Powell-Eyring fluid. Riaz et al. [12,13] mainly researched on the MHD Oldroyd-B fluid. Nowadays there are many different studies that aim to solve the inverse problem in multiphase porous media flow. Hazra et al. [14] solved the direct and inverse modeling in multiphase porous media flow using numerical simulation techniques, and Wang and Zabaras [15] identified the contamination source in multiphase porous media flow. Nilssen et al. [16] recovered the diffusion parameters in multiphase porous media flow with the augmented Lagrangian method. With the development of artificial intelligence, several recent studies focused on the application of deep learning models to the inverse problem in multiphase porous media flow [17,18,19,20].
This paper considers the identification of space-dependent permeability k for the nonlinear convection-diffusion equation in the multiphase porous media flow:
u t + · ( ϕ , φ ) · ( k · N ( u ) u ) = f ( x , t ) , ( x , t ) Ω × ( 0 , T ) ,
under the boundary and initial conditions
u ( x , t ) = χ ( x , t ) , ( x , t ) Ω × ( 0 , T ) ,
u ( x , 0 ) = ψ ( x ) , x Ω .
An additional condition is
u ( x s , t ) = γ ( x s , t ) , s = 1 , 2 , , S , t ( 0 , T ) ,
where ϕ and φ are both S-shaped Buckley–Leverett flux functions, N is a nonlinear diffusion function, f is a source function that is piecewise smooth, Ω is set to be an unit square for simplicity, and γ is the measured data (seepage velocity data).
Equation (1) correlates with the multiphase porous media flow. The partial differential equation (PDE) system describing the immiscible displacement of oil by water in porous media under no gravity condition is as follows [6]:
· V = f 1 ( x , t ) , V = k · ς ( x , u ) ( p ϑ ( u ) h ) , β ( x ) u t + · ( ϕ ( u ) V + k · φ ( u ) h ) · ( k · N ( u ) u ) = f 2 ( x , t ) ,
where ϕ = ς r a ς , ς r a = k r a μ , φ = ( ϑ ρ ) ϕ ϵ , and the definitions of model parameters are listed in Table 1. Units for k , V , p, h, ρ are, respectively, m 2 , m / s , N, m, k g / m 3 , and other parameters are dimensionless.
Equation (1) closely resembles Equation (5) if the time derivative and convection terms are not considered. The coefficients of time derivative terms in Equations (1) and (5) are respectively constant 1 and β ( x ) ; the convection term in Equation (5) has varying coefficient and permeability dependence, and the one in Equation (1) does not.
Oil reservoir simulation on the basis of this inverse problem is an effective tool, which can provide help for petroleum reservoir management, such as the choices of the fluid production and injection rates, well locations, imaging method (see Figure 1). Generally, the permeability model has problems with equivalence, non-uniqueness, hidden or suppressed layers and lack of model resolution, because the measured data are insufficient, inconsistent and inaccurate [21]. Moreover, the primary difficulties of the traditional methods (e.g., Levenberg–Marquardt, Gauss–Newton) are the local convergence property and that the cost function has numerous local minima. To cope with these problems, a homotopy method is developed associated with the constraints of well logs to identify the permeability coefficient.
The homotopy method has global convergence under certain weak assumptions [22], and so has been successfully used to solve various nonlinear problems [23,24,25,26,27,28]. Recently, some authors also used it to solve inverse problems. For instance, Mallick et al. [29,30] used this method to estimate the thermal parameters within annular fins. Biswal et al. [31] applied the homotopy method for the inverse analysis of Jeffery–Hamel flow problem. Sattari Shajari and Shidfar [32] studied the homotopy solution of the wave equation inverse source problem. Liu [33,34] combined the homotopy method with multiscale ideas to develop the hybrid algorithm, and applied it to nonlinear inverse problems. Hu et al. [35] identified the parameters of a cracked beam using the homotopy algorithm. Courbot and Colicchio [36] analyzed the solution of the gridless sparse recovery problem using the homotopy method. The homotopy method has also been applied to the fractional inverse Stefan problem [37], the backward heat conduction problem [38], the porosity reconstruction on the basis of Biot elastic model [39], etc.
Usually, the measured data have a low signal-to-noise ratio. With the aim of restraining noise and improving identified model quality, the constraint condition has a wide application in many inversion fields [40,41,42,43]. The reason lies in that the constraint data are recorded from the interior of the model to be identified, and are less noisy than the measured data.
In previous works [44,45], we have verified the effectiveness of the homotopy and multigrid-homotopy methods for the inverse problem of the nonlinear diffusion equation:
u t · ( k · N ( u , u ) u ) = f ( x , t ) ,
which is an intermediate step of permeability identification in multiphase porous media flow. Different from [44,45], this paper not only considers the permeability identification based on the nonlinear convection-diffusion Equation (1), which can more accurately describe the multiphase flow process in porous media than the nonlinear diffusion Equation (6), but also introduces a well-log constraint to this inverse problem. The resulted constrained inverse problem can be transformed into a nonlinear constrained optimization problem. Because of the ill-posedness of problem and the local convergence property of traditional methods, we first impose Tikhonov regularization, and then use a widely convergent homotopy method to solve the normal equation of the regularized cost function. The final part of the paper consists of a report of numerical experiments that demonstrates the performance of the method.

2. Discretization

Equations (1)–(4) can be discretized by using the finite-difference scheme as follows:
u i , j k u i , j k 1 h t + · ( ϕ ( u i , j k 1 ) , φ ( u i , j k 1 ) ) · ( k i , j N i , j k u i , j k ) = f i , j k , i = 1 , , I 1 ; j = 1 , , J 1 ; k = 1 , , K , u 0 , j k = χ 0 , j k , j = 0 , , J ; k = 1 , , K , u 1 , j k = χ 1 , j k , j = 0 , , J ; k = 1 , , K , u i , 0 k = χ i , 0 k , i = 0 , , I ; k = 1 , , K , u i , 1 k = χ i , 1 k , i = 0 , , I ; k = 1 , , K , u i , j 0 = ψ i , j , i = 0 , , I ; j = 0 , , J , u x s k = γ x s k , s = 1 , , S ; k = 1 , , K ,
where
u i , j k = u ( i h x , j h y , k h t ) , k i , j = k ( i h x , j h y ) , f i , j k = f ( i h x , j h y , k h t ) , ψ i , j = ψ ( i h x , j h y ) , χ i , j k = χ ( i h x , j h y , k h t ) , γ x s k = γ ( x s , k h t ) , I = 1 / h x , J = 1 / h y , K = T / h t ,
h t , h x , h y are, respectively, the time and spatial step lengths. · ( ϕ ( u i , j k 1 ) , φ ( u i , j k 1 ) ) and · ( k i , j N i , j k u i , j k ) are, respectively, the discretizations of convection and diffusion terms, which will be described in the Appendix A and Appendix B.
By Equation (7), this inverse problem can be formulated as the following nonlinear operator equation:
R ( K ) = Γ ,
where
K = ( k 1 , 1 , k 1 , 2 , , k 1 , J , k 2 , 1 , k 2 , 2 , , k 2 , J , , k I , 1 , k I , 2 , , k I , J ) , Γ = ( γ x 1 1 , γ x 2 1 , , γ x S 1 , γ x 1 2 , γ x 2 2 , , γ x S 2 , , γ x 1 K , γ x 2 K , , γ x S K ) .
The measured data are denoted by γ ¯ x s k , which can form the vector Γ ¯ in the same sequence as Γ :
Γ ¯ = ( γ ¯ x 1 1 , γ ¯ x 2 1 , , γ ¯ x S 1 , γ ¯ x 1 2 , γ ¯ x 2 2 , , γ ¯ x S 2 , , γ ¯ x 1 K , γ ¯ x 2 K , , γ ¯ x S K ) ,
and the permeability, known from the well logs of a well located at point i 0 in the x -direction, is denoted by
K ¯ i 0 = ( k ¯ i 0 , 1 , k ¯ i 0 , 2 , , k ¯ i 0 , J ) .
Let Y denote a matrix
Y = 0 0 0 1 0 0 0 0 0 0 0 0 0 1 0 0 0 0 0 0 0 0 0 1 0 0 0 J × ( I × J ) ,
such that Y K = ( k i 0 , 1 , k i 0 , 2 , , k i 0 , J ) , then Equation (8) turns into a nonlinear constrained optimization problem
min K R ( K ) Γ ¯ 2 , subject   to   Y K = K ¯ i 0 .
By combining the objective function and constraint into a penalty function, we can attack Equation (9) by solving an unconstrained problem. For instance, a penalty function can be defined as
R ( K ) Γ ¯ 2 + ω Y K K ¯ i 0 2 ,
where ω is the penalty parameter. After that, the solution of Equation (9) can be obtained by minimizing this unconstrained function Equation (10) for a large value of ω .

3. Identification Method

In the first part of this section, a basic iterative method is given by the successive linearization technique. Then, in the second part, the strategy used for the development of a homotopy method is presented.

3.1. Basic Iterative Method

It is common knowledge that the inverse problem is ill-posed, therefore the regularization technique has to be used for the improvement of numerical stability of the algorithm:
min K { R ( K ) Γ ¯ 2 + ω Y K K ¯ i 0 2 + ϖ 1 W 1 K 2 + ϖ 2 W 2 K 2 } ,
where ϖ 1 , ϖ 2 denote the regularization parameters, W 1 , W 2 denote, respectively, the second-order smooth matrices in the x - and y -direction.
W 1 K = 0 k 1 , 1 2 k 2 , 1 + k 3 , 1 k 2 , 1 2 k 3 , 1 + k 4 , 1 k I 2 , 1 2 k I 1 , 1 + k I , 1 0 0 k 1 , 2 2 k 2 , 2 + k 3 , 2 k 2 , 2 2 k 3 , 2 + k 4 , 2 k I 2 , 2 2 k I 1 , 2 + k I , 2 0 0 k 1 , J 2 k 2 , J + k 3 , J k 2 , J 2 k 3 , J + k 4 , J k I 2 , J 2 k I 1 , J + k I , J 0 , W 2 K = 0 k 1 , 1 2 k 1 , 2 + k 1 , 3 k 1 , 2 2 k 1 , 3 + k 1 , 4 k 1 , J 2 2 k 1 , J 1 + k 1 , J 0 0 k 2 , 1 2 k 2 , 2 + k 2 , 3 k 2 , 2 2 k 2 , 3 + k 2 , 4 k 2 , J 2 2 k 2 , J 1 + k 2 , J 0 0 k I , 1 2 k I , 2 + k I , 3 k I , 2 2 k I , 3 + k I , 4 k I , J 2 2 k I , J 1 + k I , J 0 .
It is easy to show that Equation (11) is equivalent to (its normal equation):
R ( K ) ( R ( K ) Γ ¯ ) + ω Y ( Y K K ¯ i 0 ) + ( ϖ 1 W 1 W 1 + ϖ 2 W 2 W 2 ) K = 0 .
For the sake of avoiding the influence of second derivative, a successive linearization method can be introduced to construct a basic iterative method.
We first assume that the kth approximation K k of Equation (12) has been obtained, and use the linear function
E k ( K ) = R ( K k ) ( K K k ) + R ( K k ) ,
instead of R ( K ) in Equation (12):
R ( K k ) ( R ( K k ) ( K K k ) + R ( K k ) Γ ¯ ) + ω Y ( Y K K ¯ i 0 ) + ( ϖ 1 W 1 W 1 + ϖ 2 W 2 W 2 ) K = 0 .
Then, by denoting the solution of Equation (13) as K k + 1 , we obtain the basic iterative scheme:
K k + 1 = K k + Δ K k , k = 0 , 1 , 2 , [ R ( K k ) R ( K k ) + ω Y Y + ϖ 1 W 1 W 1 + ϖ 2 W 2 W 2 ] Δ K k = [ R ( K k ) ( R ( K k ) Γ ¯ ) + ω Y ( Y K k K ¯ i 0 ) + ( ϖ 1 W 1 W 1 + ϖ 2 W 2 W 2 ) K k ] .
This method is fast and stable, but only has local convergence.

3.2. Homotopy Method

In order to expand the domain of convergence, the homotopy method is introduced to solve Equation (12) by taking into account the fixed-point homotopy equation
A ( K , a ) = a [ R ( K ) ( R ( K ) Γ ¯ ) + ω Y ( Y K K ¯ i 0 ) + ( ϖ 1 W 1 W 1 + ϖ 2 W 2 W 2 ) K ] + ( 1 a ) ( K K 0 ) = 0 ,
where a [ 0 , 1 ] and K 0 are, respectively, the homotopy parameter and arbitrary initial value. From Equation (15), it is easy to see
A ( K , 0 ) = K K 0 , A ( K , 1 ) = R ( K ) ( R ( K ) Γ ¯ ) + ω Y ( Y K K ¯ i 0 ) + ( ϖ 1 W 1 W 1 + ϖ 2 W 2 W 2 ) K .
As the homotopy parameter a changed continuously from 0 to 1, the trivial problem A ( K , 0 ) = 0 is continuously deformed to the original problem A ( K , 1 ) = 0 , that is, K is changed from K 0 to the solution of Equation (12). In topology, this is called deformation.
To achieve the specific numerical algorithm, we first divide the interval [ 0 , 1 ] into 0 = a 0 < a 1 < < a B = 1 , and then solve A ( K , a b ) = 0 ( b = 1 , 2 , , B ) sequentially by some iterative scheme, whose initial value is chosen as the solution K b 1 of the previous equation A ( K , a b 1 ) = 0 .
For A ( K , a b ) = 0 , in the same manner as the construction of Equation (14), we have
K b k + 1 = K b k + Δ K b k , k = 0 , 1 , , b m , [ a b R ( K b k ) R ( K b k ) + a b ω Y Y + a b ϖ 1 W 1 W 1 + a b ϖ 2 W 2 W 2 + ( 1 a b ) I ] Δ K b k = [ a b R ( K b k ) ( R ( K b k ) Γ ¯ ) + a b ω Y ( Y K b k K ¯ i 0 ) + ( a b ϖ 1 W 1 W 1 + a b ϖ 2 W 2 W 2 ) K b k + ( 1 a b ) ( K b k K 0 ) ] , K b 0 = K b 1 , K b = K b b m + 1 , b = 1 , 2 , , B ,
where I refers to the unit matrix. After K B is obtained, Equation (14) may be used to make correction. Therefore, Equations (14) and (17) are combined into a stabilized method which has a wider domain of convergence for the constrained permeability identification of the nonlinear convection-diffusion equation.
By choosing a b = b B and b m = 0 , Equation (17) can be simplified as
K b + 1 = K b + Δ K b , b = 0 , 1 , , B 1 , [ b B R ( K b ) R ( K b ) + b B ω Y Y + b B ϖ 1 W 1 W 1 + b B ϖ 2 W 2 W 2 + ( 1 b B ) I ] Δ K b = [ b B R ( K b ) ( R ( K b ) Γ ¯ ) + b B ω Y ( Y K b K ¯ i 0 ) + ( b B ϖ 1 W 1 W 1 + b B ϖ 2 W 2 W 2 ) K b + ( 1 b B ) ( K b K 0 ) ] .
To test the necessity of introduction of constraints, we also give the homotopy method for the ordinary permeability identification for the nonlinear convection-diffusion equation in the multiphase porous media flow, and compare Equations (14) and (18) with it. Specifically, let ω = 0 , then Equations (14) and (18) turn into the homotopy method for the ordinary permeability identification problem
K b + 1 = K b + Δ K b , b = 0 , 1 , , B 1 , [ b B R ( K b ) R ( K b ) + b B ϖ 1 W 1 W 1 + b B ϖ 2 W 2 W 2 + ( 1 b B ) I ] Δ K b = [ b B R ( K b ) ( R ( K b ) Γ ¯ ) + ( b B ϖ 1 W 1 W 1 + b B ϖ 2 W 2 W 2 ) K b + ( 1 b B ) ( K b K 0 ) ] ,
and
K k + 1 = K k + Δ K k , k = 0 , 1 , 2 , [ R ( K k ) R ( K k ) + ϖ 1 W 1 W 1 + ϖ 2 W 2 W 2 ] Δ K k = [ R ( K k ) ( R ( K k ) Γ ¯ ) + ( ϖ 1 W 1 W 1 + ϖ 2 W 2 W 2 ) K k ] .

4. Numerical Experiments

This section tested two synthetic examples to show the performance of our proposed method. The parameters required in the identification process are given by
ϕ ( u ) = u 2 ( 1 5 ( 1 u 2 ) ) u 2 + ( 1 u ) 2 , φ ( u ) = u 2 u 2 + ( 1 u ) 2 , N ( u ) = u 2 u + 1 , ψ ( x ) = s i n ( π x ) s i n ( π y ) , χ ( x , t ) = 0 , ω = 10 4 , ϖ 1 = ϖ 2 = 10 5 , K 0 5 , h x = h y = 1 24 , i 0 = 12 24 , T = 0.06 , h t = 0.002 , B = 5 .
Example 1.
The exact permeability k in this example is shown in Figure 2. In order to illustrate the noise sensitivity, we add Gaussian noise to the measured data Γ ¯ . With 5%, 10%, 15%, and 20% Gaussian noises, the identified permeability results by the homotopy method with constraints are shown in Figure 3, Figure 4, Figure 5 and Figure 6, respectively.
Example 2.
This example selects a model of two anomalous bodies in a homogeneous medium with a permeability of 4.25, and the anomalous bodies have the permeability of 6.18 and 7.55, see Figure 7. Figure 8, Figure 9, Figure 10 and Figure 11, respectively, show the identified permeability results of the homotopy method with constraints, with 5%, 10%, 15%, and 20% Gaussian noises.
To test the necessity of introduction of constraints and the wide convergence of the homotopy strategy, the homotopy method without constraints and basic iterative method with constraints are used for the above measured data, and the relative errors of all the identified permeability results by the three methods are listed in Table 2, where × means that there is no convergent result.
Table 2 shows that
(1)
The stability of the homotopy method with constraints is better than the homotopy method without constraints;
(2)
The region of convergence of the homotopy strategy is wider than the basic iterative method with constraints;
(3)
The homotopy method with constraints has wide convergence region and strong anti-noise ability.

5. Conclusions

To restrain the noise and improve identified model quality, the constraint condition is introduced to the inverse problem of the nonlinear convection-diffusion equation in the multiphase porous media flow. Then, Tikhonov regularization is applied to the constrained optimization problem obtained by the discretization to overcome the ill-posed property. We developed a homotopy method that is widely convergent. The results from the numerical experiments show that this method is feasible and effective. Compared with the homotopy method without constraints, our approach has better stability; compared with the basic iterative method with constraints, our approach has a wider convergence region. So far, there is no literature on the use of fractional operators for constraint inverse problems, which will be an interesting future work.

Author Contributions

Conceptualization, T.L. and Y.Z.; methodology, T.L. and Y.Z.; software, T.L. and K.X.; validation, T.L. and K.X.; formal analysis, T.L. and K.X.; investigation, T.L. and K.X.; resources, Y.Y., R.Q. and Y.Q.; data curation, T.L.; writing—original draft preparation, T.L., K.X. and Y.Z.; writing—review and editing, T.L., K.X. and Y.Z.; visualization, T.L.; supervision, Y.Z.; project administration, T.L., Y.Q. and C.L.; funding acquisition, T.L., Y.Q. and C.L. All authors have read and agreed to the published version of the manuscript.

Funding

This research was funded by the Natural Science Foundation of He-bei Province of China (A2020501005, A2020501007), the Fundamental Research Funds for the Central Universities (N2123015).

Institutional Review Board Statement

Not applicable.

Informed Consent Statement

Not applicable.

Data Availability Statement

Not applicable.

Acknowledgments

On behalf of my co-authors, we would like to express our great appreciation to the editor and reviewers.

Conflicts of Interest

The authors declare no conflict of interest.

Appendix A. Discretization of the Diffusion Term

Let
D x u i , j k = u i , j k u i 1 , j k h x , D + x u i , j k = u i + 1 , j k u i , j k h x ,
be the discrete derivatives of x, and a corresponding notation is used for the the discrete derivatives of y.
The mean values for the discretized permeability are defined as
k ¯ i , j + 1 2 x = 1 2 ( k i + 1 2 , j + 1 2 + k i 1 2 , j + 1 2 ) , k ¯ i + 1 2 , j y = 1 2 ( k i + 1 2 , j + 1 2 + k i + 1 2 , j 1 2 ) ,
and the mean values for the nonlinear diffusion function are defined as
( N ¯ x ) i + 1 2 , j k = 1 2 ( N i + 1 , j k + N i , j k ) , ( N ¯ y ) i , j + 1 2 k = 1 2 ( N i , j + 1 k + N i , j k ) ,
where N i , j k = N ( u i , j k ) .
Then, the diffusion term can be discretizated as
· ( k i , j N i , j k u i , j k ) = D x ( k ¯ i + 1 2 , j y ( N ¯ x ) i + 1 2 , j k D + x u i , j k ) + D y ( k ¯ i , j + 1 2 x ( N ¯ y ) i , j + 1 2 k D + y u i , j k ) .

Appendix B. Discretization of the Convection Term

The convection term can be discretizated by the Engquist–Osher upwind scheme [46] as follows:
· ( ϕ ( u i , j k 1 ) , φ ( u i , j k 1 ) ) = D x ϕ E O ( u i , j k 1 , u i + 1 , j k 1 ) + D y φ E O ( u i , j k 1 , u i , j + 1 k 1 ) ,
where the Engquist–Osher numerical flux functions ϕ E O ( u i , j k 1 , u i + 1 , j k 1 ) and φ E O ( u i , j k 1 , u i , j + 1 k 1 ) are defined by
ϕ E O ( u i , j k 1 , u i + 1 , j k 1 ) = 1 2 ( ϕ ( u i , j k 1 ) + ϕ ( u i + 1 , j k 1 ) ) 1 2 u i , j k 1 u i + 1 , j k 1 | ϕ ( ξ ) | d ξ , φ E O ( u i , j k 1 , u i , j + 1 k 1 ) = 1 2 ( φ ( u i , j k 1 ) + φ ( u i , j + 1 k 1 ) ) 1 2 u i , j k 1 u i , j + 1 k 1 | φ ( ξ ) | d ξ .
Finally, we give the explicit formulas for ϕ E O and φ E O , for examples of ϕ ( u ) = u 2 ( 1 5 ( 1 u 2 ) ) u 2 + ( 1 u ) 2 and φ ( u ) = u 2 u 2 + ( 1 u ) 2 . For the sake of simplicity, the following calculations only use one subscript index.

Appendix B.1. Buckley–Leverett Flux Function

φ ( u ) is an S-shaped flux function of Buckley–Leverett type, and
φ ( u ) 0 , u ( 0 , 1 ) ,
therefore
φ E O ( u i , u i + 1 ) = φ ( u i ) .

Appendix B.2. Buckley–Leverett Flux Function with Gravitational Effects

ϕ ( u ) is an S-shaped flux function of Buckley–Leverett type under gravity effects, and for the calculation of u i u i + 1 | ϕ ( ξ ) | d ξ , we first discuss the sign of ϕ in ( 0 , 1 ) . Actually, ϕ has only one zero point in ( 0 , 1 ) , that is
x z e r o 0.37 .
Therefore, we can get
ϕ ( ξ ) < 0 , ξ ( 0 , x z e r o ) , = 0 , ξ = x z e r o , > 0 , ξ ( x z e r o , 1 ) .
Additionally, we have the formula of the integral
u i u i + 1 | ϕ ( ξ ) | d ξ = ϕ ( u i + 1 ) ϕ ( u i ) , u i , u i + 1 x z e r o , ϕ ( u i ) ϕ ( u i + 1 ) , u i , u i + 1 < x z e r o , ϕ ( u i ) + ϕ ( u i + 1 ) 2 ϕ ( x z e r o ) , u i < x z e r o , u i + 1 x z e r o , 2 ϕ ( x z e r o ) ϕ ( u i ) ϕ ( u i + 1 ) , u i x z e r o , u i + 1 < x z e r o ,
and according to the definition of ϕ E O , then we have
ϕ E O ( u i , u i + 1 ) = ϕ ( u i ) , u i , u i + 1 x z e r o , ϕ ( u i + 1 ) , u i , u i + 1 < x z e r o , ϕ ( x z e r o ) , u i < x z e r o , u i + 1 x z e r o , ϕ ( u i ) + ϕ ( u i + 1 ) ϕ ( x z e r o ) , u i x z e r o , u i + 1 < x z e r o .

References

  1. Rana, B.M.J.; Arifuzzaman, S.M.; Islam, S.; Reza-E-Rabbi, S.; Al-Mamun, A.; Mazumder, M.; Roy, K.C.; Khan, M.S. Swimming of microbes in blood flow of nano-bioconvective Williamson fluid. Therm. Sci. Eng. Progr. 2021, 25, 101018. [Google Scholar] [CrossRef]
  2. Reza-E-Rabbi, S.; Ahmmed, S.F.; Arifuzzaman, S.M.; Sarkar, T.; Khan, M.S. Computational modelling of multiphase fluid flow behaviour over a stretching sheet in the presence of nanoparticles. Eng. Sci. Technol. Int. J. 2020, 23, 605–617. [Google Scholar] [CrossRef]
  3. Moosavian, N. Pipe network modeling for analysis of flow in porous media. Can. J. Civil Eng. 2019, 46, 1151–1159. [Google Scholar] [CrossRef]
  4. Shapiro, A.A. Mechanics of the separating surface for a two-phase co-current flow in a porous medium. Transp. Porous Med. 2016, 112, 489–517. [Google Scholar] [CrossRef] [Green Version]
  5. Hunt, A.G.; Sahimi, M. Flow, transport, and reaction in porous media: Percolation scaling, critical-path analysis, and effective medium approximation. Rev. Geoph. 2017, 55, 993–1078. [Google Scholar] [CrossRef]
  6. Espedal, M.S.; Karlsen, K.H. Numerical solution of reservoir flow models based on large time step operator splitting algorithm. In Filtration in Porous Media and Industrial Applications: Lecture Notes in Mathematics; Fasano, A., Ed.; Springer: Berlin/Heidelberg, Germany, 2000; Volume 1734, pp. 9–77. [Google Scholar]
  7. Saeed, S.T.; Riaz, M.B.; Baleanu, D.; Akgul, A.; Husnine, S.M. Exact analysis of second grade fluid with generalized boundary conditions. Intell. Autom. Soft Comput. 2021, 28, 547–559. [Google Scholar] [CrossRef]
  8. Saeed, S.T.; Abro, K.A.; Almani, S. Role of single slip assumption on the viscoelastic liquid subject to non-integer differentiable operators. Math. Meth. Appl. Sci. 2021, 44, 6005–6020. [Google Scholar] [CrossRef]
  9. Saeed, S.T.; Riaz, M.B.; Baleanu, D.; Abro, K.A. A mathematical study of natural convection flow through a channel with non-singular kernels: An application to transport phenomena. Alex. Eng. J. 2020, 59, 2269–2281. [Google Scholar] [CrossRef]
  10. Abdeljawad, T.; Riaz, M.B.; Saeed, S.T.; Iftikhar, N.I. MHD Maxwell fluid with heat transfer analysis under ramp velocity and ramp temperature subject to non-integer differentiable operators. Comput. Model. Eng. Sci. 2021, 126, 821–841. [Google Scholar] [CrossRef]
  11. Firdous, H.; Saeed, S.T.; Ahmad, H.; Askar, S. Using non-Fourier’s heat flux and non-Fick’s mass flux theory in the radiative and chemically reactive flow of Powell-Eyring fluid. Energies 2021, 14, 6882. [Google Scholar] [CrossRef]
  12. Riaz, M.B.; Siddique, I.; Saeed, S.T.; Atangana, A. MHD Oldroyd-B fluid with slip condition in view of local and nonlocal kernels. J. Appl. Comput. Mech. 2021, 7, 116–127. [Google Scholar]
  13. Riaz, M.B.; Saeed, S.T. Comprehensive analysis of integer order, Caputo-Fabrizio and Atangana-Baleanu fractional time derivative for MHD Oldroyd-B fluid with slip effect and time dependent bounday conditions. Discr. Contin. Dynam. Syst. 2021, 14, 3719–3746. [Google Scholar] [CrossRef]
  14. Hazra, S.; Class, H.; Helmig, R.; Schulz, V. Forward and inverse problems in modeling of multiphase flow and transport through porous media. Computat. Geosci. 2004, 8, 21–47. [Google Scholar] [CrossRef]
  15. Wang, J.; Zabaras, N. A Markov random field model of contamination source identification in porous media flow. Int. J. Heat Mass Transf. 2006, 49, 939–950. [Google Scholar] [CrossRef]
  16. Nilssen, T.K.; Karlsen, K.H.; Mannseth, T.; Tai, X.C. Identification of diffusion parameters in a nonlinear convection-diffusion equation using the augmented lagrangian method. Comput. Geosci. 2009, 13, 317–329. [Google Scholar] [CrossRef] [Green Version]
  17. Yan, B.; Harp, D.R.; Chen, B.; Pawar, R. A physics-constrained deep learning model for simulating multiphase flow in 3D heterogeneous porous media. Fuel 2022, 313, 122693. [Google Scholar] [CrossRef]
  18. Yan, B.; Chen, B.; Harp, D.R.; Jia, W.; Pawar, R.J. A robust deep learning workflow to predict multiphase flow behavior during geological CO2 sequestration injection and Post-Injection periods. J. Hydrol. 2022, 607, 127542. [Google Scholar] [CrossRef]
  19. Magzymov, D.; Ratnakar, R.R.; Dindoruk, B.; Johns, R.T. Evaluation of machine learning methodologies using simple physics based conceptual models for flow in porous media. In Proceedings of the SPE Annual Technical Conference and Exhibition, Dubai, United Arab Emirates, 21–23 September 2021. [Google Scholar]
  20. Almajid, M.M.; Abu-Al-Saud, M.O. Prediction of porous media fluid flow using physics informed neural networks. J. Petrol. Sci. Eng. 2022, 208, 109205. [Google Scholar] [CrossRef]
  21. Jackson, D.D. Interpretation of inaccurate, insufficient and inconsistent data. Geophys. J. R. Astron. Soc. 1972, 28, 97–109. [Google Scholar] [CrossRef] [Green Version]
  22. Watson, L.T. Globally convergent homotopy methods: A tutorial. Appl. Math. Comput. 1989, 31, 369–396. [Google Scholar]
  23. Mousa, M.M.; Alsharari, F. Convergence and error estimation of a new formulation of homotopy perturbation method for classes of nonlinear integral/integro-differential equations. Mathematics 2021, 9, 2244. [Google Scholar] [CrossRef]
  24. Agarwal, P.; Akbar, M.; Nawaz, R.; Jleli, M. Solutions of system of Volterra integro-differential equations using optimal homotopy asymptotic method. Math. Meth. Appl. Sci. 2021, 44, 2671–2681. [Google Scholar] [CrossRef]
  25. Mousa, M.M.; Kaltayev, A. Homotopy perturbation method for solving nonlinear differential-difference equations. Z. Naturforsch. A 2010, 65, 511–517. [Google Scholar] [CrossRef]
  26. Hammad, H.A.; Agarwal, P.; Guirao, J.L.G. Applications to boundary value problems and homotopy theory via tripled fixed point techniques in partially metric spaces. Mathematics 2021, 9, 2012. [Google Scholar] [CrossRef]
  27. Mousa, M.M.; Kaltayev, A. Application of the homotopy perturbation method to a magneto-elastico-viscous fluid along a semi-infinite plate. Int. J. Nonlinear Sci. Numer. Simul. 2009, 10, 1113–1120. [Google Scholar] [CrossRef]
  28. Saad, K.M.; Iyiola, O.S.; Agarwal, P. An effective homotopy analysis method to solve the cubic isothermal auto-catalytic chemical system. AIMS Math. 2018, 3, 183–194. [Google Scholar] [CrossRef]
  29. Mallick, A.; Ranjan, R.; Das, R. Application of homotopy perturbation method and inverse prediction of thermal parameters for an annular fin subjected to thermal load. J. Therm. Stress. 2016, 39, 298–313. [Google Scholar] [CrossRef]
  30. Mallick, A.; Ranjan, R.; Prasad, D.K.; Das, R. Inverse prediction and application of homotopy perturbation method for efficient design of an annular fin with variable thermal conductivity and heat generation. Math. Model. Anal. 2016, 21, 699–717. [Google Scholar] [CrossRef]
  31. Biswal, U.; Chakraverty, S.; Ojha, B.K. Application of homotopy perturbation method in inverse analysis of Jeffery-Hamel flow problem. Eur. J. Mech. B Fluid. 2021, 86, 107–112. [Google Scholar] [CrossRef]
  32. Sattari Shajari, P.; Shidfar, A. Application of weighted homotopy analysis method to solve an inverse source problem for wave equation. Inverse Probl. Sci. Eng. 2019, 27, 61–88. [Google Scholar] [CrossRef]
  33. Liu, T. A multigrid-homotopy method for nonlinear inverse problems. Comput. Math. Appl. 2020, 79, 1706–1717. [Google Scholar] [CrossRef]
  34. Liu, T. A wavelet multiscale-homotopy method for the parameter identification problem of partial differential equations. Comput. Math. Appl. 2016, 71, 1519–1523. [Google Scholar] [CrossRef]
  35. Hu, L.; Huang, L.; Lu, Z.R. Crack identification of beam structures using homotopy continuation algorithm. Inverse Probl. Sci. Eng. 2017, 25, 169–187. [Google Scholar] [CrossRef]
  36. Courbot, J.B.; Colicchio, B. A fast homotopy algorithm for gridless sparse recovery. Inverse Probl. 2021, 37, 025002. [Google Scholar] [CrossRef]
  37. Słota, D.; Chmielowska, A.; Brociek, R.; Szczygieł, M. Application of the homotopy method for fractional inverse Stefan problem. Energies 2020, 13, 5474. [Google Scholar] [CrossRef]
  38. Liu, J.; Wang, B. Solving the backward heat conduction problem by homotopy analysis method. Appl. Numer. Math. 2018, 128, 84–97. [Google Scholar] [CrossRef]
  39. Liu, T. Porosity reconstruction based on Biot elastic model of porous media by homotopy perturbation method. Chaos Soliton. Fract. 2022, 158, 112007. [Google Scholar] [CrossRef]
  40. Enting, I.G.; Pearman, G.I. Description of a one-dimensional carbon cycle model calibrated using techniques of constrained inversion. Tellus B 1987, 39, 459–476. [Google Scholar] [CrossRef] [Green Version]
  41. Lambert, M.; Lesselier, D. Binary-constrained inversion of a buried cylindrical obstacle from complete and phaseless magnetic fields. Inverse Probl. 2000, 16, 563–576. [Google Scholar] [CrossRef]
  42. Atzberger, C.; Richter, K. Spatially constrained inversion of radiative transfer models for improved LAI mapping from future Sentinel-2 imagery. Remote Sens. Environ. 2012, 120, 208–218. [Google Scholar] [CrossRef]
  43. Auken, E.; Christiansen, A.V. Layered and laterally constrained 2D inversion of resistivity data. Geophysics 2004, 69, 752–761. [Google Scholar] [CrossRef]
  44. Zhao, J.; Liu, T.; Liu, S. An adaptive homotopy method for permeability estimation of a nonlinear diffusion equation. Inverse Probl. Sci. Eng. 2013, 21, 585–604. [Google Scholar] [CrossRef]
  45. Liu, T. Parameter estimation with the multigrid-homotopy method for a nonlinear diffusion equation. J. Comput. Appl. Math. 2022, 413, 114393. [Google Scholar] [CrossRef]
  46. Engquist, B.; Osher, S. One-sided difference approximations for nonlinear conservation laws. Math. Comp. 1981, 36, 321–351. [Google Scholar] [CrossRef]
Figure 1. Practical description of the problem.
Figure 1. Practical description of the problem.
Processes 10 01143 g001
Figure 2. The exact permeability k in Example 1. (a) Three-dimensional model. (b) Two-dimensional profile.
Figure 2. The exact permeability k in Example 1. (a) Three-dimensional model. (b) Two-dimensional profile.
Processes 10 01143 g002
Figure 3. The identified permeability k with 5% Gaussian noise in Example 1. (a) Three-dimensional model. (b) Two-dimensional profile.
Figure 3. The identified permeability k with 5% Gaussian noise in Example 1. (a) Three-dimensional model. (b) Two-dimensional profile.
Processes 10 01143 g003
Figure 4. The identified permeability k with 10% Gaussian noise in Example 1. (a) Three-dimensional model. (b) Two-dimensional profile.
Figure 4. The identified permeability k with 10% Gaussian noise in Example 1. (a) Three-dimensional model. (b) Two-dimensional profile.
Processes 10 01143 g004
Figure 5. The identified permeability k with 15% Gaussian noise in Example 1. (a) Three-dimensional model. (b) Two-dimensional profile.
Figure 5. The identified permeability k with 15% Gaussian noise in Example 1. (a) Three-dimensional model. (b) Two-dimensional profile.
Processes 10 01143 g005
Figure 6. The identified permeability k with 20% Gaussian noise in Example 1. (a) Three-dimensional model. (b) Two-dimensional profile.
Figure 6. The identified permeability k with 20% Gaussian noise in Example 1. (a) Three-dimensional model. (b) Two-dimensional profile.
Processes 10 01143 g006
Figure 7. The exact permeability k in Example 2. (a) Three-dimensional model. (b) Two-dimensional profile.
Figure 7. The exact permeability k in Example 2. (a) Three-dimensional model. (b) Two-dimensional profile.
Processes 10 01143 g007
Figure 8. The identified permeability k with 5% Gaussian noise in Example 2. (a) Three-dimensional model. (b) Two-dimensional profile.
Figure 8. The identified permeability k with 5% Gaussian noise in Example 2. (a) Three-dimensional model. (b) Two-dimensional profile.
Processes 10 01143 g008
Figure 9. The identified permeability k with 10% Gaussian noise in Example 2. (a) Three-dimensional model. (b) Two-dimensional profile.
Figure 9. The identified permeability k with 10% Gaussian noise in Example 2. (a) Three-dimensional model. (b) Two-dimensional profile.
Processes 10 01143 g009
Figure 10. The identified permeability k with 15% Gaussian noise in Example 2. (a) Three-dimensional model. (b) Two-dimensional profile.
Figure 10. The identified permeability k with 15% Gaussian noise in Example 2. (a) Three-dimensional model. (b) Two-dimensional profile.
Processes 10 01143 g010
Figure 11. The identified permeability k with 20% Gaussian noise in Example 2. (a) Three-dimensional model. (b) Two-dimensional profile.
Figure 11. The identified permeability k with 20% Gaussian noise in Example 2. (a) Three-dimensional model. (b) Two-dimensional profile.
Processes 10 01143 g011
Table 1. Parameter definitions.
Table 1. Parameter definitions.
ParameterDefinition
k absolute permeability
V total Darcy velocity
ς total mobility of phases
pglobal pressure
ϑ density of wetting phase
hheight
β porosity
f 1 injection well
f 2 production well
ς r a mobility ratio
k r a relative permeability
μ viscosity
ρ density
ϵ phase mobility of nonwetting phase
Table 2. Relative errors of identified permeability k by three different methods. HM—homotopy method with constraints, HM—homotopy method without constraints, BIMC—basic iterative method with constraints.
Table 2. Relative errors of identified permeability k by three different methods. HM—homotopy method with constraints, HM—homotopy method without constraints, BIMC—basic iterative method with constraints.
Example NumberNoise LevelHMCHMBIMC
4.15%6.22%6.90%×
10%6.38%7.40%×
15%7.13%××
20%7.44%××
4.25%5.71%6.47%×
10%5.81%7.39%×
15%6.06%××
20%7.19%××
Publisher’s Note: MDPI stays neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Share and Cite

MDPI and ACS Style

Liu, T.; Xia, K.; Zheng, Y.; Yang, Y.; Qiu, R.; Qi, Y.; Liu, C. A Homotopy Method for the Constrained Inverse Problem in the Multiphase Porous Media Flow. Processes 2022, 10, 1143. https://doi.org/10.3390/pr10061143

AMA Style

Liu T, Xia K, Zheng Y, Yang Y, Qiu R, Qi Y, Liu C. A Homotopy Method for the Constrained Inverse Problem in the Multiphase Porous Media Flow. Processes. 2022; 10(6):1143. https://doi.org/10.3390/pr10061143

Chicago/Turabian Style

Liu, Tao, Kaiwen Xia, Yuanjin Zheng, Yanxiong Yang, Ruofeng Qiu, Yunfei Qi, and Chao Liu. 2022. "A Homotopy Method for the Constrained Inverse Problem in the Multiphase Porous Media Flow" Processes 10, no. 6: 1143. https://doi.org/10.3390/pr10061143

APA Style

Liu, T., Xia, K., Zheng, Y., Yang, Y., Qiu, R., Qi, Y., & Liu, C. (2022). A Homotopy Method for the Constrained Inverse Problem in the Multiphase Porous Media Flow. Processes, 10(6), 1143. https://doi.org/10.3390/pr10061143

Note that from the first issue of 2016, this journal uses article numbers instead of page numbers. See further details here.

Article Metrics

Back to TopTop