Theoretical Estimate of the Total Shelf Length in a Gas Fields Model

• Alexander K. Skiba and Nadezhda K. Skiba
From the journal Open Computer Science

Abstract

The authors consider a continuous dynamic approximating model of a gas fields group and, on its basis, set maximum and minimum issues. The tasks proposed for research are optimal control problems with mixed constraints with free-final-time and moving right end. We analytically solve the rapid-action problem. The central mathematical apparatus is Pontryagin maximum principle in Arrow form, using Lagrange multipliers. The theoretically obtained results of the analysis are of particular interest.

1 Introduction

Natural gas refers to non-renewable minerals located under the surface of the earth at high pressure and temperature. It is located in the earth's layer in a gaseous state in the form of separate accumulations (gas deposits). The composition of natural gas includes many useful substances used in industry.

We extract natural gas from the earth using wells. We try to place the wells evenly throughout the field. Thus, we ensure the same pressure drop across an entire deposit [1].

Natural gas is a valuable product with an enormous demand. Russia is one of the wealthiest countries in the world in terms of natural gas reserves. Proven natural gas reserves in Russia for 2014 amount to 24.6% of world reserves [2].

However, the locations of gas production and its consumption are at a considerable distance. One of the most effective ways to deliver gas to consumers is by transporting it through a gas pipeline.

Gas pipeline capacity Q¯ limits the volume of gas flowing through it. There are many problems for gas workers in the fields. The most important task is the execution of the current gas production equal to the value Q¯ .

We introduce the following notation:

• by i we denote an ordinal number of the deposit in a group consisting of n elements;

• by N¯i we denote a total stock of producing wells in the corresponding field;

• by Ni(t) we denote an operating wells stock at time t;

• by qi(t) we denote an average flow rate of producing wells at time t.

If at time t the following strict inequality holds

(1) i=1nqi(t)N¯i>Q¯,

then there exists such an N^i=Ni(t)[0,N¯i] at which the equality

(2) i=1nqi(t)N^i=Q¯,

is satisfied.

There are many such solutions. Therefore, on this set, it is possible to solve other both static and dynamic problems. As the flow rate of producing wells decreases, there comes a moment T at which the equality

(3) i=1nqi(T)N¯i=Q¯.

is fulfilled.

Further, we can no longer fulfil the plan for gas production. Therefore, the value T is the shelf length of the gas fields group. The question arises: within what limits can be the shelf length of the gas deposits group be? Thus, it is necessary to solve two optimal control problems with a mixed restriction on the maximum and minimum. We can analytically solve these problems and obtain theoretically substantiated estimates of the gas fields shelf length. We devote this paper to this statement. As the central mathematical apparatus, we choose the Pontryagin maximum principle in the Arrow form.

2 Pontryagin maximum principle in Arrow form

About 50 years ago, K. Arrow, the 1972 Nobel Prize winner in economics, published an article [3], where he, taking Pontryagin maximum principle [4] as a basis, modifies it and formulates propositions that allow solving problems optimal control with mixed constraints. The proposals additionally include some elements of nonlinear programming, such as Lagrangian function, Lagrange multipliers, and the complementary slackness conditions. The author of this article used a modified principle in his papers [5, 6].

Among the Russian works on optimal control, we can distinguish monographs by A. M. Ter-Krikorov [7], S. M. Aseev, and A.V. Kryazhemsky [8]. An interesting international paper is an article by E. Balder [9], in which the author proves the existence theorem for optimal control problems with an infinite horizon.

We formulate the Arrow proposition.

Proposition 1

Let v*(t) be a choice of instruments (0 ≤ tT) which maximizes 0TU[x(t),v(t),t]dt subject to the condition,

• (a) x˙=S[x(t),v(t),t],0tT ,

and subject to a set of constraints,
• (b) F[x(t), v(t), t] ≥ 0

At least one constraint in (b) contains the instruments v(t) at t ∈ [0, T]. The constraint set (b) may or may not include the state variables x(t) at t ∈ [0, T], the initial conditions x(t) at t = 0, and the non-negative terminal conditions x(t) at t = T.

If the Constraint Qualification (Regularity Conditions) holds, then there exist auxiliary variables p(t), such that, for each t,

• (c) v*(t) maximizes H[x(t), v(t), p(t), t] subject to the constaints (b), where H(x, v, p, t) = U(x, v, t) + pS(x, v, t);

• (d) p˙i=L/xi , when evaluated at x = x(t), v = v*(t), p = p(t), where

• (e) L(x, v, p, t) = H(x, v, p, t) + qF(x, v, t), and the Lagrange multipliers q are such that

• (f) ∂L/∂vk = 0, for x = x(t), v = v*(t), p = p(t), q(t) ≥ 0, q(t)F[x(t), v*(t), t] = 0, and

• (g) p(T) ≥ 0, p(T)x(T) = 0.

Note that if a set of constraints (b) consists only of inequalities

• (b) F[v(t)] ≥ 0 at t ∈ [0, T],

then the optimal control problem posed in Proposition 1 can be solved using the maximum principle formulated by L.S. Pontryagin.

3 Model Description and Problems Posing

Let us consider a model for the functioning of a gas fields group with interacting wells [10, 11]. By Vi(t) we denote the recoverable gas reserve at time t. Between variables, we establish a differential relationship, which we describe as a system of differential equations:

(4) V˙i=Niqi,q˙i=qi0Vi0qi(t)Ni(t)wheni=1,n¯,

under the initial conditions Vi0>0,qio>0 . We impose a limit of 0NiNi¯ on the stock of operating wells. Here Ni¯>0 . We assume an even distribution of production wells over the entire area of each field. We suppose the fulfilment of equalities for all values of t:

(5) qi(t)=qi0Vi0Vi(t)wheni=1,n¯,

and limt→⧜ Vi(t) = 0 when i=1,n¯ . By differentiating the equalities in (5), we get the differential equations (4). We also suppose different values for qi0Vi0N¯i .

There is a general limitation on the total shelf length of the gas deposits given by

i=1nqi(t)Ni(t)Q¯.

We calculate cumulative gas production using the following formula:

i=1n0Tqi(t)Ni(t)dt.

The paper [6] considered two questions: the problem of maximizing the accumulated extraction, and the problem of maximizing the total shelf length for the gas fields group. We present their mathematical formulations.

Problem 1

We maximize the functional

(6) i=1n0Tqi(t)Ni(t)dt

on the fixed interval [0, T] under the differential connections

(7) q˙i=αi0qi(t)Ni(t),0tT,i=1,n¯

with the initial conditions

(8) qi0>0,Vi0>0,i=1,n¯,

where

(9) αi0=qi0Vi0.

We impose the following restrictions on the controls:

(10) 0Ni(t)Ni¯,0tT,i=1,n¯,

(11) i=1nqi(t)Ni(t)Q¯,0tT.

We also assume that

(12) i=1nqi0N¯i>Q¯.

Here qi(t) for i=1,n¯ are phase variables, and Q¯ is a constant value.

The controls Ni(t) for i=1,n¯ belong to a set of measurable functions. The right end of the phase trajectory is free. The quantities qi0Vi0N¯i differ between themselves.

Problem 2

We maximize the functional T under the differential connections (7), with initial conditions (8) and constraints (10) and

(13) i=1nqi(t)Ni(t)=Q¯,t[0,T].

We assume that the initial conditions satisfy the strict inequality (12). As in Problem 1, the controls Ni(t) when i=1,n¯ , belong to a set of measurable functions. The following statements hold under these assumptions.

Proposition 2

For any admissible trajectory from sets (7), (8), and (10) the following strict inequalities hold: qi(T) > 0, i=1,n¯ .

Proposition 2 follows from the differential equations in (7), the initial conditions (8) and restrictions on controls (10).

Proposition 3

For any admissible trajectory from the sets (7), (8), (10), and (13), and under the initial constraint (12), there exists T > 0 for which the following equality holds:

(14) i=1nqi(T)N¯i=Q¯.

The point q(T) is the right end of an admissible trajectory q(t). Further development of the system (7) is impossible due to existing limitations.

Proof

Multiply both sides of each differential equation in (7) by the corresponding value N¯i , then sum up and integrate both sides. After simple transformations, we get

(15) i=1nqi(t)N¯i=i=1nqi0N¯ii=1nαi0N¯i0tqi(t)Ni(t)dt.

Taking into account (13), we continue the transformation to yield

(16) i=1nqi(t)N¯i=i=1nqi0N¯iα10N¯1Q¯ti=2n(αi0N¯iα10N¯1)0tqi(t)Ni(t)dt.

Let α10N¯1 have the smallest value among the other values αi0N¯i for i=2,n¯ , then

(17) i=1nqi(t)N¯ii=1nqi0N¯iα10N¯1Q¯t.

The right-hand side of the inequality (17) reaches Q¯ in a finite time t′. Therefore, its left-hand side achieves the value Q¯ in the time T ∈ (0, t′). Thus, at time T, equality (14) holds.

Proposition 4

Let the equality (14) holds on some trajectory at time T, then

(18) i=1nqi(t)N¯i<Q¯

at T < t.

The following proposition follows from Propositions 3 and 4.

Proposition 5

For any admissible trajectory from the sets (7), (8), (10), and (13), and under the initial constraint (12), the control on the right end is strictly definite, i.e. Ni(T)=N¯i when i=1,n¯ .

The exact upper bound for the total length of the shelf can be obtained by solving the Problem 2. We have described the solution in [6]. The essence of the solution is as follows. We sort all deposits in ascending order of values αi0N¯i . Further, we consistently bring into development all the existing fields from the group. The solution to the next problem leads us to obtain a lower bound. Below we will describe in detail its complete solution.

Problem 3

It is required to minimize the functional T under the differential connections in (7), with the initial conditions (8) and with the constraints (10), (13).

We assume that the initial conditions satisfy the strict inequality (12). As in Problems 1 and 2, controls Ni(t) when i=1,n¯ , belong to a set of measurable functions.

In [6], the search for a solution to Problem 2 is carried out indirectly. First, we solve Problem 1. Then, by taking limits, we find the solution to Problem 2. However, the search for a solution to Problem 2 can be carried out directly using the Arrow proposition and the corresponding transversality conditions.

The transversality conditions applied below are not mentioned in the paper [3] along with the Arrow proposition. Nevertheless, they are used in solving similar problems of optimal control together with the Pontryagin maximum principle in the classical formulation [12].

We now reformulate Problems 2 and 3 as follows.

Problem 2a

It is required to find optimal controls Ñi when i=1,n¯ of the dynamic system (7) from the class of admissible controls, which translates the dynamic system from the given initial state (8) to the hyperplane (14) for the maximum time T under the restrictions (10) and (13) on controls and phase variables along the trajectories.

Problem 3a

It is required to find optimal controls Ñi for i=1,n¯ of the dynamic system (7) from the class of admissible controls, which translates the dynamic system from the given initial state (8) to the hyperplane (14) for the shortest time T under the restrictions (10) and (13) on phase variables and controls along the trajectories.

Problems 2a and 3a are optimal control problems with a moving right end and free-final-time. Moreover, Problem 3a refers to the optimal rapid action problem. The optimal solution for Problems 2a and 3a exists. This fact, for example, follows from the monograph [13]. In search of decisions to Problems 2a and 3a, there is a logical relationship. The existing differences are insignificant.

In [6] we have already found the optimal solution to Problem 2. Therefore, in this paper, we restrict ourselves only to searching for the optimal decision to Problem 3a.

4 Investigation of Problem 3a

We recall

(19) αi0=qi0Vi0,i=1,n¯.

Next, we order the phase variables in descending order of αi0N¯i . We note that this order is unique since the values αi0N¯i are different. According to Proposition 1, we write out the Hamiltonian and the Lagrangian as follows

(20) H(q,N,ψ)=ci=1nψiαi0qiNi,

(21) L(q,N,ψ,γ,δ,β)=c+i=1n[ψiαi0qiNi+γi(Ni¯Ni)+δiNi]+β[Q¯i=1n(qiNi)],

where

(22) c=1.

Next, we describe a set of admissible controls G(q), on which we maximize Hamiltonian (20) as follows

(23) G(q)={Nn|0NN¯,i=1nqiNi=Q¯}.

The main objective of the study is to search a continuous vector function ψ(t) and control Ñ(t), which satisfy the adjoint system of differential equations (27) and transversality conditions:

(24) H(q˜(T),N˜(T),ψ(T))=ci=1nψi(T)αi0q˜i(T)N¯i=0;

(25) ψi(T)=νN¯i,i=1,n¯,

where ν is some number. For each t ∈ [0, T] the control Ñ(t) maximizes the Hamiltonian

(26) H(q˜(t),N˜(t),ψ(t))=maxNG(q˜(t))H(q˜(t),N,ψ(t))=maxNG(q˜(t))i=1n[cαi0ψi(t)q˜i(t)Ni].

We represent the system of equations as

ψ˙i(t)=L/qi=αi0ψi(t)N˜i(t)+β(t)N˜i(t),i=1,n¯,

or

(27) ψ˙i(t)=[αi0ψi(t)+β(t)]N˜i(t),i=1,n¯.

According to Proposition 1, the Lagrange multipliers β(t), γ(t), δ(t) must satisfy the equalities:

L/Ni=αi0ψi(t)q˜i(t)β(t)q˜i(t)γi(t)+δi(t)=0,i=1,n¯,

or

(28) q˜i(t)[αi0ψi(t)+β(t)]=γi(t)δi(t),i=1,n¯;

(29) β(t)[Q¯i=1nq˜i(t)N˜i(t)]=0;

(30) γi(t)[Ni¯N˜i(t)]=0,γi(t)0,i=1,n¯;

(31) δi(t)N˜i(t)=0,δi(t)0,i=1,n¯.

In nonlinear programming, the above relations (29), (30) and (31) are complementary slackness conditions. We introduce the notation φi(t)=αi0ψi(t) for i=1,n¯ . In our new notation, we represent the maximization of the Hamiltonian (26), the system of adjoint equations (27), equality (28), transversality conditions (24) and (25) in the following form:

(32) H(q˜(t),N˜(t),φ(t))=maxNG(q˜(t))H(q˜(t),N,φ(t))=c+maxNG(q˜(t))i=1nφi(t)q˜i(t)Ni;

(33) φ˙i(t)=[φi(t)β(t)]αi0N˜i(t),i=1,n¯;!

(34) q˜i(t)[φi(t)β(t)]=γi(t)δi(t),i=1,n¯;

(35) H(q˜(T),N˜(T),φ(T))=c+i=1nφi(T)q˜i(T)N¯i=0;

(36) φi(T)=ναi0N¯i,i=1,n¯.

Let us consider the transversality conditions in more detail. It follows from (24) and (25) that the right end of the adjoint trajectory ψ(t) lies in the negative region. However, if we solved Problem 2a, then ψ(T) would belong to the positive domain. Using the transversality conditions (24) and (25), and taking into account (22), (35), and (36), we obtain

(37) ν=1i=1nαi0q˜i(T)N¯i2,

(38) φi(T)=αi0N¯ii=1nαi0q˜i(T)N¯i2,i=1,n¯.

This implies the following proposition.

Proposition 6

The right end of the “adjoint” vector function φ(t) lies in the positive domain; that is, φi(T) > 0 when i=1,n¯ . The values of vector components φi(T)when i=1,n¯ decrease.

We do not know the value of the vector q˜i(T) . However, if we knew it, then the coefficient ν and, accordingly, the values of the vectors φ(T) and ψ(T), we could determine. Therefore, we can only get their estimates, given as follows

(39) 1α10N¯1Q¯ν1αn0N¯nQ¯;

(40) αi0N¯iα10N¯1Q¯φi(T)αi0N¯iαn0N¯nQ¯,i=1,n¯;

(41) N¯iαn0N¯nQ¯ψi(T)N¯iα10N¯1Q¯,i=1,n¯.

Theorem 1

Let phase variables qi(t) be in descending order αi0N¯i , then:

1. we determine the optimal trajectory by the formulas:

(42) q˜i(t)=qi0exp[αi00tN˜i(θ)dθ],i=1,n¯;

2. the following rule uniquely determines the vector of optimal controls Ñ(t). For each t there exists an integer k ∈{1, . . . , n} such that the following double inequality holds i=1k1q˜i(t)N¯i<Q¯i=1kq˜i(t)N¯i , and

(43) N˜i(t)=Ni¯fori=1,2,,k1,

(44) N˜i(t)=0fori=k+1,,n,

(45) N˜k(t)=[Q¯i=1k1q˜i(t)N¯i]/q˜k(t);

3. we combine a set of values t related to one integer k into half-intervals (τk−1, τk]. Under the condition q10N¯1>Q¯ , the number of such half-intervals in the whole considered period is n.

Proof

Further, we will adhere to the condition formulated in paragraph (3) of Theorem 1. We install

(46) φ1(0)>φ2(0)>>φn(0).

Below, in Proposition 7, we show the same order of functions φi(t); that is, for any t ∈ [0, T] the strict inequalities hold: φ1(t) > φ2(t) > . . . > φn(t).

Maximizing the Hamiltonian (32) reduces to function maximization

(47) maxNG(q˜(t))i=1nφi(t)q˜i(t)Ni.

Now we introduce the following notation: u=q˜N and u¯=q˜N¯ . As a result, the maximization of the Hamiltonian (32) reduces to the following linear programming problem:

(48) i=1nφiuimax

under the conditions

(49) i=1nui=Q¯,

(50) 0uu¯.

At the same time, the coefficients φi decrease.

Taking (49) into account, we transform the linear function (48) and obtain the following linear programming problem at maximum

(51) i=1n1(φiφn)uimax.

Note that the coefficients φi=φiφn in the description of the linear function (51) are positive and decrease as we increase the sequence number components of the vector u.

The following procedure allows us to find such an optimal vector ũ, which gives the maximum of the linear function (51) and, accordingly, the maximum of the linear function (48).

At the first step, we choose u¯1 . If u¯1Q¯ , then u˜1=Q¯ and ũi = 0 when i=2,n¯ . In this case, we determine all the optimal values of the vector ũ. If u¯1<Q¯ , then u˜1=u¯1 and we replace Q¯=Q¯u¯1 .

Then we go to the second component u2, and repeat the first step's procedure, but with the changed value of Q¯ . At the second step, we change the adjusted magnitude of Q¯ again.

If we stop at the k-th step, then u˜i=u¯i when i=1,k1¯ , u˜k=Q¯ , ũi = 0 when i=k+1,n¯ . If we go through all the steps from one to n − 1, then u˜i=u¯ when i=1,n1¯ , u˜n=Q¯ .

From the above proof, it follows that for any k-th component of the control vector N(t), the following statements hold

1. if u˜k=u¯k , then N˜k(t)=N¯k ;

2. if ũk = 0, then Ñk(t) = 0;

3. if u˜k(0,u¯) , then N˜k(t)(0,N¯) .

Let's consider the optimal trajectory's dynamic behavior on the segment [0, T]. Suppose that at the initial time q10N¯1>Q¯ , then N˜1(0)(0,N¯1) , Ñi(0) = 0 for i=2,n¯ . In connection with a decrease of the phase variable q˜1(t) , there is a moment τ1 such that:

1. the following relations hold N˜1(t)(0,N¯1) and Ñi(t) = 0, i=2,n¯ on the interval t ∈ (0, τ1);

2. N˜1(τ1)=N¯1 and Ñi(τ1) = 0 for i=2,n¯ .

From the system of “adjoint” equations (33) we take into account (34), (30) and (31) at t ∈ (0, τ1) to get γ1(t) = 0; δ1(t) = 0; β(t) = φ1(t) = β1 = const; φi(t) = const for i=2,n¯ . In this case, the order relations established among the components of the vector function φ(t) at the initial moment t = 0 do not change on the segment [0, τ1].

Then there is a moment τ2 such that:

1. the following relations hold N˜1(t)=N¯1 , N˜2(t)(0,N¯2) and Ñi(t) = 0, i=3,n¯ on the interval (τ1, τ2);

2. N˜1(τ2)=N¯1 , N˜2(τ2)=N¯2 and Ñi(τ2) = 0 for i=3,n¯ .

From the system of “adjoint” equations (33), and taking into account (34), (30) and (31) at t ∈ (τ1, τ2), then we get γ2(t) = δ2(t) = 0; β(t) = φ2(t) = β2 = const; φi(t) = const for i=3,n¯ ; the function φ1(t) increases on the segment [τ1, τ2]. The latter follows from the inequality β1 = φ1(τ1) = φ1(0) > β2 = φ2(τ1) = φ2(0).

On the segment [0, τ1], the function φ1(t) − φ2(t) is constant, and on [τ1, τ2] this function strictly increases. The order relations established among the components of the vector functions φ(t) at the time t = 0, do not change on the segment [0, τ2]. Thus, we divide the segment [0, T] into n parts [0, τ1], [τ1, τ2], ..., [τn−1, τn]. Here τn = T. At the end of each kth segment the equality holds

i=1kqk(τk)N¯k=Q¯,k=1,n¯.

Proposition 7

Let t ∈ (0, τ1), k = 1, then:

β(t)=φ1(t)=φ1(0)=β1;N˜1(t)=Q¯q˜1(t)andN˜i(t)=0fori=2,n¯;φi(t)=φi(0)=βifori=2,n¯.

All functions φi(t) are constant on the segment [0, τ1]. On it, the strict inequalities hold φ1(t) > φ2(t) > . . . > φn(t).

Let t ∈ (τ1, τ2), k = 2, then:

β(t)=φ2(t)=φ2(0)=β2<β1;N˜1(t)=N¯1,N˜2(t)=Q¯q˜1(t)N¯1q˜2(t)andN˜i(t)=0fori=3,n¯;φ1(t)=φ2(τ1)+[φ1(τ1)φ2(τ1)]exp[α10N¯1(tτ1)]=φ2(0)+[φ1(0)φ2(0)]exp[α10N¯1(tτ1)],φi(t)=φi(0)=βifori=3,n¯.

The function φ1(t) strictly increases on the segment [τ1, τ2]. The remaining functions φi(t) for i=2,n¯ , are constant. On it, the following strict inequalities hold φ1(t) > φ2(t) > . . . > φn(t).

Let t ∈ (τk−1, τk), k=3,n¯ , then:

β(t)=φk(t)=φk(0)=βk<βk1;N˜i(t)=N¯ifori=1,k1¯;N˜k(t)=Q¯i=1i=k1q˜i(t)N¯iq˜k(t)andN˜i(t)=0fori=k+1,n¯;φi(t)=φk(τk1)+[φi(τk1)φk(τk1)]exp[αi0N¯i(tτk1)]=φk(0)+[φi(τk1)φk(0)]exp[αi0N¯i(tτk1)]fori=1,k1¯;φi(t)=φi(0)=βifori=k+1,n¯.

The functions φi(t), where i=1,k1¯ strictly increases on the segment [τk−1, τk], and the remaining functions φi(t), where i=k,n¯ , is constant. On it, the following strict inequalities hold φ1(t) > φ2(t) > . . . > φn(t).

Proof

Using the relations (28)–(31), (33), (42) and taking into account the continuity of the vector function φ(t), we write out formulae for calculating the variables β(t), N(t), φ(t) on each of the n segments.

Taking into account the earlier established orders φ1(0) > φ2(0) > . . . > φn(0), we prove the strict inequalities φ1(t) > φ2(t) > . . . > φn(t) on each time segment.

For the dimension n = 3 in Figure 1, we schematically show the behaviour dynamics of optimal controls Ñ1(t), Ñ2(t), Ñ3(t), the functions φ1(t), φ2(t), φ3(t) and the Lagrange multiplier β(t) obtained by solving Problem 3a.

Figure 1

Behavior dynamics of optimal controls Ñ1(t), Ñ2(t), Ñ3(t), the functions φ1(t), φ2(t), φ3(t) and the Lagrange multiplier β(t) by solving Problem 3a of dimension n = 3.

Let us describe a scheme for determining the continuous vector function φ(t). First, we find τ1, τ2, . . . , τn = T and qi(T) for i=1,n¯ .

Next, from the formula in (38), we calculate the exact values of the vector φ(T), whose components, according to Proposition 6, are in descending order.

The results in Proposition 7 are the basis for further description of the algorithm. We consider the segments sequentially starting with [τn−1, T] and ending with [0, τ1]. We know the Lagrange multiplier β(t) = βn on the n-th interval (τn−1, T), the function φn(t) = φn(0) = φn(T) = βn on the segment [0, T] and the values φi(T), i=1,n¯ , under which the following strict inequalities hold φ1(T) > φ2(T) > . . . > φn(T). The following formula calculates all other components of the vector function φ(t) on the nth segment [τn−1, T]:

(52) φi(t)=φn(T)+[φi(T)φn(T)]exp[αi0N¯i(Tt)],

where i=1,n1¯ .

The following inequalities hold φi(t) > φn(t) for i=1,n1¯ on the n-th segment of [τn−1, T]. Thus, we find φi(τn−1) when i=1,n¯ .

We proceed to the sequential study of functions φk(t) on the k-th segment of [τk−1, τk] for k=1,n1¯ . We know all the values of φi(τk) when i=1,n¯ . We also know the functions φi(t) = φi(τk) = βi, i=k+1,n¯ on the segment [τk−1, τk].We determine the Lagrange multiplier β(t) = βk on the k-th interval (τk−1, τk) and the function φk(t) = φk(0) = φk(τk) = βk on the segment [0, τk]. The following formula calculates all other components of the vector function φ(t) on the kth section [τk−1, τk]:

(53) φi(t)=φk(τk)+[φi(τk)φk(τk)]exp[αi0N¯i(τkt)]

when i=1,k1¯ .

Thus, we find all φi(τk−1), where i=1,n¯ for each magnitude of k=1,n1¯ . Including for the first segment, we determine all values of φi(0) when i=1,n¯ .

Therefore, on the segment [0, T], we have constructed a single-valued vector function φ(t) satisfying the transversality conditions (38) and the maximum principle in the Arrow form. Moreover, this function is unique. Indeed, any other order of the components of the vector φ(0) does not allow us to hold the necessary optimality conditions.

Taking into account the existence of an optimal trajectory for the problem under study, we come to the following conclusion. We have found the only “adjoint” vector function φ(t) and the corresponding optimal trajectory. We have proved Theorem 1.

Next, we consider the case when all indicators αi0N¯i are equal to each other and calculate the time T.

Proposition 8

For any admissible trajectories from sets (7)–(10) and (13) of the total shelf lengths under the condition

(54) α0N¯=α10N¯1==αn0N¯n

the lengths of their shelves are the same and they are calculated using the following formula

(55) T=i=1nVi0Q¯1α0N¯.

Proof

We determine the accumulated production by the formula:

(56) Q¯T=i=1n[Vi0Vi(T)].

Taking into account (5) and (54), we rewrite (56) as follows

(57) Q¯T=i=1n[Vi0qi(T)N¯iαi0N¯i]=i=1nVi0i=1nqi(T)N¯iα0N¯.

According to Propositions 3, 4, and 5, at the end of any admissible trajectory, the equality

(58) i=1nqi(T)N¯i=Q¯.

holds. From (57). taking into account (58), we get (55).

We now turn to the description of the algorithm for determining the minimum and maximum length of the total shelf.

Algorithm for calculating the minimum and maximum length of the total shelf for a group of gas fields

To minimize, we sort deposits in descending order of values αi0N¯i . To maximize, sort these values in reverse order.

First step

We put the first field into development and, using Proposition 8 with n = 1 we calculate the temporary step length:

(59) τ1=V10Q¯1α10N¯1.

Denote q10α10Q¯τ1 by q11 . Using (4) and (5), we get the following dynamics of the production wells flow rate of the first field

(60) q1(t)={q10α10Q¯tfort[0,τ1];q11exp(α10N¯1(tτ1))fortτ1.

Note that q11N¯1=Q¯ .

Second step

Denote by θ2 the temporary second step length. Denote by q12 the production wells flow rate of the first field at the end of the second step, and denote by q22 the production wells flow rate of the second field at the end of the second step.

We determine the temporary length of the second step from a solution of equation (61) with one unknown. Taking into account (4) and (5), τ2 = τ1 + θ2 and q11N¯1=q12N¯1+q22N¯2=Q¯ , we write out this equation as follows

(61) V20=Q¯[θ2+(1exp(α10N¯1θ2))(1α20N¯21α10N¯1)].

We consider equality (61) as an equation with one unknown variable θ2 We denote the right-hand side of equality (62) by F(θ2). We differentiate the function F(θ2) and analyze it. As a result, we obtain the following inequality

(62) F(θ2)=Q¯[1+exp(α10N¯1θ2)(α10N¯1α20N¯21)]>0forθ20.

Hence, the function F(θ2) is an increasing function for all positive values θ2.

The inequalities V20>F(0) and V20<F((V10+V20)/Q¯) hold. Equation (61) has no explicit solution. We can solve it numerically by dividing the segment in half. At the end of the second stage, we describe the producing wells flow rate dynamics of the second field

(63) q2(t)={q20fort[0,τ1];q20α20Q¯[tτ11exp(α10N¯1(tτ1)α10N¯1]fort[τ1,τ2];q22exp(α10N¯1(tτ2))fortτ2.

k-th step, k=3,n¯

Denote by θk the temporary k-th step length. Denote by qik when i=1,k¯ the t production wells flow rate of the i-th field at the end of the k-th step.

We determine the temporary length of the k-th step from the solution of the equation with one unknown. Taking into account (4) and (5), τk = τk−1 + θk and i=1k1qik1N¯i=i=1kqikN¯i=Q¯ , we write out this equation as follows

Vk0=Q¯θk+i=1k1qik1N¯i(1exp(αi0N¯iθk))(1αk0N¯k1αi0N¯i).

It has one solution. To prove this, we denote the right-hand side of the last equality by F(θk). It is easy to show the following properties of the function F(θk): F′(θk) > 0 for all θk > 0; Vk0>F(0) ; Vk0<F(i=1kVi0/Q¯) .

The last equation has no explicit solution. We can solve it numerically by dividing the segment in half. At the end of the k-th stage, we describe the producing wells flow rate dynamics of the k-th field

qk(t)={qk0fort[0,τk1];qk0αk0i=1k1qik1N¯i[tτk11exp(αi0N¯i(tτk1)αi0N¯i]fort[τk1,τk];qkkexp(αk0N¯k(tτk))fortτk.

This completes the description of the algorithm.

5 Conclusion

When developing a group of gas fields, the main problem at any given time is the fulfilment of the current gas production plan. In the initial period of development of the gas fields group, we can achieve this effect in different ways, including by turning on or off wells. However, there comes a time when there is not enough capacity to fulfil the current gas production plan. All available wells in all fields are involved in the development. At this point, the process of developing the gas fields group is interrupted. Each development dynamics of the gas fields group has its shelf length. It is required to determine its exact estimates.

We solve the problems of optimal control for the minimum and maximum shelf lengths of gas fields group under strict restrictions on total gas production. To solve the set tasks, we sort all fields by the value of the parameter αi0N¯i . When optimizing for a minimum, we vary them from the maximum magnitude of the indicator to the minimum one. Otherwise, we arrange parameters in the reverse order.

We divide the entire optimal period into n parts equal to the number of deposits in the group. At the beginning of each part of the optimal period, we bring into development the next field in order. All previously commissioned fields continue to operate at their maximum capacities until the end of the optimal period and cannot ensure the implementation of the entire current production plan. Therefore, we include a new deposit in development.

At the end of the current part of the optimal period, all wells of the newly developing field are in development and are operating at full capacity.

At the end of the optimal period, all wells in the group of fields operate at full capacity with maximum efficiency. Any further development of the group only leads to non-fulfilment of the current plan for the volume of gas production. This fact interrupts the dynamic process. If we set all parameters αi0N¯i equal, then any acceptable field development policy is optimal.

We have presented in this article an algorithm for calculating the minimum length of the total shelf for a group of gas fields.

1. Conflict of interest: Authors state no conflict of interest.

References

[1] Vyakhirev R.I., Korotaev Y.P., Kabanov N.I., Theory and experience of gas recovery, Nedra, Moscow, 1998 (in Russian)Search in Google Scholar

[2] OPEC. OPEC Annual Statistical Bulletin, 2015, http://large.stanford.edu/courses/2016/ph241/almajid2/docs/ASB2015.pdfSearch in Google Scholar

[3] Arrow K.J., Applications of control theory to economic growth, Mathematics of the decision sciences, 1968, 2, 85–119Search in Google Scholar

[4] Pontryagin L.S., Boltyanskii V.G., Gamkrelidze R.V., Mishechenko E.F., The mathematical theory of optimal processes, Fizmatgiz, Moscow, 1961 (in Russian)Search in Google Scholar

[5] Skiba A.K., Optimal growth with a convex-concave production function, Econometrica, 1978, 3(46), 527–53910.2307/1914229Search in Google Scholar

[6] Skiba A.K., Maximization of the accumulated extraction in a gas fields model, In: Proceedings of the 9th International Conference on Optimization and Applications, (1–5 October 2018, Petrovac, Montenegro), CCIS, 2019, 974, 453–46910.1007/978-3-030-10934-9_32Search in Google Scholar

[7] Ter-Krikorov A.M., Optimal control and mathematical economics, Nauka, Moscow, 1977 (in Russian)Search in Google Scholar

[8] Aseev S.M., Kryazhimskii A.B., The Pontryagin maximum principle and optimal economic growth problems, Proc. Steklov Inst. Math., Nauka, Moscow, 2007, 257 (in Russian)10.1134/S0081543807020010Search in Google Scholar

[9] Balder E.J., An existence result for optimal economic growth problems, Journal of Mathematical Analysis and Applications, 1983, 95, 195–21310.1016/0022-247X(83)90143-9Search in Google Scholar

[10] Margulov R.D., Khachaturov V.R, Fedoseev A.V.m System analysis in long-term planning of gas production, Nedra, Moscow, 1992 (in Russian)Search in Google Scholar

[11] Khachaturov V.R., Solomatin A.N., Zlotov, A.V., Planning and design of development of oil and gas producing regions and deposits: mathematical models, methods, application,. LENAND, Moscow, 2015 (in Russian)Search in Google Scholar

[12] Moiseev N.N., Elements of the theory of optimal systems, Nauka, Moscow, 1975 (in Russian)Search in Google Scholar

[13] Lee E. B., and Markus L., Foundations of optimal control theory, Wiley, New York, 1961Search in Google Scholar