CN104615862A - High water-cut oilfield well position determining method based on evolutionary algorithm - Google Patents

Abstract

The invention provides a high water-cut oilfield well position determining method based on an evolutionary algorithm. The method comprises the steps of building an oil reservoir model of a predetermined oil area, and dividing the oil reservoir model into grids; obtaining the oil reservoir data, and obtaining a surplus oil potential enrichment area in a history matching mode according to the oil reservoir data; computing a well position constraint value for representing the surplus oil distribution and injected water waves and volume in the surplus oil potential enrichment area to obtain the well position advantage area; carrying out iteration on an oil well in the advantage area through the evolutionary algorithm to compute a fitness function, and evaluating the fitness function through a Kriging algorithm, wherein the evolutionary algorithm is adjusted according to the fitness function, and the Kriging algorithm considers the multi-dimensional space, the nugget effect and the local weight; obtaining a maximum fitness function value according to the iteration computing of the evolutionary algorithm so as to determine the well position. The improved evolutionary algorithm and the Kringing algorithm serve as the optimization tools, a dynamic potential index is built to reduce the searching space of the evolutionary algorithm, and the optimization computing convergence precision and speed are greatly improved.

Description

Based on the method for the high water cut oil field determination well location of evolution algorithm

Technical field

The invention belongs to high water cut oil field modification scenario well location technical field, particularly relate to a kind of method of the high water cut oil field determination well location based on evolution algorithm.

Background technology

For the oil reservoir of China's waterflooding, major part is the waterflooding High water cut maturing field of decades, the complicacy of subsurface reservoir heterogeneous body and oil water relation, remaining oil distribution is caused to present the complex characteristic of " totally scattered, concentration of local ", great difficulty is caused to the adjustment of oil field later stage, directly affects waterflooding development effect and economic benefit.Determining that Reasonable adjustment well location is the significant subject technology of High water cut maturing field adjustment, is the important means improving reserves exploitation degree and recovery ratio.

Reasonable adjustment well location is by the impact of complex nonlinear relation between the uncertain factor such as geology reservoir, fluid properties and factor, and conventional Method for Numerical needs a large amount of modeling schemes, and is difficult to obtain optimum adjustment well location.Intelligent optimization algorithm is the effective means determining optimal correction well location.

The optimization method of adjustment well location conventional at present is mainly divided three classes: random optimization method, gradient optimization method, Quality Map method.Found all there is certain inadaptability for these three class methods of High water cut maturing field by investigation and analysis:

(1) random optimization method.Such algorithm is from the search of many group initial solutions, easy acquisition globally optimal solution, but many geology of the complicacy of high water-cut stage underground oil water relation and existence is uncertain, make such algorithm search space larger, computing time is long, in addition the artificial property of Optimal Parameters selection, there is larger room for promotion in the performance of such algorithm.

(2) gradient optimization method.Such algorithm, by the gradient of calculating target function, makes well location move to greatest gradient direction, thus realizes well location Optimizing Search.But such algorithm generally needs the source code of simulator, a large amount of iteration is needed to restrain, and very responsive to choosing of initial solution.Therefore, such algorithm is difficult to carry out in the maturing field developed on a large scale.

(3) Quality Map method.These class methods have two kinds of purposes: one is intuitively determine optimum adjustment well location by the high level region in qualitative map, avoids loaded down with trivial details numerical simulation; Two is the constraint condition being used as evolution algorithm, instructs evolution algorithm well spacing in reasonable region, to reduce search volume, to improve speed of convergence.Because this class methods counting yield is higher, and be build based on geologic model and Production development history because of it, computational accuracy is also higher, is therefore widely used in the well location optimization problem under high water-cut stage complex geological condition.Although current Quality Map method considers some Static and dynamic attribute of oil reservoir, do not consider reservoir drive mechanism, existing well pattern earial drainage and inject the impact of water swept volume on adjustment well location, certain error is caused to the well location optimized.

Summary of the invention

In view of Problems existing in current adjust well bit optimization, the invention provides a kind of method of the high water cut oil field determination well location based on evolution algorithm, its technical scheme is as follows:

Based on the method for the high water cut oil field determination well location of evolution algorithm, comprise the following steps:

S1: the reservoir model setting up predetermined oil district, by described reservoir model grid division;

S2: the oil reservoir data obtaining described predetermined oil district, draws potentiality of remaining oil enrichment region according to described oil reservoir data according to history matching;

S3: the well location binding occurrence calculating each grid in described potentiality of remaining oil enrichment region, described well location binding occurrence represents remaining oil distribution and injects water swept volume, draws the grid that maximum well location binding occurrence is corresponding, described grid representation well location dominant area;

S4: use evolution algorithm to carry out its fitness function of iterative computation the oil well in described dominant area, and utilize Kriging method to assess fitness function, described evolution algorithm carries out self-adaptative adjustment according to fitness function, and described Kriging method considers hyperspace, nugget effect and partial weight;

S5: according to evolution algorithm iterative computation, obtains maximum adaptation degree functional value, thus determines well location.

As above based on the method for the high water cut oil field determination well location of evolution algorithm, in S3, determine that the step of well location dominant area is:

S31: determine the particle travel-time, the described particle travel-time utilizes following formula to calculate:

In formula: the travel-time that TOF (x, y, z) needs for particle in-position (x, y, z), ψ, χ are called double-current function, represent gradient; for reservoir porosity, s is the distance along streamline;

S32: on the basis in particle travel-time, calculate movable oil bulk index MOVI, oil phase fluid ability index CFI and streamline dynamic measurement DM, above-mentioned three indexes are calculated by following formula:

Movable oil bulk index MOVI:

MOVI＝(S _ro-S _orw)×PV×NTG

Oil phase fluid ability index CFI:

CFI＝K×K _ro×DZ×NTG

Streamline dynamic measurement DM:

DM＝(TOF _p+TOF _i)

In formula: S _ro: remaining oil saturation distribution when prediction starts, %;

S _orw: residual oil saturation distribution at the end of prediction, %;

PV: volume of voids, m ³;

NTG: net-gross ratio, represents net-to-gross ratio, dimensionless;

K: absolute permeability distributes, 10 ^-3μm ²;

Kro: oil relative permeability, refers in water-oil phase system, the ratio of permeability to oil and absolute permeability, dimensionless;

DZ: the size of longitudinal grid, m;

Above-mentioned parameter all obtains by obtaining oil reservoir data;

TOF _ifor taking water injection well as the particle travel-time that starting point calculates, TOF _pfor taking producing well as the particle travel-time that starting point calculates, DM is total travel-time;

S33: movable oil bulk index MOVI, oil phase fluid ability index CFI that S12 is obtained and streamline dynamic measurement DM Regularization, the concrete steps of Regularization are as follows:

1. by the 10% and 90% maximal value X determining attribute on model attributes cumulative distribution frequency plot _maxwith minimum value X _min;

2. according to the following formula Regularization is carried out to each attribute:

$X=\left\{\begin{array}{cc}0,& X<X_{\min }\\ \frac{X-X_{\max }}{X_{\max }-X_{\min }},& X_{\min }\&le;X<X_{\max }\\ 1,& X\&GreaterEqual;X_{\max }\end{array}\right.$

S34: based on three indexes of above-mentioned Regularization, utilizes following formula to calculate dynamic potential index DOI:

$\mathrm{DOI}=\sqrt[3]{\mathrm{MVOI}\&times;\mathrm{CFI}\&times;\mathrm{DM}}$

Dynamic potential index DOI is the concentrated expression of remaining oil distribution situation and well pattern sweep conditions, DOI more close to 1 region, show that these region remaining oil reserves are larger, crude oil fluidity better, and not yet by floood conformance, thus these regions become the dominant area boring adjust well.

As above based on the method for the high water cut oil field determination well location of evolution algorithm, in S4, consider that the step of the Kriging method of nugget effect and partial weight is:

S41: consider hyperspace, utilizes following formula to calculate empirical covariance:

$\mathrm{Cov}{\left(h\right)}_1=\frac{\mathrm{\&Sigma;}\left[z\left(x_i\right)z(x_i+h)\right]}{N\left(h\right)-1}-\frac{\mathrm{\&Sigma;z}\left(x_i\right)\mathrm{\&Sigma;z}(x_i+h)}{N\left(h\right)[N\left(h\right)-1]}$

Wherein: Cov (h) ₁the empirical covariance function value of sample point within expression distance h, z (x _i) be the adaptive response function related functions value with sample point, N (h) represents the total sample point number within distance h;

S42: consider nugget effect and partial weight, utilize following formula to calculate theoretical covariance:

$\mathrm{Cov}{\left(h\right)}_2=\mathrm{\&rho;}\left(h\right){\widehat{\mathrm{\&sigma;}}}^2$

$\mathrm{\&rho;}\left(h\right)=2(1-\mathrm{nu}){\mathrm{\&rho;}}_{\mathrm{\&theta;}}\left(h\right){\left(\frac{\mathrm{bh}}{2}\right)}^{s-1}{K}_{s-1}\left(\mathrm{bh}\right)\mathrm{\&Gamma;}{(s-1)}^{-1}$

Wherein: Cov (h) ₂for theoretical covariance, ρ (h) is the spatial correlation function at a distance of h any two points, for Kriging estimation variance, nu is block gold number, and K is second-order modified Bessel function, and Г is gamma function, and parameter b, s are the real number characterizing relevance function shape, can obtain parameter b, s by solving following unconstrained optimization problem:

$\underset{b,s}{\min }\mathrm{\&Sigma;}{[2{\widehat{\mathrm{\&sigma;}}}^2(1-\mathrm{nu}){\mathrm{\&rho;}}_{\mathrm{\&theta;}}\left(h\right){\left(\frac{\mathrm{bh}}{2}\right)}^{s-1}{K}_{s-1}\left(\mathrm{bh}\right)\mathrm{\&Gamma;}{(s-1)}^{-1}-\mathrm{Cov}\left(h\right)]}^2$

ρ _θh () is partial weight function, utilize following formula to calculate:

${\mathrm{\&rho;}}_{\mathrm{\&theta;}}\left(h\right)=\left\{\begin{array}{cc}(1-\frac{h}{\mathrm{\&theta;}})\cos \left(\frac{\mathrm{\&pi;h}}{\mathrm{\&theta;}}\right)+\frac{\sin \left(\frac{\mathrm{\&pi;h}}{\mathrm{\&theta;}}\right)}{\mathrm{\&pi;}},& h<\mathrm{\&theta;}\\ 0,& h\&GreaterEqual;\mathrm{\&theta;}\end{array}\right.$

Wherein: θ is supporting degree;

S43: according to the empirical covariance calculated and theoretical covariance, utilize least square fitting covariance;

S44: according to the linear relationship of covariance and variogram, draw variogram;

S45: adopt complete pivot elimination method to solve system of linear equations, solve and obtain weight coefficient, complete the foundation of Kriging method.

As above based on the method for the high water cut oil field determination well location of evolution algorithm, the evolution algorithm described in S5 comprises genetic algorithm, and crossover probability and the mutation probability of described genetic algorithm automatically adjust according to auto-adaptive function:

Utilize following formula that crossover probability is regulated automatically according to fitness value change:

$p_x=\left\{\begin{array}{cc}\frac{a(F_{\max }-F_x)}{F_{\max }-F_{\mathrm{avg}}},& F_x\&GreaterEqual;F_{\mathrm{avg}}\\ b,& F_x<F_{\mathrm{avg}}\end{array}\right.$

In formula: F _maxfor maximum adaptation value in population, F _avgfor kind of a group mean adaptive value, F _xfor fitness value larger in two individualities will be intersected, a=0.8, b=1;

Utilize following formula that mutation probability is regulated automatically according to fitness value change:

$p_m=\left\{\begin{array}{cc}\frac{c(F_{\max }-F_m)}{F_{\max }-F_{\mathrm{avg}}},& F_m\&GreaterEqual;F_{\mathrm{avg}}\\ d,& F_m<F_{\mathrm{avg}}\end{array}\right.$

In formula: F _maxfor maximum adaptation value in population, F _avgfor kind of a group mean adaptive value, F _mfor the fitness value of individuality that will make a variation, c=0.01, d=0.05.

As above based on the method for the high water cut oil field determination well location of evolution algorithm, the evolution algorithm described in S5 comprises particle cluster algorithm, and the particle displacement of described particle cluster algorithm upgrades and flying speed automatically adjusts according to auto-adaptive function:

$x_{i,j}^{k+1}=x_{i,j}^k+v_{i,j}^{k+1}\mathrm{\&Delta;t}$

In formula: be respectively the position of i-th particle when evolving to kth+1 time step and a jth element of speed; for a jth element of the optimum solution that i-th particle itself when evolving to kth time step finds; for a jth element of the globally optimal solution that whole population when evolving to kth time step is found at present; c ₁, c ₂for positive Studying factors; r ₁, r ₂for equally distributed random number between [0,1]; Δ t is time increment, for compressibility coefficient:

ω is inertia weight:

$w=\left\{\begin{array}{cc}{w}_{\min }-\frac{(w_{\max }-w_{\min })(F-F_{\min })}{(F_{\mathrm{avg}}-F_{\min })},& F\&le;F_{\mathrm{avg}}\\ w_{\max }& F>F_{\mathrm{avg}}\end{array}\right.$

In formula: ω _min, ω _maxbe respectively maximal value and the minimum value of inertia weight, F is the current adaptive value of particle; F _avg, F _minfor average criterion adaptive value and the minimum adaptive value of current all particles.

By above technical scheme, beneficial effect of the present invention is: the method for the high water cut oil field determination well location based on evolution algorithm of the present invention is to take into account the evolution algorithm of precision and solving speed and Kriging method is main optimization tool, take cumulative oil production as objective function, and construct the search volume that the dynamic potential index considering the interference of old well, fault boundary, invalid reservoir, well spacing restriction and embody remaining oil distribution and injection water swept volume reduces evolution algorithm, thus make optimization computation process quick, save time, convergence precision and speed increase greatly.

Accompanying drawing explanation

In order to be illustrated more clearly in the technical scheme in the embodiment of the present invention, below the accompanying drawing used required in describing embodiment is briefly described, apparently, accompanying drawing in the following describes is only some embodiments of the present invention, for those of ordinary skill in the art, under the prerequisite not paying creative work, other accompanying drawing can also be obtained according to these accompanying drawings.

Fig. 1 is the process flow diagram of the method for the high water cut oil field determination well location based on evolution algorithm of embodiment of the present invention;

Fig. 2 a is the remaining oil panel map in the oily district of well location to be adjusted;

Fig. 2 b is the remaining oil Soil profile figure in the oily district of well location to be adjusted;

Fig. 2 c is the potentiality of remaining oil enrichment region in the oil district one of the main force of well location to be adjusted;

Fig. 2 d is the potentiality of remaining oil enrichment region in the oil district two of the main force of well location to be adjusted;

Fig. 3 a is particle travel-time schematic diagram;

Fig. 3 b is particle travel-time plane distribution;

Fig. 4 a is the model attributes cumulative distribution frequency plot of movable oil bulk index MVOI;

Fig. 4 b is the model attributes cumulative distribution frequency plot of oil phase fluid ability index CFI;

Fig. 4 c is the model attributes cumulative distribution frequency plot of streamline dynamic measurement DM;

Fig. 5 a is the panel map of the movable oil bulk index MVOI after Regularization;

Fig. 5 b is the panel map of the oil phase fluid ability index CFI after Regularization;

Fig. 5 c is the panel map of the streamline dynamic measurement DM after Regularization;

Fig. 6 a is the dynamic potential index DOI distribution plan after Regularization;

Fig. 6 b is that after Regularization, after Regularization, dynamic potential index DOI longitudinally superposes flat distribution map;

Fig. 7 a is the well location dominant area that dynamic potential index DOI>0.2 determines;

Fig. 7 b is the well location dominant area that dynamic potential index DOI>0.4 determines;

Fig. 7 c is the well location dominant area that dynamic potential index DOI>0.6 determines;

Fig. 7 d is the well location dominant area that dynamic potential index DOI>0.8 determines;

Fig. 8 a is the cumulative oil production evolution curve that existing GA algorithm obtains;

Fig. 8 b is the optimal well location coordinate that existing GA algorithm obtains;

Fig. 9 a is the cumulative oil production evolution curve that existing PSO algorithm obtains;

Fig. 9 b is the optimal well location coordinate that existing PSO algorithm obtains;

Figure 10 a is the cumulative oil production evolution curve that the HGA algorithm of embodiment of the present invention obtains;

Figure 10 b is the optimal well location coordinate that the HGA algorithm of embodiment of the present invention obtains;

Figure 11 a is the cumulative oil production evolution curve that the HPSO algorithm of embodiment of the present invention obtains;

Figure 11 b is the optimal well location coordinate that the HPSO algorithm of embodiment of the present invention obtains;

Figure 12 is different well location optimized algorithm effect of optimization contrast under DOI>0.6 constraint condition.

Embodiment

Below in conjunction with the accompanying drawing in the embodiment of the present invention, be clearly and completely described the technical scheme in the embodiment of the present invention, obviously, described embodiment is only the present invention's part embodiment, instead of whole embodiments.Based on the embodiment in the present invention, those of ordinary skill in the art, not making the every other embodiment obtained under creative work prerequisite, belong to the scope of protection of the invention.

Refer to Fig. 1.Based on the method for the high water cut oil field determination well location of evolution algorithm, comprising:

S1: the reservoir model setting up predetermined oil district, by described reservoir model grid division;

S2: the oil reservoir data obtaining described predetermined oil district, draws potentiality of remaining oil enrichment region according to described oil reservoir data according to history matching;

S3: the well location binding occurrence calculating each grid in described potentiality of remaining oil enrichment region, described well location binding occurrence represents remaining oil distribution and injects water swept volume, draws the grid that maximum well location binding occurrence is corresponding, described grid representation well location dominant area;

S4: use evolution algorithm to carry out its fitness function of iterative computation the oil well in described dominant area, and utilize Kriging method to assess fitness function, described evolution algorithm carries out self-adaptative adjustment according to fitness function, and described Kriging method considers hyperspace, nugget effect and partial weight;

S5: according to evolution algorithm iterative computation, obtains maximum adaptation degree functional value, thus determines well location.

By above technical scheme, the beneficial effect of embodiment of the present invention is: the method for the high water cut oil field determination well location based on evolution algorithm of the present invention is to take into account the evolution algorithm of precision and solving speed and Kriging method is main optimization tool, take cumulative oil production as objective function, and construct consider old well interference, fault boundary, invalid reservoir, the dynamic potential index of well spacing restriction and embodiment remaining oil distribution and injection water swept volume reduces the search volume of evolution algorithm, thus make optimization computation process quick, save time, convergence precision and speed increase greatly.

Be applied as example with the present invention at certain high water cut oil field block below, further illustrate the method for the high water cut oil field determination well location based on evolution algorithm of the embodiment of the present invention.

At present this block of this oil field waterflooding development 20 years, underground oil water relation is comparatively complicated.After oil field injection and extraction well pattern has carried out infilling adjustment repeatedly, water drive control (only 55%) and development degree (20% ~ 35%) all still lower, inject water individual layer and advance by leaps and bounds phenomenon still relatively seriously, oil recovery rate is low, and oil offtake is successively decreased soon.Need, in conjunction with fine geological study result, to improve well pattern and injection-production relation further, improve the development effectiveness of late high water content period.Adopt multiple method to calculate, this block rational spacing between wells is 150-200m, and average well spacing is 220m at present, therefore still there are considerable well-location adjustment potentiality at present.

The present invention is directed to above-mentioned oil field block as follows based on the method step of evolution algorithm adjustment well location:

S1: the reservoir model setting up predetermined oil district, by described reservoir model grid division;

The inventive method is chosen middle part physical property main force's development block that is better, stratigraphic dip about 16 ° and is studied.Set up reservoir model to the block chosen, and carry out stress and strain model to it, grid is 77 × 49 × 15, and total nodes is 56595.

S2: the oil reservoir data obtaining described predetermined oil district, draws potentiality of remaining oil enrichment region according to described oil reservoir data according to history matching;

Adopt history matching auxiliary with manual modification parameter, obtain the overall distribution situation of this block remaining oil as shown in Fig. 2 a, 2b, and the remaining oil distribution situation of main layer position is as shown in Fig. 2 c, 2d.

After oilfield exploitation enters high water cut rate, the scattered reason of remaining oil distribution is owing to there is a large amount of minor fault, oil water relation is very complicated, add that scale of sand bodies is little, change fast feature, after tomography cutting, determine connecting degree in reservoir plane not high, injection-production pattern is irregular, and development late stage remaining oil distribution is more scattered.These optimizations being characterized as later stage adjustment well location bring difficulty.

S3: the well location binding occurrence calculating each grid in described potentiality of remaining oil enrichment region, described well location binding occurrence represents remaining oil distribution and injects water swept volume, draws the grid that maximum well location binding occurrence is corresponding, described grid representation well location dominant area;

According to above-mentioned block potentiality of remaining oil enrichment region distribution situation, calculate the dynamic potential index DOI in study area, determine well location dominant area according to dynamic potential index DOI.Concrete grammar and step as follows:

S31, determine the particle travel-time

Streamline simulation is the tool of waterflooding development management and late high water content period Fine Reservoir Numerical, its method three-dimensional simulation is reduced to a series of one dimensional flow line model, in reduction process, calculate fluid-flow data, there is the advantage of the larger quantitative series certificate of process.Owing to having the feature keeping obvious leading displacement edge and reduce grid bearing impact in displacement process, flowing is made to calculate not by the restriction of geologic grid, effectively prevent the impact of numerical dispersion that Finite Difference Meshes causes and grid orientation effect, streamline simulation is more accurate for portraying of the complicated oil water relation of High water cut maturing field.

An important parameter in streamline simulation is particle travel-time (Time Of Flight, TOF), representing that certain particle arrives the time in space required for certain point (function of coordinate x, y and z) from producing well or water injection well, is a physical quantity relevant with space.

Theoretical according to streamline simulation, the particle travel-time may be defined as:

In formula: the travel-time that TOF (x, y, z) needs for particle in-position (x, y, z), unit is sky, ψ, χ are called double-current function, represent gradient; for reservoir porosity; S is the distance along streamline, and namely particle moves the trajectory of formation.The 3D velocity field restore calculation one dimensional flow line model particle travel-time that usual employing streamline simulation or finite-difference modeling device calculate.

The particle travel-time can portray the flowing speed of subsurface deposit different fluid phase, and visual waterflood front.The present invention by asking arithmetic mean to streamline TOF values all in Finite Difference Meshes, and using the TOF value of average T OF value as this grid, can obtain TOF distribution as described in Fig. 3 a, 3b.TOF distribution is the good index of water injection well sweep conditions and producing well draining situation under reflection well pattern.

S32, calculating movable oil bulk index MOVI, oil phase fluid ability index CFI and streamline dynamic measurement DM

On the basis in particle travel-time, introduce the movable oil bulk index MOVI characterizing movable remaining oil reserves, the oil phase fluid ability index CFI characterizing oil phase fluid ability and characterize inject ripples and, the streamline dynamic measurement DM of producing well earial drainage situation.Above-mentioned three indexes are calculated by following formula:

(1) movable oil bulk index MOVI

MOVI＝(S _ro-S _orw)×PV×NTG (2)

(2) oil phase fluid ability index CFI

CFI＝K×K _ro×DZ×NTG (3)

(3) streamline dynamic measurement DM

DM＝(TOF _p+TOF _i) (4)

In formula: S _ro: remaining oil saturation distribution when prediction starts, %;

S _orw: residual oil saturation distribution at the end of prediction, %;

PV: volume of voids, m ³;

NTG: net-gross ratio, represents net-to-gross ratio, dimensionless;

K: absolute permeability distributes, 10 ^-3μm ²;

Kro: oil relative permeability, refers in water-oil phase system, the ratio of permeability to oil and absolute permeability, dimensionless;

DZ: the size of longitudinal grid, m.

Above-mentioned parameter all obtains by obtaining oil reservoir data.

Have movable and immovable two parts in remaining oil, MOVI refers to the distribution situation of movable oil part; And CFI characterizes the flox condition of movable oil part; TOF _ifor taking water injection well as the particle travel-time that starting point calculates, embody the swept volume injecting water; TOF _pfor taking producing well as the particle travel-time that starting point calculates, reflect the earial drainage situation of producing well; DM is total travel-time, and it embodies the combined influence of producing well earial drainage volume and water injection well swept volume.

S33, to above-mentioned three index Regularizations

Due to movable oil bulk index MOVI, oil phase fluid ability index, CFI is different with the dimension DM of streamline dynamic measurement, and numerical value difference may greatly to 1 ~ 2 order of magnitude.In order to make each index have comparability, need, according to the maximin of each index, each index is carried out Regularization.The concrete steps of Regularization are as follows:

1. by the 10% and 90% maximal value X determining attribute on model attributes cumulative distribution frequency plot _maxwith minimum value X _min;

2. then according to the following formula Regularization is carried out to each attribute:

$X=\left\{\begin{array}{cc}0,& X<X_{\min }\\ \frac{X_{\max }-X}{X_{\max }-X_{\min }},& X_{\min }\&le;X<X_{\max }\\ 1,& X\&GreaterEqual;X_{\max }\end{array}\right.---\left(5\right)$

In embodiment of the present invention, according to the remaining oil distribution after main force's fault block factor of porosity, permeability, reservoir thickness distribution and history matching, calculate movable oil bulk index MVOI, oil phase fluid ability index CFI and streamline dynamic measurement DM tri-indexes, shown in model attributes cumulative distribution frequency plot Fig. 4 a, 4b, 4c of drawing movable oil bulk index MVOI, oil phase fluid ability index CFI and streamline dynamic measurement DM.To add up on model attributes cumulative distribution frequency plot the maximal value X that 10% and 90% determines attribute _maxwith minimum value X _minas shown in table 1.

The property value that each index of table 1 10% and 90% accumulated probability is corresponding

After Regularization, three index panel maps are as shown in Fig. 5 a, 5b, 5c.

S34, calculate dynamic potential index DOI

Based on three indexes of above-mentioned Regularization, following formula is utilized to calculate the dynamic potential index DOI of reflection well location dominant area:

$\mathrm{DOI}=\sqrt[3]{\mathrm{MVOI}\&times;\mathrm{CFI}\&times;\mathrm{DM}}---\left(6\right)$

Dynamic potential index DOI is the concentrated expression of remaining oil distribution situation and well pattern sweep conditions, DOI more close to 1 region, show that these region remaining oil reserves are larger, crude oil fluidity better, and not yet by floood conformance, thus these regions become the dominant area boring adjust well.Have the following advantages using dynamic potential index DOI as the dominant area tool of well location optimization: high REGION OF WATER INJECTION OILFIELD accuracy of identification is higher; Consider that the new index injecting water swept volume is more reasonable to the identification can adopting remaining oil.

According to three indexes after Regularization, calculate dynamic potential index DOI distribution plan as shown in Figure 6 a, thickness weighting process is carried out to dynamic potential index DOI, obtain longitudinally superposing flat distribution map as shown in Figure 6 b.

Dynamic potential index distribution comprehensively can embody potentiality of remaining oil, fluid ability and well pattern sweep conditions, higher to the high REGION OF WATER INJECTION OILFIELD accuracy of identification of well location dominant area.

According to calculating formula (6), the region that dynamic potential index DOI is larger, is more suitable for boring adjust well.But be difficult to the quantitative relationship determining well location dominant area and dynamic potential index DOI in practical operation, therefore the present invention makes DOI>0.2,0.4,0.6,0.8 4 level dynamic potential index distribution plan respectively as the new constraint condition of adjust well bit optimization.According to dynamic potential index distribution plan, the well location dominant area of different stage can be determined, as shown in Fig. 7 a, 7b, 7c, 7d.

The dynamic potential index distribution plan of different stage all has impact to the search volume of well location optimized algorithm and optimal well location.The Quality Map that DOI value is less, search volume is comparatively large, and the uncertainty of well location optimized algorithm is higher, even likely makes algorithm convergence occur to local solution " precocity ".By the applying of above-mentioned constraint condition, greatly reduce the search volume of algorithm.

S4: use evolution algorithm to carry out its fitness function of iterative computation the oil well in described dominant area, and utilize Kriging method to assess fitness function, described evolution algorithm carries out self-adaptative adjustment according to fitness function, and described Kriging method considers hyperspace, nugget effect and partial weight;

In order to solve evolution algorithm call numerical simulation often, computing time long problem, the common Kriging method in geostatistics is extended to multidimensional by the present invention, and introduce nugget effect, partial weight makes improvements.The Kriging method improved replaces the process of the matching variogram in existing Kriging method with matching covariance.

In common gram in geostatistics, golden equation can use following matrix representation:

[K][λ]＝[D] (7)

Wherein:

$\left[K\right]=\left[\begin{array}{ccccc}{\mathrm{\&gamma;}}_{11}& {\mathrm{\&gamma;}}_{12}& ...& {\mathrm{\&gamma;}}_{1n}& 1\\ {\mathrm{\&gamma;}}_{21}& {\mathrm{\&gamma;}}_{22}& ...& {\mathrm{\&gamma;}}_{2n}& 1\\ .& .& & .& .\\ .& .& & .& .\\ .& .& & .& .\\ {\mathrm{\&gamma;}}_{n1}& {\mathrm{\&gamma;}}_{n2}& ...& {\mathrm{\&gamma;}}_{\mathrm{nn}}& 1\\ 1& 1& ...& 1& 0\end{array}\right],\left[\mathrm{\&lambda;}\right]=\left[\begin{array}{c}{\mathrm{\&lambda;}}_1\\ {\mathrm{\&lambda;}}_2\\ .\\ .\\ .\\ {\mathrm{\&lambda;}}_u\\ \mathrm{\&mu;}\end{array}\right],\left[D\right]=\left[\begin{array}{c}\mathrm{\&gamma;}(x_1,x)\\ \mathrm{\&gamma;}(x_2,x)\\ .\\ .\\ .\\ \mathrm{\&gamma;}(x_n,x)\\ 1\end{array}\right]---\left(8\right)$

Adopt complete pivot elimination method to solve system of linear equations, solve golden coefficient lambda=K in obtaining gram ^-1d and Lagrange coefficient μ, can calculate the estimated value of common Ke Lijin calculating according to weight coefficient λ.And then can in the hope of Kriging estimation variance for:

${\widehat{\mathrm{\&sigma;}}}^2\left(x\right)=D^T\mathrm{\&lambda;}---\left(9\right)$

Matrix of coefficients in common gram in golden equation implies variogram.Therefore, building Kriging method needs to choose rational theoretical variogram, by theoretical and the matching of experience variogram, determines the variogram in matrix of coefficients.

According to geology statistical theory, covariance function and variogram linear.

S41, consideration hyperspace, calculate empirical covariance.

In well location is optimized, build the sample of Ke Lijin for the coordinate of well location that will design and the cumulative oil production of correspondence.Therefore, sample point is the point of hyperspace, and its dimension depends on the well number of optimization, generally can be expressed as: the dimension=optimization well number × 2+1 of sample point.Such as, when optimizing 3 mouthfuls of wells, the sample point of each individuality is 7 × n dimension, such as formula:

$M=\left(\begin{array}{ccccccc}{X}_{11}& Y_{11}& X_{21}& Y_{21}& X_{31}& Y_{31}& Q_{c1}\\ X_{12}& Y_{12}& X_{22}& Y_{22}& X_{32}& Y_{32}& Q_{c2}\\ ...& ...& ...& ...& ...& ...& ...\\ X_{1n}& Y_{1n}& X_{2n}& Y_{2n}& X_{3n}& Y_{3n}& Q_{\mathrm{cn}}\end{array}\right)---\left(10\right)$

Wherein: X _ij, Y _ijrepresent x, y coordinate of jth mouth well in i-th individuality respectively, i=1,2..., j=1,2,3; Q _cibe i-th individual cumulative oil production, m ³.

In view of the multi-dimensional nature that well location is optimized, the present invention introduces Euclid's distance in linear algebra and represents the actual distance in any dimension space between 2.Suppose that largest optimization well number is m, namely Kriging method needs to assess the distance in n=2m dimension space between two points.According to the formula of Euclidean distance, obtain the some A (X in n-dimensional space _1i, Y _1i, X _2i, Y _2i..., X _mi, Y _mi) to a B (X _1i+1, Y _1i+1, X _2i+1, Y _2i+1..., X _mi+1, Y _mi+1) between distance be:

$h(A,B)=\sqrt{\underset{i=1}{\overset{m}{\mathrm{\&Sigma;}}}[{(x_{i1}-x_{i2})}^2+{(y_{i1}-y_{i2})}^2]}---\left(11\right)$

Wherein: h (A, B) represents the Euclidean distance of AB 2; X _jirepresent the x coordinate of jth mouth well in i-th individuality; Y _jirepresent the y coordinate of jth mouth well in i-th individuality.

Empirical covariance function following formula calculates:

$\mathrm{Cov}{\left(h\right)}_1=\frac{\mathrm{\&Sigma;}\left[z\left(x_i\right)z(x_i+h)\right]}{N\left(h\right)-1}-\frac{\mathrm{\&Sigma;z}\left(x_i\right)\mathrm{\&Sigma;z}(x_i+h)}{N\left(h\right)[N\left(h\right)-1]}---\left(12\right)$

Wherein: Cov (h) ₁, γ (h) represents empirical covariance function value, the experience variogram value z (x of sample point within distance h _i) be the functional value of sample point, in the present invention with adaptive response function about N (h) represents total sample point number within distance h, because general reservoir model has irregular shape, therefore along with distance h change.

S42, consideration nugget effect and partial weight, calculate theoretical covariance.

Choosing the precision of Kriging method prediction of theoretical variogram is extremely important.Conventional theoretical variogram has: linear, spherical, gaussian sum exponential model, in Holes of Complicated Wells bit optimization problem, adopts these functions to be difficult to reach the matching to experience variogram.In order to address this problem, the present invention introduces consideration block gold number effect, partial weight improves Kriging method.Basic ideas are as follows:

$\mathrm{Cov}{\left(h\right)}_2=\mathrm{\&rho;}\left(h\right){\widehat{\mathrm{\&sigma;}}}^2---\left(13\right)$

$\mathrm{\&rho;}\left(h\right)=2(1-\mathrm{nu}){\mathrm{\&rho;}}_{\mathrm{\&theta;}}\left(h\right){\left(\frac{\mathrm{bh}}{2}\right)}^{s-1}{K}_{s-1}\left(\mathrm{bh}\right)\mathrm{\&Gamma;}{(s-1)}^{-1}---\left(14\right)$

Wherein: Cov (h) ₂for theoretical covariance, ρ (h) is the spatial correlation function at a distance of h any two points, nu is block gold number, K is second-order modified Bessel function, Г is gamma function, parameter b, s are the real number characterizing relevance function shape, can obtain parameter b, s by solving following unconstrained optimization problem:

$\underset{b,s}{\min }\mathrm{\&Sigma;}{[2{\widehat{\mathrm{\&sigma;}}}^2(1-\mathrm{nu}){\mathrm{\&rho;}}_{\mathrm{\&theta;}}\left(h\right){\left(\frac{\mathrm{bh}}{2}\right)}^{s-1}{K}_{s-1}\left(\mathrm{bh}\right)\mathrm{\&Gamma;}{(s-1)}^{-1}-\mathrm{Cov}\left(h\right)]}^2---\left(15\right)$

ρ _θh () is partial weight function (or compact schemes related function), can be represented by the formula:

${\mathrm{\&rho;}}_{\mathrm{\&theta;}}\left(h\right)=\left\{\begin{array}{cc}(1-\frac{h}{\mathrm{\&theta;}})\cos \left(\frac{\mathrm{\&pi;h}}{\mathrm{\&theta;}}\right)+\frac{\sin \left(\frac{\mathrm{\&pi;h}}{\mathrm{\&theta;}}\right)}{\mathrm{\&pi;}},& h<\mathrm{\&theta;}\\ 0,& h\&GreaterEqual;\mathrm{\&theta;}\end{array}\right.---\left(16\right)$

Wherein: θ is supporting degree.

The empirical covariance that S43, basis calculate and theoretical covariance, utilize least square fitting covariance.

S44, linear relationship according to covariance and variogram, draw variogram.

S45, employing complete pivot elimination method solve system of linear equations, solve and obtain weight coefficient, obtain the Kriging method considering hyperspace, nugget effect and partial weight.

Evolution algorithm is based on Darwinian evolutionism thought, by simulation biological evolution process and the self-organization of the Solve problems of mechanism, adaptive artificial intelligence technology.As the breeding in biological evolution, variation, competition and selection phenomenon, evolution algorithm realizes solving of optimization problem mainly through selection, these three kinds of operations of recombinating and make a variation.Conventional evolution algorithm has genetic algorithm, particle cluster algorithm.

Embodiment of the present invention is improved genetic algorithm and particle cluster algorithm, utilizes the genetic algorithm after improving and particle cluster algorithm to carry out iterative computation auto-adaptive function to the new potentiality of remaining oil enrichment region determined.

The auto-adaptive function that the present invention selects cumulative oil production to optimize as well location, its functional relation is:

F＝f(I,J,T,S,N) (17)

In formula: F is cumulative oil production, it is the function of well parameter to be asked; I, J are respectively the grid index of straight well X, Y-direction; T is the classification of well, water injection well or producing well; S is the state of well, opens or closes; N is total well number.Well location optimization problem is and adopts a kind of sane optimized algorithm to determine the parameter combinations (I, J, T, S, N) making F maximum.

Revised genetic algorithum

Genetic algorithm is not adopt Deterministic rules, but adopts the transition rule of probability to instruct his direction of search.

Genetic algorithm (Genetic Algorithm, GA) is the overall Stochastic Optimization Algorithms of a kind of simple general-purpose, strong robustness, applied range, has imitated the evolutionary process of the survival of the fittest on biogenetics.Coding rule, initial population, crossover probability and mutation probability are the important parameters affecting genetic algorithm convergence and robustness, the present invention proposes to consider can optimize multivariable coding rule, based on the initial population of reservoir engineering experience and the improved adaptive GA-IAGA of adaptive crossover mutation and mutation probability.

(1) individual UVR exposure

The operand of genetic algorithm represents individual symbol string, so must be a kind of symbol string well parameter coding to be asked.The present invention adopts the binary coding be easily understood, and after arranging largest optimization well number, use binary coding simultaneously and well location coordinate other to the well of well to be optimized is encoded.

For different grid systems, the coding figure place of every mouthful of well is different.Be described to optimize a bite straight well, suppose that every mouthful of well needs 10 digit codes, first 2 represent the state (drive a well/closing well) of well and the classification (producing well/water injection well) of well, the mesh coordinate of rear 8 these wells of bit representation respectively.Such as: certain binary coding for certain individuality in population is 1111011111, is explained as follows:

1: mode bit, 1 expression drives a well, and 0 represents closing well, and therefore this well is opening.

1: the other position of well, 1 represents producing well, and 0 represents water injection well, and therefore this well is producing well.

1101: convert according to scale-of-two and the decimal system, I coordinate grid number I=13.

1111: convert according to scale-of-two and the decimal system, J coordinate grid number J=15.

Then genotype is the individuality of 1111011111, and the phenotype corresponding to it is for being positioned at the producing well of the 1 mouthful of opening in (13,15) grid place.

(2) initial population is determined

Initial population has important impact to genetic algorithm converges speed and accuracy.Initial population poor quality, greatly may increase the evolutionary generation of algorithm, cause counting yield to reduce, even likely Premature Convergence causes " precocity " to locally optimal solution.During determination initial population of the present invention, provide the individuality of 10% of initial population according to reservoir engineering experience, all the other produce individual by random device.Such initial method both ensure that the counting yield of algorithm, turn improved the rationality optimizing well location.

(3) fitness juice is calculated

Evaluate the good and bad degree of each individuality in genetic algorithm with the size of ideal adaptation degree, thus determine the size of its hereditary chance.In embodiment of the present invention, the fitness of Kriging method to individuality is utilized to assess.

(4) Selecting operation

Selecting operation is also called and copies computing, and individuality higher for fitness in current group is genetic in colony of future generation by certain rule or model.The individuality that General Requirements fitness is higher will have more chance to be genetic in colony of future generation.

(5) crossing operation

Crossing operation produces new individual main operating process in genetic algorithm, and with the individuality in the population after Selecting operation for operating basis, it exchanges the chromosome dyad between certain two individuality mutually with a certain crossover probability.

Crossover probability refers to that in chromosome, certain gene is selected as the probability of point of crossing, and higher crossover probability makes to evolve more abundant, and the local optimum ability of algorithm is strengthened, but speed of convergence is slack-off.The present invention adopts following adaptive crossover mutation, automatically regulates crossover probability according to fitness value change:

$p_x=\left\{\begin{array}{cc}\frac{a(F_{\max }-F_x)}{F_{\max }-F_{\mathrm{avg}}},& F_x\&GreaterEqual;F_{\mathrm{avg}}\\ b,& F_x<F_{\mathrm{avg}}\end{array}\right.---\left(18\right)$

In formula: F _maxfor maximum adaptation value in population, the present invention is cumulative oil production value; F _avgfor kind of a group mean adaptive value; F _xfor fitness value larger in two individualities will be intersected; A, b are constant, a=0.8, b=1 in the present invention.

(6) mutation operator

Mutation operator changes by a certain less probability the genic value on the some of individuality or certain some locus, and to select individuality in the population after intersecting for operating basis, it is also produce new individual a kind of method of operating.

The diversity of high mutation probability to Evolution of Population has positive role, but may destroy existing population.Early stage in evolution, low mutation probability should be adopted; Evolving late period, high mutation probability should be selected.So both avoided because high mutation probability covers early stage individual cross action, can be absorbed in convergence in advance in evolution and no longer make moderate progress to solution again, high mutation probability makes population have more evolvability and diversity.

Consider the change of optimum mutation probability with fitness, the present invention introduces self-adaptive mutation function, and mutation probability is set automatically according to fitness.The present invention selects mutation probability 0.01 in early days in optimization, when evolving to certain algebraically, when low mutation probability can not improve solution, selects high mutation probability 0.05.Self-adaptive mutation is calculated by following formula:

$p_m=\left\{\begin{array}{cc}\frac{c(F_{\max }-F_m)}{F_{\max }-F_{\mathrm{avg}}},& F_m\&GreaterEqual;F_{\mathrm{avg}}\\ d,& F_m<F_{\mathrm{avg}}\end{array}\right.---\left(19\right)$

In formula: F _mfor the fitness value of individuality that will make a variation; C, d are constant, c=0.01, d=0.05 in the present invention; All the other parameters are the same.

Embodiment of the present invention, adopts fitness function to carry out self-adaptation to crossover probability and mutation probability, obtains to carry out adaptive improvement evolution algorithm according to fitness function to evolution algorithm.Use the evolution algorithm improved to carry out iterative computation fitness function, merge the Kriging method set up and fitness function is assessed, define genetic algorithm (HGA).

Modified particle swarm optiziation

Compared with the genetic algorithm of maturation, particle cluster algorithm (Particle Swarm Optimization, PSO) is a kind of new evolution algorithm that development in recent years is got up.There is no the crossover and mutation of genetic algorithm unlike it, but particle follows optimal particle search in solution space, it is advantageous that simple easily realization and do not have many parameters to need adjustment.In order to strengthen overall exploring ability and the local development ability of PSO, the present invention introduces compressibility factor simultaneously and nonlinear inertial weight coefficient improves existing particle cluster algorithm.

In particle cluster algorithm particle displacement more new formula be:

$x_{i,j}^{k+1}=x_{i,j}^k+v_{i,j}^{k+1}\mathrm{\&Delta;t}---\left(20\right)$

For effectively controlling the flying speed of particle, the present invention introduces the particle rapidity more new formula considering compressibility factor:

In formula: be respectively the position of i-th particle when evolving to kth+1 time step and a jth element of speed; for a jth element of the optimum solution that i-th particle itself when evolving to kth time step finds; for a jth element of the globally optimal solution that whole population when evolving to kth time step is found at present; c ₁, c ₂for positive Studying factors; r ₁, r ₂for equally distributed random number between [0,1]; Δ t is time increment, and for ensureing algorithmic stability, Δ t is set to 1. for compressibility coefficient, be defined as:

Generally, arranging k is 1, Studying factors c ₁=c ₂=2.05, now C is 4.1, and calculating compressibility coefficient is 0.7298.

In order to overall exploring ability and the local development ability of balanced particle, the present invention introduces nonlinear inertia weight coefficient formula:

$w=\left\{\begin{array}{cc}{w}_{\min }-\frac{(w_{\max }-w_{\min })(F-F_{\min })}{(F_{\mathrm{avg}}-F_{\min })},& F\&le;F_{\mathrm{avg}}\\ w_{\max }& F>F_{\mathrm{avg}}\end{array}\right.---\left(23\right)$

In formula: ω is inertia weight; ω _min, ω _maxbe respectively maximal value and the minimum value of inertia weight, the present invention arranges ω _min=0.4, ω _max=0.9; F is the current target function value of particle; F _avg, F _minfor average criterion functional value and the minimum target functional value of current all particles.In embodiment of the present invention, a particle represents one group of pattern of well, is specially the well location of well number, every mouthful of well.

Embodiment of the present invention, adopts fitness function to carry out self-adaptation to inertia weight, obtains to carry out adaptive improve PSO algorithm according to fitness function to particle cluster algorithm.Use modified particle swarm optiziation to carry out iterative computation fitness function, merge the Kriging method set up and fitness function is assessed, define genetic algorithm (HPSO).

S5: according to evolution algorithm iterative computation, obtains maximum adaptation degree functional value, thus determines well location.

In the hope of dynamic potential index Quality Map as constraint condition, based on the numerical model corrected through history matching by oily district, other and well location is optimized to main force's fault block late high water content period optimal correction well number, well for the genetic algorithm (HGA) adopting two kinds of existing well bit optimization algorithm-genetic algorithms (GA) and particle cluster algorithm (PSO) and the present invention to obtain and genetic algorithm (HPSO).Arrange maximum evolutionary generation 50, the individual scale of initial population is 50, and constraint condition adopts the dynamic potential index constraint of four ranks, and arranging maximum well number according to reservoir engineering analysis is 3 mouthfuls.New brill producing well fixed output quota amount 50m ³/ d produces, and minimum sand face pressure is set to 15MPa; Injection rate IR 30m determined by new brill water injection well ³/ d water filling, maximum water injection pressure is 50MPa.The optimization time is 5 years (in the March, 2012 after terminating from history matching calculates in March, 2017).Table 2 be algorithms of different in various boundary conditions, namely DOI>0.2,0.4,0.6, the optimum results of 0.8.

The result of calculation of optimization 3 mouthfuls of wells of algorithms of different under the different well location constraint condition of table 2

As can be seen from Table 2, optimize constraint condition using the dynamic potential index DOI that the present invention proposes as well location, improve the search efficiency of well location optimized algorithm.Different Dynamic potential index constraint condition, computing time and speed of convergence are slightly different.Due under constraint condition DOI>0.2, algorithm search space is larger, cause random search in evolutionary process to invalid well location increase, therefore counting yield is lower, computing time is longer.

The cumulative oil production evolution curve of Fig. 8 ~ 11 each algorithm optimization well location under being respectively different well location constraint condition and the coordinate of optimization well location.

As shown in figs. 8 a and 8b, existing GA algorithm optimum solution obtains under DOI>0.6 constraint, and optimal well location is (21,18) and (8,9), and largest cumulative oil offtake is 15.54 × 10 ⁴m ³; As illustrated in figures 9 a and 9b, existing PSO algorithm optimum solution also obtains under DOI>0.6 constraint, and optimal well location is (64,12) and (6,28), and largest cumulative oil offtake is 15.69 × 10 ⁴m ³.

As shown in figures 10 a and 10b, HGA algorithm optimum solution of the present invention reaches in the 28th generation under DOI>0.6 constraint, and total computing time is 5 hours, search optimum well number 2 mouthfuls, be producing well, optimal well location is (62,16) and (6,9).

As shown in figures 1 la and 1 lb, HPSO algorithm optimum solution of the present invention reaches when DOI>0.6 equally, but only need to evolve to the 25th generation algorithm can restrain, total computing time is 5 hours, search optimum well number 2 mouthfuls, be producing well, optimal well location is (60,14) and (6,7).

Relative to existing two kinds of well location optimized algorithms, two kinds of algorithm counting yielies that the present invention proposes are higher, and the cumulative oil production calculating well location is higher, and as shown in figure 12, within 5 years, cumulative oil production is high by 2.04 × 10 ⁴m ³, show that the inventive method is reliable.HPSO algorithm of the present invention is identical for computing time with HGA algorithm, and cumulative oil production difference is within the error of calculation, and HPSO is slightly high, and optimal well location is not far from one another.Therefore it is all rational for digging a well in approximate range, shows that two kinds of algorithms that the present invention proposes are that computing velocity is fast, the well location optimization method of reliable results, in view of HPSO algorithm cumulative oil production of the present invention is higher, are therefore preferably final well location optimized algorithm.

Comprehensive analysis is determined, study area optimal correction well number is 2 mouthfuls, is producing well, and the grid of optimal well location in geologic model is respectively (60,14) and (6,7), and after 5 years, cumulative oil production can reach 17.73 × 10 ⁴m ³, geologic reserve recovery percent of reserves reaches 26.1%.

Theoretical foundation of the present invention is stronger, the situation such as oil water relation for High water cut maturing field complexity provides a kind of method determining High water cut maturing field determination well location, compared with existing well bit optimization method, the method can obtain higher cumulative oil production, water drive swept volume is larger, and optimize computation process quick, save time, speed of convergence increases greatly.Shown by instance analysis, the present invention determines that adjustment well location is rationally effective, has application value, provide scientific basis and computing method for the maturing field High water cut stage excavates residual production potential further to the reasonable arrangement of High water cut maturing field adjust well number well location.

The foregoing is only several embodiments of the present invention, those skilled in the art can carry out various change or modification to the embodiment of the present invention according to content disclosed in application documents and not depart from the spirit and scope of the present invention.

Claims (5)

1., based on the method for the high water cut oil field determination well location of evolution algorithm, it is characterized in that, comprise the following steps:

S1: the reservoir model setting up predetermined oil district, by described reservoir model grid division;

S2: the oil reservoir data obtaining described predetermined oil district, draws potentiality of remaining oil enrichment region according to described oil reservoir data according to history matching;

S3: the well location binding occurrence calculating each grid in described potentiality of remaining oil enrichment region, described well location binding occurrence represents remaining oil distribution and injects water swept volume, draws the grid that maximum well location binding occurrence is corresponding, described grid representation well location dominant area;

S4: use evolution algorithm to carry out its fitness function of iterative computation the oil well in described dominant area, and utilize Kriging method to assess fitness function, described evolution algorithm carries out self-adaptative adjustment according to fitness function, and described Kriging method considers hyperspace, nugget effect and partial weight;

S5: according to evolution algorithm iterative computation, obtains maximum adaptation degree functional value, thus determines well location.

2., as claimed in claim 1 based on the method for the high water cut oil field determination well location of evolution algorithm, it is characterized in that: in S3, determine that the step of well location dominant area is:

S31: determine the particle travel-time, the described particle travel-time utilizes following formula to calculate:

In formula: the travel-time that TOF (x, y, z) needs for particle in-position (x, y, z), unit is sky, ψ, χ are called double-current function, represent gradient; for reservoir porosity, s is the distance along streamline;

S32: on the basis in particle travel-time, calculate movable oil bulk index MOVI, oil phase fluid ability index CFI and streamline dynamic measurement DM, above-mentioned three indexes are calculated by following formula:

Movable oil bulk index MOVI:

MOVI＝(S _ro-S _orw)×PV×NTG

Oil phase fluid ability index CFI:

CFI＝K×K _ro×DZ×NTG

Streamline dynamic measurement DM:

DM＝(TOF _p+TOF _i)

In formula: S _ro: remaining oil saturation distribution when prediction starts, %;

S _orw: residual oil saturation distribution at the end of prediction, %;

PV: volume of voids, m ³;

NTG: net-gross ratio, represents net-to-gross ratio, dimensionless;

K: absolute permeability distributes, 10 ^-3μm ²;

Kro: oil relative permeability, refers in water-oil phase system, the ratio of permeability to oil and absolute permeability, dimensionless;

DZ: the size of longitudinal grid, m;

Above-mentioned parameter all obtains by obtaining oil reservoir data;

TOF _ifor taking water injection well as the particle travel-time that starting point calculates, TOF _pfor taking producing well as the particle travel-time that starting point calculates, DM is total travel-time;

S33: movable oil bulk index MOVI, oil phase fluid ability index CFI that S12 is obtained and streamline dynamic measurement DM Regularization, the concrete steps of Regularization are as follows:

1. by the 10% and 90% maximal value X determining attribute on model attributes cumulative distribution frequency plot _maxwith minimum value X _min;

2. according to the following formula Regularization is carried out to each attribute:

$X=\left\{\begin{array}{cc}0,& X<X_{\min }\\ \frac{X-X_{\max }}{X_{\max }-X_{\min }},& X_{\min }\&le;X<X_{\max }\\ 1,& X\&GreaterEqual;X_{\max }\end{array}\right.$

S34: based on three indexes of above-mentioned Regularization, utilizes following formula to calculate dynamic potential index DOI:

$\mathrm{DOI}=\sqrt[3]{\mathrm{MVOI}\&times;\mathrm{CFI}\&times;\mathrm{DM}}$

Dynamic potential index DOI is the concentrated expression of remaining oil distribution situation and well pattern sweep conditions, DOI more close to 1 region, show that these region remaining oil reserves are larger, crude oil fluidity better, and not yet by floood conformance, thus these regions become the dominant area boring adjust well.

3., as claimed in claim 1 based on the method for the high water cut oil field determination well location of evolution algorithm, it is characterized in that: in S4, consider that the step of the Kriging method of nugget effect and partial weight is:

S41: consider hyperspace, utilizes following formula to calculate empirical covariance:

$\mathrm{Cov}{\left(h\right)}_1=\frac{\mathrm{\&Sigma;}\left[z\left(x_i\right)z(x_i+h)\right]}{N\left(h\right)-1}-\frac{\mathrm{\&Sigma;z}\left(x_i\right)\mathrm{\&Sigma;z}(x_i+h)}{N\left(h\right)[N\left(h\right)-1]}$

Wherein: Cov (h) ₁the empirical covariance function value of sample point within expression distance h, z (x _i) be the adaptive response function related functions value with sample point, N (h) represents the total sample point number within distance h;

S42: consider nugget effect and partial weight, utilize following formula to calculate theoretical covariance:

$\mathrm{Cov}{\left(h\right)}_2=\mathrm{\&rho;}\left(h\right){\widehat{\mathrm{\&sigma;}}}^2$

$\mathrm{\&rho;}\left(h\right)=2(1-\mathrm{nu}){\mathrm{\&rho;}}_{\mathrm{\&theta;}}\left(h\right){\left(\frac{\mathrm{bh}}{2}\right)}^{s-1}{K}_{s-1}\left(\mathrm{bh}\right)\mathrm{\&Gamma;}{(s-1)}^{-1}$

Wherein: Cov (h) ₂for theoretical covariance, ρ (h) is the spatial correlation function at a distance of h any two points, for Kriging estimation variance, nu is block gold number, and K is second-order modified Bessel function, and Г is gamma function, and parameter b, s are the real number characterizing relevance function shape, can obtain parameter b, s by solving following unconstrained optimization problem:

$\underset{b,s}{\min }\mathrm{\&Sigma;}{[2{\widehat{\mathrm{\&sigma;}}}^2(1-\mathrm{nu}){\mathrm{\&rho;}}_{\mathrm{\&theta;}}\left(h\right){\left(\frac{\mathrm{bh}}{2}\right)}^{s-1}{K}_{s-1}\left(\mathrm{bh}\right)\mathrm{\&Gamma;}{(s-1)}^{-1}-\mathrm{Cov}\left(h\right)]}^2$

ρ _θh () is partial weight function, utilize following formula to calculate:

${\mathrm{\&rho;}}_{\mathrm{\&theta;}}\left(h\right)=\left\{\begin{array}{cc}(1-\frac{h}{\mathrm{\&theta;}})\cos \left(\frac{\mathrm{\&pi;h}}{\mathrm{\&theta;}}\right)+\frac{\sin \left(\frac{\mathrm{\&pi;h}}{\mathrm{\&theta;}}\right)}{\mathrm{\&pi;}},& h<\mathrm{\&theta;}\\ 0,& h\&GreaterEqual;\mathrm{\&theta;}\end{array}\right.$

Wherein: θ is supporting degree;

S43: according to the empirical covariance calculated and theoretical covariance, utilize least square fitting covariance;

S44: according to the linear relationship of covariance and variogram, draw variogram;

S45: adopt complete pivot elimination method to solve system of linear equations, solve and obtain weight coefficient, complete the foundation of Kriging method.

4. as claimed in claim 1 based on the method for the high water cut oil field determination well location of evolution algorithm, it is characterized in that: the evolution algorithm described in S5 comprises genetic algorithm, crossover probability and the mutation probability of described genetic algorithm automatically adjust according to auto-adaptive function:

Utilize following formula that crossover probability is regulated automatically according to fitness value change:

$p_x=\left\{\begin{array}{cc}\frac{a(F_{\max }-F_x)}{F_{\max }-F_{\mathrm{avg}}},& F_x\&GreaterEqual;F_{\mathrm{avg}}\\ b,& F_x<F_{\mathrm{avg}}\end{array}\right.$

In formula: F _maxfor maximum adaptation value in population, F _avgfor kind of a group mean adaptive value, F _xfor fitness value larger in two individualities will be intersected, a=0.8, b=1;

Utilize following formula that mutation probability is regulated automatically according to fitness value change:

$p_m=\left\{\begin{array}{cc}\frac{c(F_{\max }-F_m)}{F_{\max }-F_{\mathrm{avg}}},& F_m\&GreaterEqual;F_{\mathrm{avg}}\\ d,& F_m<F_{\mathrm{avg}}\end{array}\right.$

In formula: F _maxfor maximum adaptation value in population, F _avgfor kind of a group mean adaptive value, F _mfor the fitness value of individuality that will make a variation, c=0.01, d=0.05.

5. as claimed in claim 1 based on the method for the high water cut oil field determination well location of evolution algorithm, it is characterized in that: the evolution algorithm described in S5 comprises particle cluster algorithm, the particle displacement of described particle cluster algorithm upgrades and flying speed automatically adjusts according to auto-adaptive function:

$x_{i,j}^{k+1}=x_{i,j}^k+v_{i,j}^{k+1}\mathrm{\&Delta;t}$

In formula: be respectively the position of i-th particle when evolving to kth+1 time step and a jth element of speed; for a jth element of the optimum solution that i-th particle itself when evolving to kth time step finds; for a jth element of the globally optimal solution that whole population when evolving to kth time step is found at present; c ₁, c ₂for positive Studying factors; r ₁, r ₂for equally distributed random number between [0,1]; Δ t is time increment, for compressibility coefficient:

ω is inertia weight:

$w=\left\{\begin{array}{cc}{w}_{\min }-\frac{(w_{\max }-w_{\min })(F-F_{\min })}{(F_{\mathrm{avg}}-F_{\min })},& F\&le;F_{\mathrm{avg}}\\ w_{\max },& F>F_{\mathrm{avg}}\end{array}\right.$

In formula: ω _min, ω _maxbe respectively maximal value and the minimum value of inertia weight, F is the current adaptive value of particle; F _avg, F _minfor average criterion adaptive value and the minimum adaptive value of current all particles.