Method and device for optimizing comprehensive production of oil reservoir

Abstract

A method and device for optimizing comprehensive production of an oil reservoir. The method comprises: on the basis of historical production data, performing regional level recognition, single-well level recognition, and production optimization model recognition, and updating oil reservoir dynamic recognition; on the basis of the historical production data, using the following steps to perform cyclic historical fitting on model parameters of a plurality of production optimization models: performing historical fitting on the model parameters, and performing well-to-well communication relationship calibration on the plurality of production optimization models; determining, on the basis of the updated oil reservoir dynamic recognition, a constraint condition corresponding to each production optimization model, solving an objective function corresponding to each production optimization model, obtaining an optimal decision variable corresponding to each production optimization model, and obtaining a development regulation scheme corresponding to each production optimization model; obtaining a comprehensive development regulation scheme of an oil reservoir; and updating the oil reservoir dynamic recognition on the basis of a regulation effect obtained by real-time monitoring of a development operation.

Description

Method and device for comprehensive production optimization of oil reservoirs

[0001] Related Applications

[0002] This application claims priority to and incorporates by reference the entire disclosure of Chinese Patent Application No. 202410977054.7 filed on July 19, 2024. TECHNICAL FIELD

[0003] The present application relates to the technical field of oil and gas field development engineering, and in particular to a method and device for comprehensive production optimization of oil reservoirs. BACKGROUND

[0004] This section is intended to provide background or context to the embodiments of the application recited in the claims. The description herein does not constitute admission that the prior publication, square, or subject matter described herein and / or the corresponding physical item is prior art.

[0005] Improving oil recovery and economic benefits is the ultimate goal of oilfield development, and closed-loop reservoir management is an effective way to achieve this goal. However, the production optimization process in the current closed-loop reservoir management faces serious experience dependence, low precision and efficiency, and difficulty in solving constraint optimization problems due to the large number of control parameters and long analysis or simulation time. To effectively solve these problems and improve the intelligent level of oil reservoir production optimization, it is necessary to develop efficient solution algorithms for time-consuming constraint optimization problems, to study simplified physical models to solve simulation time-consuming problems and numerical simulation models to solve precision (reliability) problems, and to cross-fuse the two, to comprehensively evaluate the advantages and disadvantages of various optimization models (production optimization methods), and then to propose a comprehensive production optimization concept and workflow that balances reliability and calculation speed, considers real-time and long-term, and gives equal weight to physical meaning and data-driven, to improve the intelligent level of reservoir management.

[0006] For a long time, reservoir management is a continuous and long-period research task that needs to integrate the processes of geological modeling, history matching, reservoir simulation, production optimization, plan implementation, dynamic monitoring and analysis, and effect feedback. Intelligent closed-loop reservoir management, which integrates these research tasks, has gradually become a research hotspot in the field of oil and gas field development due to its significant achievements in efficiency improvement and precision improvement. Production optimization is an eternal theme of closed-loop reservoir management. It applies optimization theory to the process of reservoir development and management, uses a certain physical or numerical model, and uses optimization algorithms to make the given objective function optimal under the condition of meeting various complex constraints of the reservoir system, determine the optimal decision variable, and thus give the best development control plan. The core elements of production optimization include decision variables, objective functions, constraint conditions, optimization models, and optimization algorithms. For new oil fields, the objects of production optimization often involve the development methods, layer series, well pattern, and injection-production systems specified in the development plan. For old oil fields, it may also need to optimize the infill well location, injection-production parameters (well control parameters), and adjustment parameters of measures. When the development method, layer series, and well pattern of the oil field are determined and there is no new drilling or measure workload, the adjustment of injection-production parameters (work system) is the most frequent.

[0007] The most commonly used traditional production optimization methods in current reservoir development adjustment mainly include two categories: reservoir engineering methods and numerical simulation methods. Reservoir engineering methods make heuristic adjustments to injection-production systems based on theoretical derivation, dynamic analysis, geological reservoir understanding, and empirical knowledge from past development practices. The adjustment effect is used to optimize the adjustment idea and update the empirical knowledge. The advantages of this method are: simple and practical, strong operability, low data demand and time cost, and real-time performance. However, this method also has obvious disadvantages, such as: excessive dependence on the knowledge depth of the researchers and their understanding or recognition level of the target reservoir, different researchers may draw different conclusions from the same data source; only a rough optimization adjustment direction can be given, with low precision; and when the number of wells in the reservoir is large, the analysis efficiency becomes low.

[0008] The numerical simulation method is based on the fine geological modeling, and the dynamic data is used to carry out history matching of the numerical model. Due to the limited time and cost, researchers often set several simulation schemes according to the comprehensive understanding of the reservoir, and then select the optimal solution according to the simulation results. These methods can integrate seismic, logging and production dynamic data, and have high reliability. However, the data requirement of the history matching process is large, the time-consuming is too long, and it is difficult to determine the optimal solution for the subsequent prediction and optimization. Even if the history-matched numerical model is called by optimization algorithms such as gradient-based algorithms (such as finite difference gradient, adjoint gradient), approximate gradient-based algorithms (such as simultaneous perturbation stochastic approximation (SPSA), stochastic Gaussian search direction (SGSD) or Gaussian simultaneous perturbation stochastic approximation (G-SPSA), ensemble-based optimization (EnOpt) and stochastic simplex approximate gradient (StoSAG)), or classical non-gradient algorithms (such as particle swarm optimization (PSO), genetic algorithm (GA), differential evolution (DE), pattern search (PS), simulated annealing (SA) and the like), the calculation cost is still high, and the optimization effect is still improved under the limited simulation times.

[0009] Due to the difficulties of traditional production optimization methods such as serious empiricism, long time-consuming and high calculation cost, the later production optimization research mainly focuses on improving the optimization calculation efficiency and calculation accuracy. In order to overcome the disadvantages of the empiricism of pure reservoir engineering method and the time-consuming and laborious of traditional numerical simulation method, many researchers have proposed a reduced-complexity physical model (RCPM) which is not dependent on geological modeling but only subject to the main physical laws. The model parameters can be obtained by history matching of daily production data, and the history matching and dynamic prediction efficiency is greatly improved.

[0010] The pioneering research on simplifying physical model is the development of the Capacitance Resistance Model (CRM) or Capacitance Model (CM). CRM can be divided into three types according to the control volume, which are Capacitance-Resistance Model for Total system (CRMT), Capacitance-Resistance Model for Production wells (CRMP) and Capacitance-Resistance Model for Injection-Production well pairs (CRMIP). CRM is based on the multiple linear regression model, and its essence is the material balance response of fluid injection and production data in a certain time, which can be used to represent the conversion relationship between the input and output signals of the reservoir system. CRM is widely used in the evaluation of interwell connectivity in water drive reservoirs. Later, many researchers introduced the fractional flow model (such as the power-law fractional flow model, the Brooks-Corey relative permeability index fractional flow model, the Koval model, and the Y function) to predict oil and water production and optimize water injection. Some researchers also extended CRM to multi-layer reservoirs. In addition, the uncertainty of the model can be reduced by establishing a CRM error system.Comprehensive analysis shows that current researches based on CRM have made some achievements in model development, connectivity evaluation, production prediction and production optimization, especially in the field of mature water drive reservoir management. However, there are still some problems: (1) Most of the researches focus on a specific stage of water drive development (such as the middle and high water cut stage), and there is no model applicable to the whole water drive process; (2) The changes of fluid properties and model parameters are not fully considered, and the time-varying characteristics are not explained under the constraint of clear physical meaning; (3) There is insufficient research on the situation of horizontal well injection and production or mixed injection and production of vertical and horizontal wells; (4) For strong heterogeneous reservoirs, the water cut change law of production wells is complex, and the current flow rate model coupled with CRM is difficult to accurately predict oil and water flow rate, so a simple and practical phase flow estimation method needs to be developed; (5) Connectivity evaluation and real-time production optimization under complex injection and production system are still a great challenge, and for the situation of injection and production well shutdown, bottom hole flowing pressure fluctuation and other large flow field changes, the current CRM has no better solution; (6) The research on gas reservoir development or gas injection development of oil reservoirs is still not in-depth; (7) There is no systematic research on the feedback mechanism of CRM and other production optimization methods.

[0011] Interwell Numerical Simulation Model (INSIM) is another typical representative of simplified physical model; this model is established on the basis of network model and ICNS (Interwell Connectivity Numerical Simulation) model, similar to CRMIP, its basic principle is to regard the reservoir as a simplified system composed of many injection-production control units, each control unit contains two parameters of transmissibility and control volume, combined with material balance equation, Buckley-Leverett front tracking equation and fractional flow model can calculate phase flow and saturation, using injection-production data can carry out history matching to invert model parameters, and then evaluate interwell connectivity and infer geological characteristics. After that, many researchers have explored the development and application expansion of INSIM type method. Some researchers have established a multi-layer reservoir interwell connectivity inversion model based on connectivity model and INSIM, also considered the situation of wellbore fluid diversion caused by measures such as shutting down wells and switching injection, improved the saturation tracking method of ICNS connectivity model, and proposed a production optimization method based on interwell connectivity evaluation. Some researchers have developed INSIM-FT model by using front tracking (FT) method to calculate saturation based on INSIM. Some researchers have used the concept of distribution coefficient and injection efficiency to apply INSIM to water drive history matching, interwell connectivity characterization and production optimization. Some researchers have extended INSIM-FT to three-dimensional multi-layer reservoirs to overcome the shortcomings of INSIM-FT which only considers two-dimensional flow of vertical wells in single-layer reservoirs, considering the effects of gravity and well trajectory, and then developed INSIM-FT-3D. Some researchers have developed INSIM-FPT (Inter-well Numerical Simulation Model with Flow-Path Tracking) model to more accurately calculate interwell connectivity coefficient and control volume, track flow path and make water drive prediction. After that, some researchers have extended INSIM-FPT to three-dimensional reservoir application, proposed INSIM-FPT-3D model and a water drive optimization method based on oil saturation and injection distribution coefficient. Some researchers have proposed Connection Element Method (CEM) based on INSIM, which discretizes the reservoir calculation domain into a series of connection units, uses two parameters of connection transmissibility and connection volume to reflect the percolation capacity and control range of connection units. When the nodes are only well points, CEM is equivalent to INSIM type method.To overcome the shortcomings of previous INSIM, which uses upstream saturation (or upstream weight) to obtain the conductivity of each connection unit, ignores the changes of total mobility and gravity during the propagation of waterflood front, and cannot accurately generate the bottom hole pressure, a new physics-based data-driven interwell three-dimensional (3D) numerical simulation model with front tracking (INSIM-BHP) was proposed, which realizes the calculation of harmonic average conductivity and fine-scale gravity term based on the fine-scale saturation distribution calculated by front tracking, accurately represents the flow between connection nodes, and combines water body modeling to quickly and accurately calculate the bottom hole pressure. Comprehensive analysis shows that the current research based on INSIM has achieved certain results in model development, simplified numerical simulation, connectivity evaluation, production optimization, and other aspects, especially in the field of waterflood reservoir management, and has obtained considerable success. In the future, it will gradually develop towards the simulation of complex reservoir media (such as faults, fractures, unconventional reservoirs, and multiple media) and the simulation of complex fluid flow (such as three-phase problems of gas injection, water-gas alternation, etc.) and the development of flow network model. However, the following aspects still need to be continuously explored: (1) For complex well types such as horizontal wells and deviated wells or various well types combined with injection and production, the inversion of connectivity relationship is not accurate enough, and the setting of nodes and flow paths still needs to be improved; (2) Due to the speed advantage brought by model simplification, the inversion of node pressure field and saturation field between nodes is not accurate enough under the flow rate history matching; (3) For the case of complex injection and production system and production history of multi-layer reservoirs with layer conversion and injection, the history matching has strong multi-solution property, and the uncertainty of connectivity evaluation and injection-production optimization results needs to be further reduced; (4) The interlayer fluid exchange has not been effectively considered; (5) The numerical simulation and production optimization application for gas reservoir development or oil reservoir injection and development, water-gas alternation, etc. still needs to be improved; (6) The parameter feedback and verification mechanism between other simulation methods with higher reliability has not been systematically studied.

[0012] In addition to the CRM and INSIM-type models, the Flow Network Model (FNM) is another representative simplified physical model developed in recent years; it further introduces the step of grid division on the basis of CRMIP and INSIM, simplifying the multi-dimensional reservoir space into a series of one-dimensional flow networks between each well pair, coupling each one-dimensional model at the intersecting nodes, and the flow in each connection depends on the fluid volume and the corresponding reservoir permeability in the displacement region covered by the two wells. The flow network model can be calibrated by adjusting these two parameters (fluid volume and corresponding reservoir permeability) using the output of the full-order reservoir model or field observation data, and then implementing production optimization. Since each one-dimensional network uses two-phase flow control equations for finite difference simulation, its model parameters are much fewer than those of the full-order simulation model, so the calculation speed is improved. FNM can be used to characterize interwell connectivity, estimate the displacement volume between well pairs, and predict reservoir production over a long period of time for the purpose of production optimization. Some researchers have used a one-dimensional flow path network connecting different well perforation locations to model the reservoir and establish an equivalent two-dimensional Cartesian model, and proposed a physically-based data-driven model GPSNet using a commercial simulator to solve the problem of severe calculation constraints in closed-loop reservoir management using full-physical models for history matching and optimization, which discretizes the one-dimensional flow path and uses the ESMDA (Ensemble Smoother with Multiple Data Assimilation) method to history match production data, thereby obtaining the attribute parameters of each grid block along each flow path; they tested GPSNet using water drive and steam drive examples. Some researchers added additional nodes to the FNM to represent indirect flow paths and aquifers, and combined the StoSAG (Stochastic Simplex Approximate Gradient) algorithm to apply it to well control optimization. Some researchers reconstructed the physically-based data-driven model GPSNet, and used the GPSNet flow network model combined with a commercial simulator to conduct history matching and water drive optimization research using a partial model of an actual oilfield as an example. Some researchers also proposed a flow network model (General Physics-Based Data-Driven Framework, GPDF) based on the INSIM-type model through grid division on the simplified one-dimensional connection units of the reservoir, which also uses the SPSA algorithm to solve the characteristic parameters of the connection units. Since it retains the simplified physical properties of the INSIM method and the generality of traditional grid-based simulation methods, its prediction accuracy and applicability have been enhanced.On the basis of INSIM and flow network model (FNM), some researchers proposed a simplified physical model FlowNet for history matching and production prediction of fractured shale or tight reservoirs. Similar to FNM, this model simplifies the real reservoir into a flow network composed of a series of well nodes (fracture cluster nodes, fracture nodes, SRV nodes and matrix nodes) and one-dimensional connection units between them. The flow equations in the FlowNet grid system are solved by using a fully implicit nonlinear solver, and the parameters such as pressure, phase saturation and production can be obtained. The efficient history matching program based on the FlowNet model of fractured reservoirs can determine the parameters of each connection unit, so as to realize the rapid prediction of production. Based on INSIM and FNM, some researchers defined the parameters such as conductivity, control volume, fractal mass dimension and fractal index in each one-dimensional connection unit to map the reservoir properties, considered the fractal characteristics of reservoir permeability and porosity, combined the data-driven physical model with fractal theory, and then proposed a data-driven reservoir simulation framework based on fractal physics (FlowNet-fractal). On the basis of the previous flow network model GPSNet, some researchers discretized the reservoir into two-dimensional connections (x-z plane) between a series of completion sections, developed a GPSNet-2D that can be coupled with commercial simulators, and applied it to the history matching (Ensemble Smoother with Multiple Data Assimilation, ES-MDA) and robust optimization (genetic algorithm, GA) of steam flooding in a part of the thermal recovery model. This model can better depict the plane / vertical flow path and flexibly consider the complex physical phenomena in the thermal recovery process (such as steam separation near the injection well, rapid steam breakthrough at the production well, and heat loss of overlying / underlying formations, etc.). Some researchers applied GPSNet to a water drive oilfield containing thousands of wells, used its internal simulator and ESMDA algorithm for history matching, and based on the history-matched model, used the GA algorithm to optimize the injection rate and injection-production ratio. In addition, in recent years, research on other flow network models has gradually emerged. For example, some researchers proposed a reservoir graph network (RGNet) modeling method based on discrete flight time grid division to solve the problems of traditional numerical simulation model construction and calibration being cumbersome and the cost of short-term decision-making being too high, and pure data-driven models lacking physical insights. This method simplifies the three-dimensional reservoir flow problem into a graph network representation problem, which can be solved by using commercial reservoir simulators. After model calibration, it can also be used for real-time prediction, scenario analysis and production optimization. Some researchers used ES-MDA algorithm based on RGNet to perform history matching on an oilfield example, and used evolutionary algorithm to optimize injection and production rates.Some researchers have applied the reservoir graph network (RGNet) model to history matching and dynamic prediction of unconventional gas reservoirs. Some researchers have applied the RGNet modeling method to a mechanistic example and an oilfield example, considering well interconnection generation and a more stable pressure history matching method, using ES-MDA and an adjoint-based method for history matching (ES-MDA can perform uncertainty analysis, while the adjoint method requires fewer simulation times), and the history matching parameters include well indices, well internal and external cell conductivities, and pore volumes. Some researchers have considered the pressure dependence of conductivities and pore volumes, and have modeled and predicted multi-well production in unconventional reservoirs using RGNet. Comprehensive analysis shows that the advantages of FNM are to reduce model parameters while retaining the basic physical meaning of the original full-order reservoir model, improve operational efficiency, and directly use field data for calibration. The simulation can also be coupled with commercial simulators. However, the disadvantages of FNM are that the generation of flow networks requires manual setting, it is essentially a one-dimensional finite difference simulation based on a grid, and certain grid attribute information is required in advance; and there is a multi-solution and greater uncertainty in history matching, and the application of FNM in complete closed-loop reservoir management still needs further research.

[0013] The development of the above-mentioned simplified physical models from CRM to INSIM to FNM reflects the balance between the pursuit of production optimization speed (efficiency) and accuracy (determinism, reliability) by different researchers. In summary, although simplified physical models effectively improve the time-consuming and laborious shortcomings of traditional production optimization methods, they also generally have greater uncertainty and weaker long-term prediction ability, and are more suitable for short-term real-time optimization. In view of this, many researchers have also made beneficial explorations for other production optimization methods with higher reliability (such as accelerated optimization relying on numerical simulators). Such research mainly focuses on improving the speed or efficiency of traditional numerical simulation optimization methods, and optimization methods based on machine learning agent modeling belong to one of them.

[0014] In the field of reservoir production optimization, surrogate models can be roughly divided into two categories: reservoir simulator surrogate and data-driven surrogate. For the former, its surrogate optimization algorithm involves three cores: initial experimental design, surrogate model selection, and sampling method or strategy (the criteria or basis for selecting potential candidate points). In this case, a numerical simulator is used to obtain the true response between the sample (a set of decision variables) and the production optimization objective function, so the objective function evaluation and surrogate model updating are both "online". If the objective function values of some samples are known, they are called "offline" samples, which can be directly included in the initial sampling range. In addition, since the sampling method (or strategy) will continuously select new valuable samples (potential candidate points) to guide the search process, the surrogate model is constantly updated, that is, the two processes of surrogate construction and potential candidate point determination are alternately performed. For the latter, i.e. data-driven surrogate optimization based on reservoir production data, the surrogate model should as accurately as possible reflect the input-output relationship between the reservoir production data, which needs to be used for training after quality control of the reservoir production data, and the trained surrogate model is used to implement production optimization combined with optimization algorithms. At this time, the training of the surrogate model should not only accurately represent the production data response of the underground reservoir system, but also prevent the problem of "overfitting" and poor generalization ability.

[0015] The main difference between reservoir simulator surrogate and data-driven surrogate is that: (1) the essence of the former is sampling, and the surrogate is updated with optimization iteration, which does not need to be accurate enough to approximate the true model at the beginning; while the essence of the latter is history fitting, and the key is the arrangement and control of training data, and the surrogate is not updated with optimization iteration, which needs to be accurate enough to approximate the true model. (2) The two processes of surrogate construction and production optimization of the former are coupled together; while for the latter, the two processes are separated, and production optimization is performed after the construction of the surrogate model with accurate history fitting results and prediction results.

[0016] Machine learning can be generally divided into supervised learning, unsupervised learning or a combination of the two, according to whether the data response or output is known (whether the training data is labeled). Supervised learning trains the model according to known input and output, so that the model can predict future output, which is generally used for classification and regression analysis. Unsupervised learning only groups and explains data based on input data, aiming to find hidden rules or internal structure from input data, which is generally used for data dimensionality reduction, feature extraction, clustering analysis, etc. The surrogate modeling method in reservoir production optimization belongs to supervised learning, aiming to build a surrogate model that can reasonably predict the objective function (NPV (Net Present Value), cumulative oil production, etc.) based on uncertain production data.

[0017] Common proxy models mainly include polynomial response surface, multivariate adaptive regression splines (MARS), Gaussian process regression (GPR) or Kriging model, random forest, decision tree, support vector machine regression, generalized additive model (GAM), radial basis function (RBF) interpolation model, artificial neural network (ANN, NNs), etc. Among them, there are many types of neural network models, which can be roughly divided into feedforward networks (also known as multilayer perceptron networks, such as convolutional neural network CNN, BP (Backpropagation) neural network, RBF neural network, etc.), feedback networks (such as recurrent neural network RecurrentNN, recursive neural network RecursiveNN, Hopfield network, Boltzmann machine, long short-term memory network LSTM (Long Short-Term Memory), etc.), and graph neural networks (such as graph convolution network, graph autoencoder, graph generation network, graph recurrent network, graph attention network), etc. According to the number of network layers, it can be divided into shallow neural network and deep neural network, and according to the learning method, it can be divided into supervised learning, unsupervised learning, semi-supervised learning and reinforcement learning, etc.

[0018] Machine learning surrogate modeling usually generates sample data based on traditional numerical simulator or real reservoir data, and trains various machine learning models (such as Radial Basis Function Network (RBFN), Recurrent Neural Network (RNN), Long Short-Term Memory Network (LSTM), Deep Reinforcement Learning (DRL), Convolutional Encoder-Decoder Network (CED), Convolutional Neural Network (CNN), Gaussian Process Regression (GPR), Kriging, Support Vector Machine for Regression (SVMR) / Support Vector Regression (SVR), Echo State Network (ESN), and other models or algorithms or their mutual combinations) to replace the simulator or real reservoir for objective function evaluation or to find reasonable points to be evaluated to gradually approach the real model, thereby realizing surrogate optimization.

[0019] For surrogate optimization based on reservoir simulator, in addition to surrogate selection, initial experimental design method and sampling method (or strategy) are the other two cores of surrogate optimization algorithm. The initial experimental design methods mainly include Latin hypercube design and its various variants, Halton sequence sampling design, Sobol sequence sampling design, and various adaptive sampling, and existing literature has systematically introduced various experimental design methods. However, a good surrogate optimization algorithm should not be sensitive to the initial experimental design method, and therefore the focus of surrogate optimization research falls on the exploration of subsequent sampling method (strategy), which selects the next batch of valuable function evaluation points (potential candidate points) by predicting the objective function value of the unsampled area through the surrogate model and combining the current function evaluation situation, and uses them for surrogate update. The representative sampling method (or strategy) is the utility function defined by radial basis function interpolation and the Metric Stochastic Response Surface (MSRS) algorithm. For optimization problems with very expensive function calculation and no additional information available, some researchers have proposed a method to solve continuous non-convex functions in R dThe method of global minimum on compact subsets, which uses radial basis function interpolation to define a utility function, the maximum of which corresponds to the next point at which the objective function is evaluated. For the global optimization of computationally expensive functions, a stochastic response surface(SRS) algorithm was proposed, which iteratively uses response surface models to approximate the expensive function. The method of identifying candidate points from a randomly generated set of points for the next function evaluation to update the response surface was given. The metric SRS(MSRS) algorithm was proposed, which uses two criteria to select the function evaluation point(best candidate point): the estimated function value obtained from the response surface model and the minimum distance from the previous evaluation point. The global optimization version and the multi-start local optimization version of MSRS were developed. The previous algorithms that use surrogate models to optimize expensive functions can only handle bounded constraint problems, so these surrogate models are only used to approximate the objective function and cannot effectively handle nonlinear(inequality) constraints. Based on the local metric stochastic radial basis function(LMSRBF) algorithm, a new derivative-free optimization algorithm for expensive black-box objective functions with expensive black-box inequality constraints was proposed, namely ConstrLMSRBF. This algorithm constructs RBF surrogate models for the objective function and all constraint functions in each iteration and uses these RBF models to guide the selection of the next objective function and constraint function evaluation points. Based on a new distance definition, the metric stochastic response surface(MSRS) algorithm was modified and extended to be applied to computationally expensive high-dimensional black-box optimization problems. A new method for generating candidate points was proposed and the corresponding global convergence was proved. The modified MSRS algorithm(referred to as Surrogate Optimization-Sensitivity Analysis, SO-SA) is more adaptable and more likely to perturb the most sensitive coordinates rather than all coordinates simultaneously when generating candidate points. The influence of surrogate model selection and sampling point selection methods on the quality of surrogate solutions was also explored. For expensive optimization problems(EOPs) with high evaluation cost of candidate solutions, surrogate-assisted evolutionary algorithms(SAEAs) that can effectively reduce computational cost and improve solution efficiency were systematically introduced. According to the types of objective functions and constraints, existing SAEAs were classified and discussed. The research results on SAEAs were summarized from the aspects of algorithm and application.For more valuable research on sampling strategies, see Shoemaker, Regis, Mueller, et al.

[0020] The research finds that the machine learning agent model is trained by using production data generated by a numerical simulator or field production data in combination with a sampling method and an optimization algorithm for one or more optimization targets, so as to replace the simulator or the real reservoir to realize the prediction and optimization process by using the machine learning model. For the oil reservoir production optimization problem with time-consuming function evaluation and difficult gradient calculation, the agent optimization method can reduce the calculation cost to a certain extent. The difference between the agent optimization methods based on the simulator mainly lies in the initial test design, the selection of the agent model, and the selection strategy of the subsequent function evaluation points (potential candidate points), and the difference between the pure data-driven agent optimization methods mainly lies in the construction / selection of the agent model and the optimization algorithm.

[0021] For the agent optimization method based on the simulator, the main advantages are that the agent can represent the physical connotation of the numerical simulator, the simulation cost after the initial agent is constructed is significantly reduced, and the optimization (or seeking potential candidate points) process no longer needs to call the simulator for operation, which can greatly reduce the time cost. It should be noted that this method is not suitable for the "training first and then optimizing" path in the pure data-driven agent optimization, because the data set acquisition process is too time-consuming. The main problems existing in this method include: (1) due to the sensitivity of the required initial sample to the dimension of the optimization problem and the involvement of multiple iterative function evaluations, the time consumption of the initial agent model construction and the subsequent agent update process is still generally long; (2) even if the training process and the optimization process can realize benign interaction and accelerate the optimization solution, the optimization result still depends on the quality or accuracy of the historical fitting model of the numerical simulator, and the realization of the production optimization closed loop still faces challenges; (3) the influence of the geological model uncertainty has not been effectively considered; (4) the synergistic effect with other oil reservoir production optimization methods has not been effectively exerted; (5) for high-dimensional complex constraint optimization problems with linear, nonlinear (inequality and equality) constraints, the efficient synergistic cooperation of the agent optimization in agent construction, sampling strategy, optimization algorithm, and constraint processing method still faces challenges.

[0022] For the pure data-driven agent optimization method, the main advantage is that the construction of the agent and the subsequent production optimization process are in a separated state, so that the time cost of executing the optimization process after the agent is constructed can be almost negligible compared with running the simulator. The main problems existing in the pure data-driven agent optimization method include: (1) the workload of production data collection is huge, and the data quality is difficult to guarantee, and the data noise may seriously affect the reliability of the model construction; (2) it is difficult to describe the change of the underground reservoir system and the physical meaning reflected by the fluid flow, and the interpretability of the optimization result is poor.

[0023] In addition to proxy modeling production optimization, streamline-based production optimization is another promising method for reservoir production optimization, which mainly relies on the speed advantage of streamline simulation and the unique physical meaning of the seepage field endowed by streamline to carry out injection-production flow allocation optimization research. In the field of streamline simulation, the first representative production optimization method is injection efficiency analysis. Some researchers first proposed the concept of injection efficiency of injection wells or injection-production well pairs according to the characteristics of streamline simulation that can determine the dynamic allocation coefficient between injection and production wells and the flow corresponding relationship. According to the relative size of injection efficiency and average injection efficiency, some principles can be designed to predict the appropriate flow target of injection and production wells, and to redistribute the injected water from low-efficiency injection wells to high-efficiency injection wells, thereby improving water drive management. Later, many researchers developed this concept and optimization idea. For example, some researchers established a simulation model of naturally fractured reservoirs based on the injection efficiency derived from streamline simulation to improve the water drive management of naturally fractured reservoirs in view of the complex fluid exchange behavior between fractures and matrix in naturally fractured reservoirs. Some researchers explained the concepts of well pair efficiency and generalized injection-production well pair efficiency based on net present value (the product of three efficiencies) in view of the flow allocation optimization problem in the process of water / gas injection and production, and proposed a fast and stable derivative-free workflow based on streamline simulation, which can improve economic value by optimizing the flow allocation of water drive and gas drive; they mainly consider reservoir performance and economic value in the objective function, and use static, dynamic and economic parameters (such as price and discount rate and flight time) to evaluate the expected net present value of each injection-production well pair under given future business decisions (such as the time of next infill); then, according to the expected net present value, the performance of each well is sorted, and the flow of each well is redistributed, so as to maximize future economic benefits.

[0024] Another representative streamline-based production optimization method is the time-of-flight analysis. For example, a method for determining the optimal injection-production rates was proposed for large-scale flow optimization problems under actual field conditions based on streamline tracking and time-of-flight (TOF) calculation. The waterflood front was controlled by adjusting the injection-production rates to improve sweep efficiency and delay the water breakthrough time. The basic idea is to make the waterflood front reach all production wells at the same time in the selected subzone of the waterflood project. The optimization of the TOF has good quasi-linear properties, and the optimization process can proceed smoothly even if the initial conditions are far from the solution. In addition, the sensitivity of the TOF to injection and production rates can be calculated analytically using a single flow simulation. Later, a method for optimizing injection-production rates was proposed based on the same idea of maximizing sweep efficiency by balancing the waterflood front arrival time at all production wells. Streamlines were used to efficiently and analytically calculate the sensitivity of the TOF to well rates, and a stochastic optimization framework considering multiple realizations was used to account for geological uncertainty. They derived the analytical form of the gradient and Hessian matrix of the objective function and implemented the optimization process under operating and facility constraints using a sequential quadratic programming method.

[0025] Other streamline-based production optimization methods include using other physical quantities derived from flow field diagnostic data to optimize and adjust flow distribution. For example, a waterflood optimization strategy based on streamline simulation to maximize volumetric sweep efficiency was proposed for waterflood management of injection cycles and inefficient displacement in oilfields. They defined the oil and gas F-Φ curves and oil and gas Lorenz coefficients (L C-HC ) by analogy with the flow-reservoir capacity curves (F-C curves) in reservoir engineering, which can be derived from streamline simulation results and used for waterflood optimization. L C-HC can be used as a unique measure of flow or dynamic heterogeneity, and minimizing L C-HC can maximize volumetric sweep efficiency. Based on this principle, they evaluated the sensitivity of L C-HC to changes in operating conditions in a Design of Experiments (DoE) study and used a response surface method to describe L C-HC as a function of these operating conditions, and then selected the operating conditions that minimize L C-HCThe operating conditions at which the minimum (or maximum sweep efficiency) is achieved. Some researchers have optimized the number of fracturing stages in tight gas reservoirs by investigating the change in drainage volume derived from gas streamlines and taking the drainage volume as a function of the number of fracturing stages. Some researchers have investigated the evolution of drainage volume from the well derived from gas streamlines and optimized the well location in naturally fractured tight gas reservoirs. Some researchers have proposed a simple and practical waterflood rate optimization workflow using streamline-based displacement efficiency maps to address the problem that conventional waterflood management methods are difficult to achieve reasonable optimization at the field scale. Displacement efficiency maps can show the flow rate distribution and the time of flight distribution of production wells. The optimization idea is to balance the average time of flight (TOF) between production wells on a sub-area basis, thereby maximizing the waterflood sweep efficiency; analytical calculations are used to determine the weighting factors for injection and production rates to minimize the variance in TOF between production wells. Some researchers have combined the flow field diagnostic information with a fuzzy logic system and proposed an adaptive method for optimizing water injection strategies in oilfields.

[0026] In addition, for the production optimization problem under the conditions of complex reservoir model (such as multiple media, faults, fractures and other complex geological features), complex fluid (such as water-gas alternation, chemical flooding, steam flooding and other processes), considering geological uncertainty, streamline simulation is not mature yet, and still needs continuous research. It is quite representative that some researchers based on Embedded Discrete Fracture Model (EDFM) proposed a robust streamline tracking framework using the boundary layer method, and introduced the application of this framework in flow visualization, flow diagnosis and flow allocation optimization. This method can reflect the flow between matrix and fracture or between fractures, and use streamlines to establish the workflow of water drive flow allocation optimization in natural fractured reservoirs. Some researchers proposed a streamline conversion method based on flow exchange calculation probability and a dual porosity dual permeability (DPDK) reservoir streamline tracking algorithm to solve the problem that the current streamline method is difficult to describe the interaction between matrix and fracture in DPDK medium, which can be applied to streamline tracking and water drive flow allocation optimization in DPDK model. Some researchers used streamline flight time and principal component analysis to extract the flow characteristics of the geological model, used K-means algorithm to generate subsets of uncertain models, and then proposed a flow allocation optimization workflow considering the influence of uncertainty, realizing the combination of streamline simulation and machine learning technology. Some researchers considered geological uncertainty, and on the basis of the streamline-based flow allocation optimization method for a single numerical model, used the optimal interpolation method to update the "joint well pair flow multiplier" by synthesizing the production data of multiple geological realizations, thereby realizing the robust optimization of injection-production flow rate. However, the intelligent flow field optimization theory combining streamline simulation and machine learning algorithm is still in the initial stage of research and development, and still needs to be vigorously researched.

[0027] The research found that although the current development of mature commercial streamline simulator can carry out three-dimensional three-phase black oil or component simulation, the complex model considering phase change, gravity and capillary force, strong compressibility fluid and other effects may reduce the speed and stability of streamline simulation, so streamline simulation is mostly applied to the study of oil-water two-phase flow model. At present, the production optimization based on streamline simulation mainly optimizes the flow distribution according to the flow field diagnostic indicators such as injection efficiency and flight time. The former is to calculate the injection efficiency of single well or well pair according to the flow distribution coefficient, set the corresponding objective function (such as maximum net present value, maximum cumulative oil production, balanced injection efficiency, etc.), and sequentially solve the optimal water injection rate and oil production according to certain flow updating criteria. This method is suitable for the middle and late stages of water drive development and is beneficial to stable oil control and water control. The latter is to make the injection water arrival time (see water time) of each production well as close as possible, so as to realize the effect of approximate balanced displacement and improve the sweep efficiency. However, this method is often more effective in the early stage of development. In addition, some researchers combine streamline simulation with mathematical methods (such as fuzzy logic, optimization algorithm, other derived information or self-defined diagnostic indicators) to optimize injection and production flow. The reliability of these methods depends on the accuracy of the historical fitting of the early streamline model to a certain extent. However, the production optimization problem at the oilfield scale usually involves very complex reservoir models, production and facility constraints, and a large number of unknown factors, which makes the production optimization method based on streamline simulation (flow field diagnosis) still face the following challenges: (1) the streamline tracking technology in complex oil and gas reservoirs (multiple medium reservoirs such as fractured reservoirs, fracture-cave reservoirs and triple medium reservoirs) is still not perfect; (2) there are uncertain parameters in the flow updating criteria, which may also need to be optimized, and the optimization algorithm can be further coupled in the flow field diagnosis optimization; (3) the optimization process still cannot effectively consider the uncertainty of the geological model, and the reliability of the result depends on the accuracy of the historical fitting of the early model to a certain extent, which can be partially solved by combining machine learning and optimization theory, and the intelligent closed-loop production optimization method combined with automatic history matching needs to be developed; (4) at present, there is no mature flow field diagnosis method based on the derived streamline information of streamline simulator or finite difference simulator and production optimization platform based on this, and the streamline simulation optimization method has not been widely applied; (5) the collaborative mechanism with simplified physical models has not been systematically studied.

[0028] The foregoing reviews the existing reservoir production optimization methods, and based on the research progress and advantages and disadvantages analysis of various production optimization methods, the existing methods have defects such as time-consuming and laborious, high cost, low efficiency, low reliability, insufficient accuracy, low intelligence, and difficulty in balancing the contradiction between real-time regulation and long-term prediction demand. SUMMARY

[0029] In a first aspect, embodiments of the present application provide an oil reservoir comprehensive production optimization method, which is used to meet injection-production adjustment requirements of oil reservoir management in each development stage and improve efficiency and accuracy of oil reservoir production optimization, and comprises the following steps:

[0030] obtaining production history data of a target oil reservoir;

[0031] performing regional level understanding, single well level understanding and production optimization model understanding according to the production history data, and updating oil reservoir dynamic understanding;

[0032] based on the production history data, performing cyclic history matching on model parameters of multiple production optimization models until optimal model parameters are obtained, the multiple production optimization models comprising a first type of prediction model and a second type of prediction model, the history matching on the model parameters of each production optimization model, and the inter-well connectivity relationship calibration of the multiple production optimization models, so that the first type of prediction model provides an initial estimation of the inter-well connectivity relationship for the second type of prediction model, and the second type of prediction model calibrates the model parameters of the first type of prediction model;

[0033] determining constraint conditions corresponding to each production optimization model according to the updated oil reservoir dynamic understanding, solving a target function corresponding to each production optimization model, obtaining optimal decision variables corresponding to each production optimization model, and obtaining a development regulation scheme corresponding to each production optimization model;

[0034] obtaining an oil reservoir comprehensive development regulation scheme according to the development regulation scheme corresponding to each production optimization model, the oil reservoir comprehensive development regulation scheme being used to guide development operation of the target oil reservoir;

[0035] updating the oil reservoir dynamic understanding according to a regulation effect obtained by real-time monitoring of the development operation.

[0036] In a second aspect, embodiments of the present application further provide an oil reservoir comprehensive production optimization device, which is used to meet injection-production adjustment requirements of oil reservoir management in each development stage and improve efficiency and accuracy of oil reservoir production optimization, and comprises the following steps:

[0037] a production history data obtaining module, configured to obtain production history data of a target oil reservoir;

[0038] an oil reservoir dynamic understanding updating module, configured to perform regional level understanding, single well level understanding and production optimization model understanding according to the production history data, and update oil reservoir dynamic understanding;

[0039] The history matching module is configured to perform cyclic history matching on model parameters of a plurality of production optimization models, including a first type of prediction model and a second type of prediction model, based on production history data until optimal model parameters are obtained, by performing history matching on the model parameters of each production optimization model and calibrating interwell connectivity for the plurality of production optimization models, such that the first type of prediction model provides an initial estimation of interwell connectivity for the second type of prediction model and the second type of prediction model calibrates the model parameters of the first type of prediction model.

[0040] The comprehensive production optimization module is implemented to determine constraint conditions corresponding to each production optimization model according to the updated reservoir dynamic understanding, solve an objective function corresponding to each production optimization model, obtain optimal decision variables corresponding to each production optimization model, and obtain a development control scheme corresponding to each production optimization model; and a comprehensive development control scheme of the reservoir is obtained according to the development control schemes corresponding to each production optimization model, and the comprehensive development control scheme of the reservoir is used to guide development operations of the target reservoir.

[0041] The adjustment effect feedback and updating module is configured to update the reservoir dynamic understanding according to a control effect obtained by real-time monitoring of the development operations.

[0042] In a third aspect, an embodiment of the present application further provides a computer device, which comprises a memory, a processor, and a computer program stored in the memory and executable on the processor, and the processor implements the reservoir comprehensive production optimization method when executing the computer program.

[0043] In a fourth aspect, an embodiment of the present application further provides a computer readable storage medium, which stores a computer program, and the computer program is executable on a processor to implement the reservoir comprehensive production optimization method.

[0044] In a fifth aspect, an embodiment of the present application further provides a computer program product, which comprises a computer program, and the computer program is executable on a processor to implement the reservoir comprehensive production optimization method.

[0045] In the embodiments of the present application, the advantages of reservoir dynamic understanding, the first type of prediction model and the second type of prediction model are fused, the production optimization based on the reservoir dynamic understanding (semi-empirical) is not dependent on the physical prediction model or the numerical prediction model, mainly uses the advantage of convenient operation to provide optimization suggestions in a very short time or high frequency (such as daily); the production optimization based on the first type of prediction model mainly uses the speed efficiency advantage to provide short-term optimization suggestions, and the production optimization based on the second type of prediction model mainly uses the accuracy (reliability) advantage to provide long-term decision reference; the reservoir dynamic understanding summarized from the main control factor analysis, the periodic dynamic analysis, the node system analysis and the like is taken as a constraint to make the result of the comprehensive production optimization more reasonable; in turn, the production optimization results of the first type of prediction model and the second type of prediction model together with the periodic dynamic analysis understanding form new reservoir dynamic understanding; according to the comprehensive production optimization concept proposed in the present application, various production optimization methods meeting different optimization requirements can be gradually developed and formed, the intelligent reservoir closed-loop production optimization method with the advantages of complementary, real-time and long-term consideration, physical meaning and data-driven, and significant optimization efficiency and accuracy improvement. BRIEF DESCRIPTION OF DRAWINGS

[0046] In order to more clearly illustrate the technical solutions in the embodiments of the present application or the prior art, the drawings needed to be used in the embodiments or the prior art description will be briefly introduced. Obviously, the drawings in the following description are only some embodiments of the present application, and other drawings can be obtained by those skilled in the art without any creative effort. In the drawings:

[0047] Fig. 1 is a flow chart of the comprehensive production optimization method of the reservoir in the embodiments of the present application;

[0048] Fig. 2 is an intelligent closed-loop reservoir management framework in the embodiments of the present application;

[0049] Fig. 3 is a classification of the reservoir production optimization model in the embodiments of the present application;

[0050] Fig. 4 is a comprehensive production optimization workflow in the embodiments of the present application;

[0051] Fig. 5 is a reservoir dynamic understanding updating process in the embodiments of the present application;

[0052] Fig. 6 is an injection-production optimization process based on the reservoir dynamic understanding updating in the embodiments of the present application;

[0053] Fig. 7 is a history fitting process of the time-varying resistance capacitance model TV-CRM in the embodiments of the present application;

[0054] Fig. 8 is a production optimization process of the time-varying resistance capacitance model TV-CRM in the embodiments of the present application;

[0055] Figure 9 is the interwell connectivity relationship of TV-CRM inversion in the embodiment of the present application;

[0056] Figure 10 is the liquid production of TV-CRM history matching in the embodiment of the present application;

[0057] Figure 11 is the pore volume of each production well control volume of TV-CRM inversion in the embodiment of the present application;

[0058] Figure 12 is the productivity index of each production well of TV-CRM inversion in the embodiment of the present application;

[0059] Figure 13 is the initial water saturation of each production well control volume of TV-CRM inversion in the embodiment of the present application;

[0060] Figure 14 is the water cut history matching of TV-CRM in the embodiment of the present application;

[0061] Figure 15 is the development control scheme and predicted production given by TV-CRM in the embodiment of the present application;

[0062] Figure 16 is a schematic diagram of a plane connection unit in the embodiment of the present application;

[0063] Figure 17 is a schematic diagram of a multi-layer connection unit in the embodiment of the present application;

[0064] Figure 18 is a schematic diagram of initial condition approximation processing in the process of one-dimensional saturation tracking in the embodiment of the present application;

[0065] Figure 19 is the initial information of the connection unit simulation model corresponding to Example 2 in the embodiment of the present application;

[0066] Figure 20 is the fitting effect of part of the oil production and water cut of single well in Example 2 in the embodiment of the present application;

[0067] Figure 21 is the characteristic parameter of the connection unit simulation model after history matching in the embodiment of the present application;

[0068] Figure 22 is the injection distribution coefficient of the connection unit simulation model inversion in the embodiment of the present application;

[0069] Figure 23 is the injection volume optimization result of part of the injection wells in the embodiment of the present application;

[0070] Figure 24 is the change of liquid production of part of the production wells before and after optimization in the embodiment of the present application;

[0071] Figure 25 is the general process of pure data-driven proxy model optimization in the embodiment of the present application;

[0072] Figure 26 is a schematic diagram of the pore volume of the kth flow line on the injection-production well pair (i, j) in the embodiment of the present application;

[0073] Figure 27 shows the curves of the flow update weight W as a function of the normalized flow field diagnostic index β under different weight indices in the unbounded weight criterion of this application (β is taken as β). ave (For example, 0.5);

[0074] Figure 28A shows the sensitivity of the relationship curve between W and β in the bounded weight criterion of this application to α (R = 0.85, W lim,min =-1, W lim,max =1);

[0075] Figure 28B) shows the sensitivity of the relationship curve between W and β in the bounded weight criterion of this application to R (α = 2, W lim,min =-1, W lim,max =1);

[0076] Figure 28C) shows the relationship curves between W and β in the bounded weight criterion of this application. lim,min Sensitivity (α=2, R=0.85, W) lim,max =1);

[0077] Figure 28D) shows the relationship curves between W and β in the bounded weight criterion of this application. lim,max Sensitivity (α=2, R=0.85, W) lim,min =-1);

[0078] Figure 29 shows the optimization process of the streamline simulation model for a single time step in the embodiments of this application;

[0079] Figure 30 shows the development and control scheme of the streamline simulation model based on well-pair injection efficiency in the embodiments of this application;

[0080] Figure 31 shows the development and control scheme of the streamline simulation model based on the displacement efficiency of movable residual oil in the well pair in the embodiments of this application;

[0081] Figure 32 shows the general optimization process based on the numerical simulator proxy model in the embodiments of this application;

[0082] Figure 33 is a flowchart illustrating the proxy optimization algorithm in an embodiment of this application;

[0083] Figure 34 shows the permeability distribution of the heterogeneous reservoir model with 5 injections and 4 productions in the embodiments of this application;

[0084] Figure 35 shows the injection and extraction history of Embodiment 3 in this application;

[0085] Figure 36 shows the development control scheme obtained by proxy optimization in Embodiment 3 of this application;

[0086] Figure 37 shows the change of net present value with the number of numerical simulations during the agent optimization process in Example 3 of this application;

[0087] Figure 38 is an optimization process of the streamline simulation proxy model in the embodiment of the present application;

[0088] Figure 39 is a basic element of the optimization process of the streamline simulation proxy model in the embodiment of the present application;

[0089] Figure 40 is an optimization process of the streamline simulation proxy model based on the injection efficiency of the well pair in the embodiment of the present application (the proxy optimization decision variable is a, R, W lim,min = -1, W lim,max = 1) ;

[0090] Figure 41 is a development control scheme obtained by optimizing the streamline simulation proxy model based on the injection efficiency of the well pair in the embodiment of the present application;

[0091] Figure 42 is an optimization calculation process of the streamline simulation proxy model based on the movable residual oil displacement efficiency of the well pair in the embodiment of the present application (the decision variable is a, R, W lim,min , W lim,max ) ;

[0092] Figure 43 is a development control scheme obtained by optimizing the streamline simulation proxy model based on the movable residual oil displacement efficiency of the well pair in the embodiment of the present application;

[0093] Figure 44 is a coupling interaction relationship of the intelligent comprehensive production optimization method of the reservoir in the embodiment of the present application;

[0094] Figure 45A) is an interwell connectivity coefficient inverted by the TV-CRM in the embodiment of the present application;

[0095] Figure 45B) is an injection flow rate distribution coefficient derived from the streamline simulation at the end of the history matching period in the embodiment of the present application;

[0096] Figure 46 is an optimization process of the streamline simulation proxy model based on the injection efficiency of the well pair in the embodiment of the present application (the decision variable is a, R, W lim,min = -1, W lim,max = 1) ;

[0097] Figure 47 is an optimization process of the streamline simulation proxy model based on the movable residual oil displacement efficiency of the well pair in the embodiment of the present application (the decision variable is a, R, W lim,min = -1, W lim,max = 1) ;

[0098] Figure 48 is an injection-production setting of the baseline scheme of Example 4 in the embodiment of the present application;

[0099] Figure 49 is a development control scheme of Example 4 obtained by optimizing the streamline simulation proxy model based on the injection efficiency of the well pair in the embodiment of the present application;

[0100] Fig. 50 is a development regulation scheme of Example 4 obtained by optimizing the flow line simulation proxy model based on well pairs movable residual oil displacement efficiency in the embodiment of the present application;

[0101] Fig. 51 is a schematic diagram of the oil reservoir comprehensive production optimization device in the embodiment of the present application;

[0102] Fig. 52 is a schematic diagram of the computer device in the embodiment of the present application. DETAILED DESCRIPTION

[0103] In order to make the purpose, technical scheme and advantages of the embodiments of the present application more clear, the embodiments of the present application are further described in detail below with reference to the drawings. Herein, the illustrative embodiments of the present application and their descriptions are used to explain the present application, but not as a limitation to the present application.

[0104] Improving oil recovery and economic benefits is the ultimate goal of oilfield development, and closed-loop reservoir management is an effective way to achieve this goal. However, the single production optimization method in reservoir management currently has a series of problems such as being empirical, poor transplantability, low efficiency, insufficient precision, high computational cost, difficulty in solving constraint optimization problems, and low intelligence level. In order to meet the injection-production adjustment needs of reservoir management at each development stage and improve the efficiency and precision of reservoir production optimization, it is necessary to improve and comprehensively utilize various optimization models or production optimization methods, to develop efficient solving algorithms for time-consuming constraint optimization problems, to propose the comprehensive production optimization concept, workflow and implementation method of cross-fusion of various production optimization methods, and to apply the method to improve the intelligent level and economic benefits of reservoir management.

[0105] The present application is directed to the many shortcomings of existing reservoir production optimization methods (for example, the oil reservoir engineering method is strong in experience, poor in portability, and low in precision, the traditional numerical simulation method has high data requirements, high computational cost, time-consuming and laborious, and it is difficult to obtain an optimal solution, the simplified physical model (such as CRM, INSIM and FNM, etc.) has strong multi-solution, weak long-term prediction ability, low reliability, limited application range, the machine learning agent modeling optimization method still takes a long time in the agent construction stage, and the optimization efficiency is low for high-dimensional problems, and the streamline simulation optimization method needs to rely on the geological model history fitting result, and the flow field diagnosis method suitable for the whole life period of the reservoir has not been effectively defined, and the flow updating criterion has not been proposed. From the basic definition of production optimization, the main principles and advantages and disadvantages of several important production optimization (models) such as production optimization based on reservoir dynamic understanding, time-varying capacitance-resistance model (TV-CRM), linkage unit simulation model (LUSM), numerical simulator agent model, streamline simulation model, and streamline simulation agent model are briefly discussed, the implementation process is improved and explained, and future research directions are proposed.

[0106] In particular, the scheme proposed by the embodiments of the present application intends to combine production optimization based on simplified physical models and production optimization based on numerical reservoir models under the constraint of reservoir dynamic understanding in the framework of intelligent closed-loop reservoir management. The production optimization based on reservoir dynamic understanding (not dependent on physical or numerical prediction models) mainly utilizes its convenient operation to provide optimization suggestions in a very short time or at a high frequency (such as daily), the production optimization based on simplified physical models mainly utilizes its speed efficiency to provide short-term optimization suggestions, and the production optimization based on numerical reservoir models (surrogate models and streamline simulation) mainly utilizes its accuracy (reliability) to provide long-term decision reference. The reservoir dynamic understanding (and development technical strategies, injection-production adjustment strategies, etc.) summarized by the main factor analysis, periodic dynamic analysis, and node system analysis is used as a constraint to make the result of the comprehensive production optimization more reasonable. In return, the results of the production optimization based on simplified physical models and the production optimization based on numerical reservoir models form new reservoir dynamic understanding together with the periodic dynamic analysis understanding. Thus, an intelligent reservoir closed-loop production optimization method is gradually developed and formed, which is based on specific optimization requirements, and realizes the advantages of various production optimization methods, such as complementary advantages, real-time and long-term consideration, and physical meaning and data-driven.

[0107] The present scheme will be described in detail below.

[0108] FIG. 1 is a flowchart of the reservoir comprehensive production optimization method in the embodiments of the present application, which includes:

[0109] Step 101, obtaining production history data of a target reservoir;

[0110] Step 102, performing regional horizontal understanding, single-well horizontal understanding, and production optimization model understanding according to the production history data, and updating the reservoir dynamic understanding;

[0111] Step 103, based on the production history data, performing cyclic history matching on model parameters of multiple production optimization models until optimal model parameters are obtained, the multiple production optimization models including a first type of prediction model and a second type of prediction model: performing history matching on the model parameters of each production optimization model, and calibrating the inter-well connectivity relationship of the multiple production optimization models, so that the first type of prediction model provides an initial estimation of the inter-well connectivity relationship for the second type of prediction model, and the second type of prediction model calibrates the model parameters of the first type of prediction model;

[0112] Step 104, according to the updated reservoir dynamic understanding, determine the constraint condition corresponding to each production optimization model, solve the objective function corresponding to each production optimization model, obtain the optimal decision variable corresponding to each production optimization model, and obtain the development control scheme corresponding to each production optimization model;

[0113] Step 105, according to the development control scheme corresponding to each production optimization model, obtain the comprehensive development control scheme of the reservoir, and the comprehensive development control scheme of the reservoir is used to guide the development operation of the target reservoir;

[0114] Step 106, according to the control effect obtained by real-time monitoring of the development operation, update the reservoir dynamic understanding.

[0115] Figure 2 is an intelligent closed-loop reservoir management (closed-loop production optimization) framework corresponding to Figure 1 in the embodiments of the present application, which comprehensively uses a first type of prediction model (simplified physical model) and a second type of prediction model (pure data-driven proxy model or non-proxy model in the reservoir numerical model) under the constraint of reservoir dynamic understanding. When the development control scheme is implemented on the target reservoir or region, fluid migration occurs in the reservoir and causes oil and water redistribution. Based on the reservoir dynamic monitoring data and production data, the reservoir dynamic understanding is summarized by main factor analysis, periodic dynamic analysis, node system analysis, etc. The development control scheme formulated according to the reservoir dynamic understanding itself can be used as a production optimization method independent of the prediction model, thereby proposing a development control scheme for the selected problem well in a semi-empirical manner. This production optimization method is more suitable for very short time or high frequency adjustments (such as daily adjustments), while other optimization methods are difficult to complete the history matching process and give the development control scheme within a few minutes. On the one hand, based on daily production data (such as weekly, monthly data), the model parameters of the simplified physical model can be determined by the history matching algorithm, such as the connectivity coefficients and control volumes of TV-CRM, the initial connection unit conductivity and initial control volume of LUSM, etc. The relative sizes of these connectivity coefficients and conductivities can provide a reference for the history matching of the reservoir numerical model. The simplified physical model after history matching can take advantage of its rapid calculation to carry out real-time injection-production optimization through optimization algorithms, giving a short-term development control scheme. On the other hand, based on the reservoir dynamic understanding, it is determined whether the previous geological model needs to be corrected or re-geological modeling is needed, and then the automatic history matching technology is used to correct and update the reservoir numerical model to make it consistent with the production response of the new production data. Automatic history matching is generally completed by machine learning proxies (such as convolutional neural networks combined with long short-term memory networks) and history matching algorithms (such as ensemble Kalman filter EnKF and multiple data assimilation ensemble smoothing algorithm ES-MDA). The non-proxy model in the reservoir numerical model after history matching can be used for the prediction of future production dynamics and as the basis for numerical simulator proxy model optimization (here, proxy optimization mainly refers to numerical simulator-based proxy optimization; pure data-driven proxy optimization is directly based on production dynamic data) and streamline simulation model optimization methods. Numerical simulator proxy model optimization, streamline simulation model optimization, and streamline simulation proxy model optimization combining the two can provide medium and long-term development control schemes for reservoir development. The short-term and medium and long-term development control schemes obtained by the simplified physical model and the numerical simulator proxy model under the constraint of reservoir dynamic understanding are collectively used as a reference for actual production operations, guiding production adjustments of different frequency requirements in the field. At the same time, the production optimization results based on the simplified physical model and the reservoir numerical model, together with the periodic dynamic analysis, form new reservoir dynamic understanding.After the development regulation is implemented, dynamic monitoring is continued, thereby entering a new round of production optimization cycle.

[0116] From the perspective of generalized optimization problem, reservoir production optimization refers to the process of applying optimization theory to solve various optimization problems in reservoir development and management. Like general optimization problems, it specifically refers to using a certain physical or numerical model to evaluate the production optimization objective function, under the condition of meeting various complex constraints of the reservoir system, using optimization algorithms to achieve the optimal value of the given objective function through automatic or heuristic optimization strategies, obtaining the optimal decision variable, and thus giving the best development regulation scheme. The core elements of production optimization include: construction of optimization problem (decision variable, objective function, constraint condition), optimization model and optimization algorithm. In this embodiment, the decision variable, objective function and constraint condition depend on the optimization problem constructed to meet the production demand, such as the coordinates of well location in well location optimization, the flow or bottom hole pressure of injection and production wells in injection-production optimization, the layers that need to be opened or closed in layer optimization, and the combination of them in joint optimization; the objective function of reservoir production optimization is generally the maximum cumulative oil production or the highest net present value, or multiple objectives such as the maximum cumulative oil production, the minimum cumulative water production, the lowest water cut, the highest net present value, etc. The constraint conditions generally include single well flow, pressure, water cut constraints, well group injection capacity, liquid production capacity, produced water treatment capacity constraints, and block or oilfield level total flow, pressure, water cut constraints.

[0117] The production optimization model refers to a physical model, a numerical model, a strategy or a means relied on for evaluating the objective function in the production optimization process and then guiding the search direction of the decision variable, such as a simplified physical model (CRM, INSIM, FNM, etc.), a finite difference numerical simulation model (i.e., a finite difference model), a streamline simulation model, a streamline simulation proxy model, etc. Of course, the optimization model can also be a heuristic optimization strategy, such as a periodic water injection strategy for adjusting injection and production flow rates according to pressure and water cut indicators, in which case it does not evaluate the objective function but only determines the trend of the objective function. FIG. 3 classifies the production optimization model into a prediction model class and a non-prediction model class according to whether the prediction model is needed to evaluate the objective function. The non-prediction model class mainly refers to a production optimization method based on reservoir dynamic understanding. The prediction model class is further classified into a simplified physical model, a reservoir numerical model and a pure data-driven proxy model according to whether a geological model needs to be established. In this embodiment, both the simplified physical model and the pure data-driven proxy do not need fine geological modeling, and thus have fast operation speed and high efficiency. The production optimization model relying on reservoir numerical simulation can be classified into a proxy model and a non-proxy model. The non-proxy model includes a finite difference model and a streamline simulation model. The proxy model includes a numerical simulator proxy model and a streamline simulator proxy model. They need to perform fine geological modeling on the reservoir system according to dynamic and static data, and perform history matching using production data to calibrate the numerical model. Such models have high data requirements and high computational cost, but are more reliable and have stronger long-term prediction capability. The input-output response data generated by the numerical simulator is used for training of the numerical simulator proxy model, and the training and optimization processes are alternately performed. The traditional numerical simulator generally uses a finite difference model, while the streamline simulation model can be converted from the finite difference model, can be generated by a special streamline simulator, and can also be researched by streamline tracing and flow field diagnosis according to the simulation results output by the finite difference model.

[0118] The production optimization algorithm refers to an algorithm used for solving the optimization problem constructed by the production optimization. According to whether the gradient of the objective function with respect to the decision variable is involved, it can be roughly classified into a gradient class, an approximate gradient class and a non-gradient class. The optimization algorithm often also involves a constraint handling method (linear / nonlinear equation and linear / nonlinear inequality) and a framework for algorithm execution, and sometimes the optimization algorithm, the constraint handling method, the solving framework or the strategy are collectively referred to as an optimization method or algorithm.

[0119] The reservoir production optimization problem is generally divided into layer series optimization, well pattern optimization, well location optimization, injection-production parameter optimization, etc. according to different optimization objects (decision variables), or combination optimization thereof. In this application, injection-production optimization is taken as an example. Reservoir injection-production optimization refers to a process of adjusting injection-production parameters (injection rate / injection pressure, liquid production rate / production well BHP (Bottom Hole Pressure), pump frequency of electric pump well, choke, etc.) of injection and production wells to maximize or minimize the objective function. The decision variable is generally set as the injection rate of the injection well and the liquid production rate of the production well. The constraint condition should be the reasonable development technical strategy limit, and the optimization objective is generally to maximize oil production, minimize water cut, maximize pressure recovery, maximize net present value, etc. First, several methods (or models) for solving the reservoir production optimization problem are introduced, and then a comprehensive production optimization method is proposed, and the specific connotation and implementation steps are explained. In the intelligent closed-loop reservoir management framework shown in FIG. 2, with the aid of three main means of reservoir dynamic understanding, the first type of prediction model and the second type of prediction model, the internal relationship among them is established through the identification and calibration of the interwell connectivity relationship, and an intelligent comprehensive production optimization workflow of the reservoir is proposed, which is a benign coupling interaction between each optimization model.

[0120] FIG. 4 is a detailed flowchart of the reservoir comprehensive production optimization method corresponding to FIG. 1 in the embodiments of the present application. As described above, the production optimization based on reservoir dynamic understanding belongs to the non-prediction model type of production optimization method. This type of method reflects the experience-based expression derived from the reservoir understanding of researchers in various disciplines. Therefore, in addition to being used for high-frequency (such as daily) adjustment, it is more important as a constraint condition for other production optimization methods, so that the results of comprehensive production optimization reflect the actual reservoir and are more operable. For the prediction model type of production optimization method, the physical quantities obtained through historical fitting or model calibration can be related through interwell connectivity. Therefore, the comprehensive production optimization workflow shown in FIG. 4 is designed, which embodies the benign interaction between the simplified physical model and the numerical optimization model, and the interaction between the reservoir dynamic understanding and them.

[0121] In step 101, production history data of a target reservoir is obtained.

[0122] The reservoir production history data includes all the historical data required for production optimization of the target reservoir. Taking injection-production optimization as an example, it includes injection-production flow rate data (water injection rate, oil production rate, water production rate, gas production rate, etc.), pressure data (wellhead pressure, bottom hole flowing pressure, static pressure, average formation pressure, etc. of single wells), and physical parameters of reservoir rock and fluid (compressibility of rock and fluid, change of oil, gas and water viscosity and volume coefficient with pressure, oil, gas and water relative permeability curve and capillary pressure curve, etc.). It should be noted that the production history data may, for example, be data obtained from the target reservoir based on sensors in the past.

[0123] In step 102, according to production history data, regional level understanding, single well level understanding and production optimization model understanding are performed, and reservoir dynamic understanding is updated;

[0124] Reservoir dynamic understanding refers to the understanding of reservoir geological characteristics and development rules formed by different researchers based on different discipline knowledge through analysis of seismic, logging, various tests and production data. The "dynamic" understanding refers to the object of the understanding and knowledge, which is the production dynamic of the reservoir (including the region, well group and single well), and refers to the fact that the understanding is changed with time and gradually deepened with the development process.

[0125] FIG. 5 is a reservoir dynamic understanding updating flow in the embodiment of the application. The formation of reservoir dynamic understanding has three main sources, i.e. regional level understanding obtained from main control factor analysis, single well level understanding obtained from periodic dynamic analysis and node system analysis, and theoretical model understanding obtained from production optimization model.

[0126] In an embodiment, according to production history data, regional level understanding, single well level understanding and production optimization model understanding are performed, and reservoir dynamic understanding is updated, including:

[0127] In the reservoir development process, the development stage and the regional development characteristics are analyzed from multiple development aspects to analyze multiple main control factors affecting the development effect of the reservoir and the region (even well group / single well) and producing the observed production dynamic (production, pressure, water cut change rule), and the multiple development aspects include at least one of geology, reservoir development and development engineering;

[0128] According to different main control factors, the reservoir is classified;

[0129] For the reservoir classification conforming to the overall law of the region, regional level understanding is obtained;

[0130] For the reservoir classification not conforming to the overall law of the region (i.e. single well to be managed separately, including well group), single well level understanding is updated by periodic dynamic analysis and node system analysis (at this time, the "one well one strategy" method can be used);

[0131] The regional level understanding is corrected through periodic reservoir dynamic analysis, and the single well level understanding is verified through the regional level understanding;

[0132] The regional or single well level understanding is updated through theoretical model understanding, and the theoretical model understanding is constrained through the single well level understanding;

[0133] According to regional level understanding, single well level understanding and production optimization model understanding, reservoir dynamic understanding is updated;

[0134] According to the updated reservoir dynamic understanding, a reasonable development technical strategy is determined;

[0135] With the updated reservoir dynamic understanding as the constraint condition, including:

[0136] With the reasonable development technical strategy as the constraint condition.

[0137] Among them, the main control factor analysis embodies the reservoir management logic of first grasping the main contradiction and then making the whole into parts.

[0138] The main control factor analysis is to summarize the medium and long-term influencing factors of the reservoir / region (partial integrity) from the perspective of the development stage and regional development characteristics. In fact, it can also be regarded as a kind of phased dynamic analysis, but the single well level understanding still needs to rely on the short-term periodic reservoir dynamic analysis, i.e. periodic reservoir dynamic analysis. For example, through real-time monitoring of the abnormality of daily and weekly production indicators (liquid volume, pressure, water cut) and the feedback of work system changes to screen the wells that may need to be adjusted or assist in determining the direction of future injection-production adjustment. The role of periodic reservoir dynamic analysis is to continuously verify and correct the regional level understanding under the framework of main control factor analysis, and to develop individual development technical strategies for those "special wells" that do not conform to the overall law of the region, so that each single well (which can be a well group) can extract the reservoir understanding basis that constitutes the reasonable development technical strategy to guide its rapid adjustment.

[0139] Although the update of reservoir recognition mainly comes from the regional level recognition obtained by the analysis of main controlling factors and the single well level recognition obtained by the periodic reservoir dynamic analysis, in the present embodiment, the single well level recognition can be referred to as the individual recognition at the single well level, but in order to avoid the experience of the researchers to some extent, the optimization results obtained by the node system analysis and other production optimization models relying on prediction models (such as the simplified physical models focusing on short-term optimization: CRM, INSIM, FNM, etc.; the proxy models in the numerical reservoir model focusing on medium and long-term optimization: numerical simulator proxy model, pure data-driven proxy model; and the non-proxy models in the numerical reservoir model: finite difference model, streamline simulation model, etc.) should be referred to as supplementary recognition to assist the dynamic update of reservoir recognition. The node system analysis (taking the production well as an example) refers to selecting a node in the migration path of the selected fluid from the reservoir to the wellhead (or the gathering station) as the solution node, analyzing the inflow dynamic relationship (IPR) from the reservoir to the solution node and the outflow dynamic relationship (OPR) from the solution node to the wellhead, determining the production coordination point through the intersection of the IPR curve and the OPR curve, thereby obtaining the optimal pressure and flow rate setting, determining the reasonable limits of the production index through the change law of the production coordination point with the sensitivity parameters, thereby assisting the update of the reservoir dynamic recognition at the single well level. For example, the connectivity coefficients obtained by the CRM historical fitting, the inter-node connection conductivities derived by the INSIM model, and the injection-production flow distribution coefficients obtained by the streamline simulation can be used to more accurately identify the injection direction of the injection well and the liquid source of the production well, estimate the injection-production ratio at the single well level, and calculate the diagnostic indexes such as injection-production efficiency in a manner similar to the flow field diagnosis of the streamline simulation, evaluate the pros and cons of injection-production effect, and then correct the reservoir recognition and assist in the formulation of reasonable development technical strategies.

[0140] The reasonable development technical strategy refers to the reasonable value or reasonable range of the content specified by the oilfield development or development adjustment, such as development strategy, development mode, layer division, well pattern and spacing, injection mode, injection-production parameter, and measure selection, which is applicable to the target oilfield (or reservoir / well group / single well, etc. at the regional level) under the current conditions. Similarly, the reasonable water injection development technical strategy refers to the reasonable range of the development index (and adjustment strategy) set at the target reservoir, region, well group, single well, etc. level during the water drive development, such as the reasonable limits of the injection-production flow rate, injection-production ratio, pressure maintenance level, water cut, etc. (and the reasonable adjustment direction when the corresponding index violates the constraint).

[0141] The reservoir dynamic understanding determines the determination of the reasonable development technical strategy, and the reasonable development technical strategy constitutes the basis of the subsequent reservoir state diagnosis and development regulation scheme. The purpose of the reservoir state diagnosis is to screen out potential "problem wells" according to the reasonable development technical strategy, and then the injection-production adjustment strategy gives the development regulation scheme for the problem wells according to the established procedure or rule. Finally, through the real-time monitoring and effect feedback of the adjustment operation, the reservoir dynamic understanding and the injection-production adjustment strategy are updated. In this way, an injection-production optimization closed loop based on the reservoir dynamic understanding update is formed.

[0142] It can be seen that in the intelligent comprehensive production optimization process, the reasonable development technical strategy needs to be used as a constraint condition, but as described before, the reservoir dynamic understanding update itself is a production optimization method, therefore, the embodiment of the present application further proposes a production optimization process based on the reservoir dynamic understanding update, in an embodiment, the method further comprises:

[0143] According to the reasonable development technical strategy, the reservoir state diagnosis is performed to determine the problem wells;

[0144] The development regulation scheme is generated for the problem wells, and the development regulation scheme is used for the adjustment operation of the problem wells;

[0145] According to the real-time monitoring and effect feedback of the adjustment operation, the reservoir dynamic understanding and the development regulation scheme are updated.

[0146] Taking one of the production optimizations, namely injection-production optimization, as an example, Fig. 6 is a general process of implementing the injection-production optimization based on the reservoir dynamic understanding.

[0147] The most critical source of the above-mentioned closed-loop reservoir dynamic understanding update is still the regional understanding obtained by the main controlling factor analysis and the individual understanding obtained by the regular dynamic analysis, because other production optimization methods are limited by the optimization frequency, and the optimization conclusions obtained by them are more of a phased supplement to the reservoir dynamic understanding, which makes the "reasonable development technical strategy" derived therefrom sometimes experience-based; in order to further weaken this "empirical expression", the comprehensive production optimization method proposed in the application regards "production optimization based on reservoir dynamic understanding" as both an independent production optimization means and a source of partial constraint conditions for other production optimization methods relying on prediction models. As shown in FIG. 5, the theoretical model understanding can assist the update of the single-well horizontal understanding, and the single-well horizontal understanding in turn can be used as a constraint condition for other methods when implementing optimization operations; for example, according to the flow field diagnosis result (or the remaining oil distribution of numerical simulation) of the streamline simulation model, it is considered that there is more movable remaining oil near the production well A and the recovery efficiency is low, and liquid production should be increased; then, according to this conclusion, the reservoir dynamic understanding is updated as: the flow rate through the well near the well can be appropriately increased, and then the reasonable limits of the production indexes (such as water cut and liquid production) of the related injection wells and the production well are determined through the recent dynamic analysis and the node system analysis, and when the production optimization method of the streamline simulation model is used again to give the development control scheme, these constraint conditions should be applied.

[0148] In step 103, based on the production history data, the model parameters of the plurality of production optimization models are cyclically history-fitted until the optimal model parameters are obtained, by the following steps: the model parameters of each production optimization model are history-fitted, and the interwell connectivity relationship of the plurality of production optimization models is calibrated, so that the first type of prediction model provides an initial estimation of the interwell connectivity relationship for the second type of prediction model, and the second type of prediction model calibrates the model parameters of the first type of prediction model.

[0149] The prediction model class production optimization method needs to use the production history data of the target reservoir to perform history matching on the production optimization model parameters. The simplified physical model and the pure data-driven agent model do not need geological modeling, and the history matching only needs daily production data of the production history data, has low calculation cost and high efficiency, but has weak long-term prediction ability; the reservoir numerical model needs geological modeling, and the history matching process is relatively complicated and time-consuming, but has strong long-term prediction ability. The history matching parameters of the capacitance resistance model (CRM) in the simplified physical model mainly include the connectivity coefficient, the time constant, and the control volume, and the history matching parameters of the connected unit simulation model (LUSM) mainly include the initial conductivity and the initial control volume of the connected unit. The finite difference (FD) model and the streamline (SL) simulation model generally adjust the grid attribute parameters (porosity, permeability), the relative permeability curve (irreducible water saturation, residual oil saturation, and the corresponding relative permeability), the capillary force curve, the water body type, and the energy to history match the water cut and the bottom hole pressure change under the fixed liquid amount, and the history matching parameters are huge, and the single simulation takes a long time, so the process is generally divided into regions and steps to gradually improve the history matching effect. In the history matching process of the reservoir numerical model (FD or SL), the connectivity coefficient of the CRM and the initial conductivity of the LUSM can be referred to, and the permeability or conductivity of the grid between the wells with large CRM connectivity coefficients and high LUSM initial conductivities is appropriately increased, and the permeability or conductivity of the grid between the wells with small CRM connectivity coefficients and low LUSM initial conductivities is appropriately reduced, so as to speed up the entire history matching process. As described above, machine learning agents (such as CNN-LSTM) and history matching algorithms (such as EnKF and ES-MDA) can also be used for automatic history matching.

[0150] Referring to FIG. 3, the prediction class production optimization model in the embodiment of the application is mainly for the simplified physical model and the reservoir numerical model. The simplified physical model can have the capacitance resistance model (CRM), the connected unit model, and the flow network model. The reservoir numerical model can have the finite difference model and the streamline simulation model, the numerical simulator agent model, and the streamline simulation agent model. However, the current production optimization model has different problems.

[0151] Therefore, the embodiment of the application improves the existing production optimization model and obtains a plurality of innovative production optimization models. The following will be introduced respectively.

[0152] In an embodiment, the first type of prediction model is a capacitance resistance model, a connected unit model, or a flow network model.

[0153] The first type is a time-varying capacitance resistance model.

[0154] The time-varying capacitance resistance model (TV-CRM) belongs to the simplified physical model, and the method further includes:

[0155] On the basis of the approximation of the oil-water two-phase material balance equation and the deliverability equation, the flow control equation and the saturation control equation are determined;

[0156] According to the flow control equation and the saturation control equation, the time-varying capacitance resistance model is established.

[0157] The above embodiment gives the steps of constructing the time-varying capacitance resistance model.

[0158] Firstly, on the basis of the approximation of the oil-water two-phase material balance equation and the deliverability equation, that is, considering the influence of the oil and water volume coefficients, the continuity principle is applied to the CRM control body for the oil and water phases respectively, and then the following material balance relationship exists:

[0159] In the formula, q o is the surface oil production, m 3 / s; B o is the oil volume coefficient, m 3 / m 3 ; V f is the apparent volume of the control body, m 3 ; φ is the porosity, decimal; S o is the oil saturation of the control body; C o is the oil compressibility, Pa -1 ; C φ is the rock compressibility, Pa -1 ; p ave is the average formation pressure of the control body, Pa; t is time, s; w is the surface injection, m 3 / s; q w is the surface water production, m 3 / s; B w is the water volume coefficient, m 3 / m 3 ; S w is the water saturation of the control body; C w is the water compressibility, Pa -1 .

[0160] The CRM basic differential equation is obtained by adding equation (1) and equation (2):

[0161] In the formula, C t and q t are respectively: C t =C φ +S o C o +S w C w (4) q t =q o Bo +q w B w (5)

[0162] In the formula: V p To control the volumetric pore volume, m 3 C t Pa is the system / comprehensive compression factor. -1 ;q t The underground liquid production is expressed in m. 3 / s.

[0163] Based on the pressure and flow rate relationship in the bounded formation boundary control flow stage, the production capacity equation of CRM under oil-water two-phase flow conditions can be approximated as:

[0164] In the formula: J is the production capacity index or liquid production index, m 3 / (Pa·s); p wf Here is the bottom hole pressure, in Pa; K ro K represents the relative permeability of the oil phase. rw The relative permeability of the aqueous phase; μ o Oil phase viscosity, Pa·s; μ w The viscosity of the aqueous phase is expressed in Pa·s; r w K is the well diameter, in meters; K is the effective permeability, in meters. 2 h is the effective reservoir thickness, in meters; A e To control the equivalent cross-sectional area of ​​the volume, m 2 γ is Euler's constant, 0.577215664901532; C A is the shape factor, which is dimensionless.

[0165] The production capacity index J here is a function of the saturation level within the CRM control system, and therefore changes over time; substituting equation (6) into equation (3), and separating the variables to integrate, we can obtain:

[0166] In the formula: τ is the time constant, s; t0 is the reference time, s.

[0167] In the interval [t] k ,t k+1 If we consider τ(ξ) in the integrand of equation (7) as a constant, then the equation can be simplified to:

[0168] Where, q t (t k+1 ) and q t (t k ) are respectively t k+1 Time and t k The underground liquid production at time t, J(t)k+1 ) and J(t k ) are the productivity indices at t k+1 and t k , respectively, τ(t k+1 ) and τ(t k ) are the time constants at t k+1 and t k , respectively, w(t k+1 ) is the surface injection rate at t k+1 , B w is the water phase volume factor, p wf (t k+1 ) and p wf (t k ) are the bottomhole pressures at t k+1 and t k , respectively.

[0169] Equation (9) is the flow control equation of the time-varying capacitance-resistance model (TV-CRM).

[0170] It can be directly applied for one injection and one production or for the whole reservoir region. If there are N inj injection wells and N pr production wells, the flow control equations of the control volume centered at production well j (TV-CRMP) and the control volume based on injection-production well pair (i,j) (TV-CRMIP) can be expressed as:

[0171] where subscript j represents the physical quantity corresponding to production well j; subscript i represents the physical quantity corresponding to injection well i; f ij is the connectivity factor (dimensionless) representing the connectivity degree of production well pair (i,j); q t,ji is the subsurface liquid production rate of injection-production well pair (i,j), m 3 / s; J ji is the productivity index of injection-production well pair (i,j), m 3 / (Pa·s); τ ji is the time constant of injection-production well pair (i,j), s; t k and t k+1 are the start and end times of the kth time period, s.

[0172] From the basic differential equation (3) of the capacitance-resistance model

[0173] Substituting equation (12) into oil-water material balance equations (1) and (2) gives

[0174] where S w(t k+1 ) and S w (t k ) are water saturation of the control volume of the production well at t k+1 and t k , respectively, V p is the pore volume of the control volume of the production well, w is the surface injection rate, q w is the surface water production rate, C w is the water phase compressibility, C φ is the rock compressibility, C t is the overall compressibility, q t is the subsurface liquid production rate, S o (t k+1 ) and S o (t k ) are oil saturation of the control volume of the production well at t k+1 and t k , respectively, q o is the surface oil production rate, B o is the oil phase volume factor, C o is the oil compressibility.

[0175] Equations (13) and (14) are the saturation control equations of the TV-CRM.

[0176] The following is an example based on the control volume of the production well j, then equation (13) can be converted to:

[0177] Equations (10) and (15) constitute the basic control equations of the TV-CRM, and the change of saturation S w represents the evolution of the corresponding productivity index J of the control volume, the compressibility C t , the time constant τ, the volume factor B o and B w , etc. The TV-CRM model parameters f ij , V p,j (t0), J j (t0), S w,j (t0) can be solved by using the production history data through the history matching process, and the connectivity coefficient f ij inverted reflects the connectivity strength of the injection well i to the production well j.

[0178] The history matching process of the time-varying resistance capacitance model TV-CRM is shown in FIG. 7. In an embodiment, based on the production history data, the model parameters of a plurality of production optimization models are history matched, including:

[0179] When the production optimization model is a time-varying capacitance resistance model, the model parameters are determined as the connectivity coefficient, the control volume pore volume, the productivity index and the saturation.

[0180] determining initial values of the model parameters;

[0181] determining an evolution process of the model parameters;

[0182] solving iteratively a constrained optimization problem constituted by the model errors based on the initial values and the evolution process of the model parameters, and determining the optimized model parameters when the model errors are less than a preset threshold.

[0183] wherein the model error is an estimated deviation of the underground liquid production, as follows:

[0184] wherein: N pro is the number of production wells; N inj is the number of injection wells; error is the model estimated error of the underground liquid production of any moment of the production well, (m 3 / s) 2 ; N tim is the number of production history data points (time steps); E sum is the total sum of the estimated errors of the underground liquid production of all production wells at all time sequences, (m 3 / s) 2 ; is the calculated value of the underground liquid production of the jth production well at the k+1 moment, is the true value of the underground liquid production of the jth production well at the k+1 moment;

[0185] The essence of the TV-CRM history fitting is to constantly update the to-be-solved model parameters f ij , V p,j (t0), J j (t0), S w,j (t0) by combining the optimization algorithm to minimize the error, i.e., to solve the following constrained optimization problem:

[0186] wherein S wc is the irreducible water saturation; S or is the residual oil saturation; V pt is the total pore volume of the production system (upper limit of the pore volume of all production wells), f ij is the connectivity coefficient of the injection well i to the production well j, is the average water saturation in the control body of the production well j at the initial moment, V p,j is the pore volume of the production well j.

[0187] The production optimization process of the time-varying resistance capacitance model TV-CRM is shown in FIG. 8. In an embodiment, according to the updated dynamic understanding of the reservoir, the constraint condition corresponding to each production optimization model is determined, the objective function corresponding to each production optimization model is solved, the optimal decision variable corresponding to each production optimization model is obtained, and the development control scheme corresponding to each production optimization model is obtained, including:

[0188] When the production optimization model is the time-varying resistance capacitance model, the injection rate of the injection well and the bottom hole pressure of the production well are taken as the decision variable, and the initial value of the decision variable is determined;

[0189] According to the initial value of the decision variable, the underground liquid production is calculated through the fractional flow model, and the surface water production and the water saturation are updated;

[0190] Based on the surface water production and the water saturation, the surface oil production and the underground liquid production are predicted through the fractional flow model;

[0191] According to the reasonable development technical strategy, the constraint condition corresponding to the time-varying resistance capacitance model is determined, the production optimization objective function corresponding to the time-varying resistance capacitance model is iteratively solved, and the decision variable when the objective function reaches the optimal value is determined as the optimal decision variable;

[0192] According to the optimal decision variable, the development control scheme is obtained.

[0193] In an embodiment, the fractional flow model represents the relationship between the surface oil production and the underground liquid production, and the relationship between the surface water production and the underground liquid production through the saturation inverted by the historical fitting of the time-varying resistance capacitance model. The relationship between the surface oil production and the underground liquid production is as follows:

[0194] The relationship between the surface water production and the underground liquid production is as follows:

[0195] Wherein, F o,j is the ratio of the surface oil production and the underground liquid production of the production well j, dimensionless; q o,j is the surface oil production of the production well j, m 3 / s; q t,j is the underground liquid production of the production well j, m 3 / s; a i is the constant coefficient for predicting the oil production; S w,j is the water saturation in the control body of the production well j, decimal; N is the polynomial degree; F w,j is the ratio of the surface water production and the underground liquid production of the production well j, dimensionless; q w,j is the surface water production of the production well j, m 3 / s; b i is the constant coefficient for predicting the water production.

[0196] The essence of TV-CRM production optimization is similar to its history fitting, i.e. combining the optimization algorithm to update the decision variable w i wf,j to maximize the objective function (such as net present value).

[0197] In an embodiment, the target function corresponding to the time-varying capacitance resistance model is the net present value, and the net present value is as follows:

[0198] The constraint condition corresponding to the time-varying capacitance resistance model is as follows:

[0199] Where w lim is the injection limit of a single well, m 3 / s; w sum is the upper limit of the total injection, m 3 / s; q sum is the upper limit of the total liquid production, m 3 / s; NPV is the net present value, yuan; p wf,min is the lower limit of the bottom hole pressure of the production well, Pa; p wf,max is the upper limit of the bottom hole pressure of the production well, Pa; N t is the total number of optimization time steps; N pro is the number of production wells; N inj is the number of injection wells; r o (t k ) is the oil price at time t k ; b is the annual discount rate, q o,j (t k ) is the surface oil production of the production well j at time t k , r w (t k ) is the water treatment cost at time t k , q w,j (t k ) is the surface water production of the production well j at time t k , r wi (t k ) is the injection cost at time t k , w i (t k ) is the surface injection of the injection well i at time t k ; t opt is the initial time of the optimization period.

[0200] The time-varying capacitance resistance model in the embodiments of the present application can be used for oil reservoir comprehensive production optimization, and can also be used for oil reservoir production optimization alone. The following gives a specific embodiment to illustrate the process of production optimization alone by the time-varying capacitance resistance model.

[0201] ​Take the production optimization as an example of injection-production optimization, a water drive oilfield production area (Example 1) has 69 production wells and 57 injection wells. The block has a production history of nearly ten years, the yield and pressure change law is complex, and it is difficult to identify the injection-production corresponding relationship artificially and cannot be quantified. The irreducible water saturation S wc = 0.28695, the residual oil saturation S or = 0.19737. Take the rock compression coefficient C φ = 5 x 10 -4 MPa -1 , the water compression coefficient C w = 4 x 10 -4 MPa -1 , the oil compression coefficient C o = 4.5 x 10 -4 MPa - 1 , the water phase viscosity μ w = 0.51 mPa·s, the oil phase viscosity μ o = 1.70 mPa·s.

[0202] The injection-production history data of Example 1 is fitted by TV-CRM, the interwell connectivity relationship is shown in Figure 9, the history matching effect of underground liquid production is shown in Figure 10, the estimated initial control volume, initial productivity index and initial water cut are shown in Figures 11, 12 and 13 respectively. Figure 14 shows the history matching effect of the water cut of each production well.

[0203] After completing the TV-CRM history matching, production optimization is carried out, the lower and upper bounds of the injection rate are set to 0 and 344 m 3 / d respectively, and the upper bound of the total injection rate is 7300 m 3 / d, thus obtaining the injection-production optimization plan (as a development control scheme) for the next 12 optimization steps (step length 30d) as shown in Figure 15, A) of Figure 15 is the injection rate setting of injection wells, B) of Figure 15 is the bottom hole pressure setting of production wells, C) and D) of Figure 15 are the liquid production and oil production predicted by the TV-CRM under the injection-production setting respectively.

[0204] Second: connected unit simulation model (LUSM)

[0205] The interwell numerical simulation model (INSIM) discretizes the reservoir system into a series of interwell connected units characterized by interwell transmissibility and connected volume, and the essence is the balance of flow rate, subsurface injection / production (source / sink) and system compression term in the connected unit, and the implicit solution of pressure and explicit solution of saturation are performed by using material balance equation, and the saturation is approximately estimated by using Buckley-Leverett equation, and although the liquid flow diversion such as well shut-in / switching injection is considered, the reasonable arrangement of non-well nodes is not considered at the beginning and the saturation calculation theory is not strict. The extension of the INSIM considering the arrangement of non-well nodes on the multi-layer reservoir is referred to as the connected unit simulation model (LUSM), that is, the multi-layer non-well node flow simulation model is extended, and it is believed that the "front tracking algorithm" in streamline simulation should be used for saturation estimation.

[0206] Figure 16 shows the connected unit related to node i; for the connected unit composed of two nodes i and j, without paying attention to the specific shape and control range of the real form, the flow characteristics of the connected unit can be simplified by using the inter-node transmissibility T ij and control (porosity) volume V p,ij .

[0207] Figure 17 shows the connected unit diagram of the multi-layer reservoir; if the two wells are far apart or need to be simulated in detail in the area of interest, non-well nodes can also be added, and there is no fluid injection and production at the non-well nodes, that is, the nodes in the embodiments of the present application include well nodes and non-well nodes.

[0208] In an embodiment, based on the production history data, the model parameters of a plurality of production optimization models are history matched, including:

[0209] When the production optimization model is a connected unit simulation model, the model parameters are determined as transmissibility, porosity volume and productivity index, the nodes in the connected unit simulation model include well nodes and non-well nodes, and the non-well nodes are arranged between well nodes with a distance exceeding a preset range or at positions requiring detailed simulation in the area of interest;

[0210] Determine the initial value of the model parameter;

[0211] Determine the node saturation by using the front tracking algorithm in streamline simulation;

[0212] According to the node saturation, the water cut of each node is determined;

[0213] According to the water cut of each node, the oil production of each node is calculated;

[0214] The constraint optimization problem composed of the oil production and the bottom hole flowing pressure of each node is calculated as the constraint optimization problem of the connected unit simulation model;

[0215] Based on the initial value of the model parameters, the constraint optimization problem of the connection unit simulation model is iteratively solved to determine the model parameters when the model error is less than the preset threshold as the optimized model parameters.

[0216] In an embodiment, the node saturation is determined using a front tracking algorithm in streamline simulation, including:

[0217] The pressure value at each node is calculated;

[0218] According to the pressure value at each node, the flow rate on the connection unit is calculated, wherein two nodes constitute a connection unit;

[0219] The node saturation is determined using a front tracking algorithm in streamline simulation.

[0220] Taking injection-production optimization as an example, consider N L In the case of multi-well injection-production of a water drive reservoir, without considering the change of oil and water viscosity caused by temperature and pressure change during fluid flow, ignoring interlayer channeling, gravity term and capillary force, the material balance relationship of the ith node can be expressed as the balance of the injection / production term, the adjacent node flow term and the system compression term:

[0221] In the formula: i and j are the serial numbers of the nodes; k is the serial number of the reservoir; t is the production time, s; T ijk is the average conductivity between nodes i and j in the kth layer, m 3 / (s·Pa); N w is the number of nodes related to i; N L is the number of layers of the reservoir; p ik , p jk is the formation pressure of the kth layer represented within the control range of nodes i and j, Pa; q i is the subsurface flow at node i (source-sink term: positive for inflow, negative for outflow), m 3 / s; V p,ik is the pore volume represented by node i in the kth layer, m 3 ; V p,ijk is the pore volume of the connection unit composed of nodes i and j in the kth layer, m 3 ; A ijk is the longitudinal equivalent cross section of the simplified control unit between nodes i and j in the kth layer, m 2 ; L ijk is the equivalent length of the simplified control unit between nodes i and j in the kth layer, m; is the average effective porosity of the connection unit (i,j) in the kth layer; C t,ik is the system (integrated) compression coefficient of node i in the kth layer, Pa -1 ; Cφ Kw is the compressibility of water in the formation, Pa -1 ; C w Kw is the compressibility of water in the formation, Pa -1 ; S o,ik Sok is the average oil saturation in the kth layer controlled by node i; Sok is the average oil saturation in the kth layer controlled by node i;

[0222] The difference discrete of equation (22) can be obtained as

[0223] In the equation, the superscript n represents the time t n

[0224] According to Darcy's law, the flow rate between nodes i and j in the kth layer can be expressed as:

[0225] In the equation, q ijk is the flow rate between nodes i and j in the kth layer, m 3 / s; q o,ijk is the oil phase flow rate between nodes i and j in the kth layer, m 3 / s; q w,ijk is the water phase flow rate between nodes i and j in the kth layer, m 3 / s; K ijk is the effective permeability of the formation between nodes i and j in the kth layer, m 2 ; K ro is the relative permeability of oil; K rw is the relative permeability of water; μ o,k is the viscosity of oil in the kth layer, Pa·s; μ w,k is the viscosity of water in the kth layer, Pa·s.

[0226] According to the definition of conductivity, the T ijk of the simplified control unit between nodes i and j in the kth layer can be written as:

[0227] In the equation, λ ijk is the total oil-water mobility of the connection unit formed by nodes i and j in the kth layer, m 2 / (Pa·s).

[0228] Since the average water saturation of each connection unit and the node saturation change with time, the physical quantities related to saturation (such as K ro , K rw , λ ijk , T ijk ) also evolve with time. In order to further simplify the calculation, the conductivity is explicitly estimated using the value of the previous time, i.e. ​

[0229] Where t n-1 The total oil-water fluidity at any given time is calculated according to the following rules:

[0230] According to the definition of rock pore compressibility coefficient, the pore volume represented by node i in the k-th layer is:

[0231] Equation (25) can be written in a form that includes interlayer flow.

[0232] In the formula: q i,k Let m be the source and sink flow of node i in layer k. 3 / s.

[0233] Arrange equation (32) to obtain

[0234] Equation (34) can be simplified as follows:

[0235] When i takes values ​​from 1 to N w When (j≠i), equation (36) can be written in matrix form:

[0236] p n =C q -1 ·(p n-1 +ξ) (39)

[0237] In the formula:

[0238] Equation (39) establishes t n-1 The pressure vector p at time t n-1 and the next moment t n pressure vector p n The connection between them is necessary, but this requires traffic from each node. For known (or flow rates that can be split between layers according to the productivity index) quantities, all physical quantities except pressure and flow rate can use the values ​​from the previous time step. If the bottomhole pressure is known, and the production rate is the response of the bottomhole flowing pressure, then the flow rate needs to be converted to pressure using the productivity index, i.e. J i,j,k =α ijk (h ijk A e,ijk C A,ijk ,r w )·λ ijk (45)

[0239] In the formula: p wf,iB i,k is the bottom-hole pressure of node i, and its value is 0 if i is a non-well node. 3 is the productivity index of node i in the kth layer, m i,j,k / s; its value is 0 if i is a non-well node. 3 is the productivity index of node i in the kth layer due to the effect of node j on node i in the connecting unit (i, j), m ijk / s, its value is 0 if i is a non-well node. ijk is a constant related to the effective thickness h ijk , control area A A,ijk , shape factor C w , and wellbore radius r ik .

[0240] Since there is only one bottom-hole pressure for a well, the well flow rate q ik corresponding to different layers is adjusted by J w , so equation (32) can be written as:

[0241] According to equations (35) and (38), equation (47) can be obtained by rearranging:

[0242] When i takes 1 to N w (j≠i), equation (48) can be easily written in matrix form:

[0243]

[0244] where

[0245] The pressure vector p n at time t n when the bottom-hole flowing pressure is known can be calculated from equation (50), and then the flow rates at the connecting units and nodes can be obtained.

[0246] where q i,j,k is the flow rate in the connecting unit (i, j) in the kth layer, m 3 / s.

[0247] From equations (45) and (29), both the conductivity and the productivity index are related to the total mobility of oil and water. To associate them, A ijk , A e,ijk , and C A,ijk must be estimated according to the initial distance L ijk , formation thickness h ijk , and other information of the nodes in the connecting unit, i.e.

[0248] Thus, J ijk By A ijk associated with the conductivity parameter , while the evolution of the saturation S ijk ) i with the saturation S w,ik . Therefore, the initial conductivity and the initial pore volume are usually set as the model parameters to be solved in the connection unit simulation model; in order to adjust more flexibly, the initial productivity index may also be included in the model parameters to be fitted, and the initial value thereof can be estimated according to equations (54) to (56), and the value at other time is:

[0249] When the pressure values at the nodes are calculated, the flow rate on the connection unit can be obtained according to equation (53), and then the fluid saturation is tracked along the flow path according to the pressure size, so as to determine the water cut of each node. The initial INSIM model uses the one-dimensional Buckley-Leverett equation to inversely calculate the node saturation through the comparison of the water cut derivative, but for the liquid flow diversion situations such as well shut-in and injection switching, and the situations where the initial saturation is not the irreducible water saturation, the description is not strict. Therefore, the front tracking algorithm in streamline simulation is used to determine the node saturation. For the saturation tracking between nodes i and j, the following saturation problem needs to be solved, and a single-layer reservoir is taken as an example:

[0250] or

[0251] In the formula, S w is the water saturation; t is the time, s; q t,i,j is the formation flow rate on the connection unit (i, j), m 3 / s; f w is the water cut; x is the one-dimensional coordinate; A ij is the longitudinal equivalent cross section of the connection unit (i, j), m 2 ; φ ij is the effective porosity of the connection unit (i, j); L ij is the equivalent length of the connection unit (i, j), m; and V p,ij is the pore volume in the control range of the connection unit (i, j), m 3 .

[0252] The initial condition in equation (60) is simplified, and first, the case that the initial saturation is discontinuous only at a position x0 is considered, and the initial saturation is obtained, that is,

[0253] In the formula: τ is the characteristic velocity, m / s; S w,left The initial saturation value to the left of x0; S w,right This represents the initial saturation value to the right of x0;

[0254] (1) If f w ′(S w,left )≥v≥f w ′(S w,right The solution to equation (61) is:

[0255] (2) If f w ′(S w,left ) < f w ′(S w,right And f w "(S w,left )·f w "(S w,right )≥0, (or if S is satisfied) w,left >S w,right >S w,ip or S w,left w,right w,ip S w,ip f w ~S w The inflection point of the curve), the solution to equation (61) is:

[0256] (3) If the water content f at both ends of x0 w If the concavity and convexity of the wave are different, the solution is a composite wave, in which case an intermediate saturation S is required. w * :

[0257] In S w,right and S w * The water saturation S between w It should meet the following requirements:

[0258] Therefore, the solution to equation (61) is:

[0259] For the saturation problem shown in equation (61), it is necessary to determine the solution type from equations (63), (65), and (68) based on the initial saturation values ​​on both sides of position x0. However, in practice, it is often a continuous initial value problem as shown in equation (59) or (60). In this case, position x can be discretized into multiple intervals, and the saturation in each interval is constant, that is, S is... w Treating it as a piecewise constant function of x simplifies the continuous initial value problem into solving a subproblem with multiple discontinuities.​​

[0260] In one embodiment, the node saturation is determined using a front tracking algorithm in a streamline simulation, including:

[0261] determining a saturation solution problem;

[0262] simplifying initial conditions in the saturation solution problem to obtain an initial saturation;

[0263] simplifying the initial saturation into a multi-segment constant function to form a plurality of saturation solution sub-problems;

[0264] determining a new saturation solution sub-problem that can occur due to a collision of saturation fronts of two adjacent sub-problems during movement of the saturation fronts;

[0265] after the collision of the fronts, solving a new saturation solution sub-problem according to the collision point, finding a next possible collision and inserting a corresponding new front, until the collision exceeds a space and time boundary, determining a saturation profile at the end of the time step, the saturation profile at the end of the time step being a piecewise constant function, and using the saturation profile at the end of the time step as initial conditions for front tracking in a next time step;

[0266] determining the node saturation according to the saturation profile at the end of the time step.

[0267] As shown in FIG. 18, the initial saturation S w (x,0) is simplified into a 5-segment constant function divided by x1, x2, x3, x4, and thus can be decomposed into four saturation solution sub-problems in the form of equation (61), and a new saturation solution sub-problem can be generated due to a "collision" of saturation fronts of two adjacent sub-problems during movement of the saturation fronts:

[0268] where S w is the water saturation, x is a one-dimensional coordinate, f w is the water content, τ is the characteristic velocity, x collision is a position where the saturation collision occurs; t collision is a time when the saturation collision occurs, S w (x, t collision ) is the water saturation at x and t collision ;

[0269] after the collision of the fronts, solving a new saturation solution sub-problem according to the collision point, finding a next possible collision and inserting a corresponding new front, until the collision exceeds a space and time boundary (0≤x≤L i,j,0≤t≤Δt), determine the saturation profile at the end of the time step, which is a piecewise constant function, and use the saturation profile at the end of the time step as the initial condition for the front tracking in the next time step; determine the node saturation according to the saturation profile at the end of the time step.

[0270] When the saturation and water cut of each node are calculated, the injection-production performance of each layer can be determined. The water cut of node i in the kth layer is:

[0271] In the formula: is the water cut of node i in the kth layer at time t; n is the water cut of node i in the kth layer at time t; n

[0272] After integrating the flow rates of all layers, the water cut and oil production of node i (when there is a sink term) are:

[0273] The proportion of the injection / production flow rate of node i in the kth layer to the total flow rate (i.e., the interlayer flow distribution coefficient) is:

[0274] The proportion of the flow rate of node i on the connecting unit (i,j) in the kth layer to the total injection / production flow rate of node i, i.e., the injection distribution coefficient or the produced fluid distribution coefficient, is:

[0275] The flow rate distribution coefficient of node j to node i on the connecting unit (i,j) in the kth layer is:

[0276] Due to the influence of the compressibility of the reservoir system, part of the injection (source term) may be absorbed by the reservoir, and part of the produced fluid (sink term) may come from the elastic action of the reservoir, so the distribution coefficient calculated by formula (75) is generally greater than the value of formula (74).

[0277] The history matching of the connecting unit simulation model requires solving the initial value of the conductivity (and the initial value of the saturation and the characteristic cross section A ijk or A e,ijk related) and the initial value of the connecting unit control volume Generally, the bottom hole pressure and water cut (or oil production / water production) are matched by fixing the liquid volume; it should be noted that when the calculated bottom hole pressure does not meet the liquid volume requirement (such as the p wf of the injection well is lower than its average formation pressure, and the p wf ​​If the reservoir pressure is negative or higher than its average reservoir pressure, the bottom hole flowing pressure should be used as a constraint, i.e. the constraint optimization problem of the oil production of each node and the bottom hole flowing pressure is:

[0278] wherein N w is the number of nodes related to i, N tim is the number of time steps of the history matching, is the calculated value of the bottom hole pressure of node i at t n , is the true value of the bottom hole pressure of node i at t n , is the calculated value of the formation oil production of node i at t n , is the true value of the formation oil production of node i at t n , is the initial value of the conductivity of the connecting unit composed of nodes i and j in the kth layer, is the initial value of the pore volume of the connecting unit composed of nodes i and j in the kth layer, S wc is the irreducible water saturation, is the initial value of the water saturation in the control range of node i in the kth layer, S or is the residual oil saturation, u is the decision variable; V pt is the total pore volume, m 3 ; p wf,inj is the bottom hole pressure of the injection well node, Pa; p inj is the average reservoir pressure of the injection well node, Pa; p wf,pro is the bottom hole pressure of the production well node, Pa; p pro is the average reservoir pressure of the production well node, Pa.

[0279] After the history matching of the connecting unit simulation model is completed, different injection and production flow rates can be set to predict the oil and water production performance of each production well and then calculate the income. The production optimization process is similar to the aforementioned TV-CRM, which will not be described here. The objective function corresponding to the connecting unit simulation model is consistent with that of the time-varying capacitance resistance model, which is formula (21).

[0280] The connecting unit simulation model in the embodiments of the present application can be used for comprehensive production optimization of an oil reservoir, and can also be used for separate production optimization of an oil reservoir. A specific embodiment is given below to illustrate the process of separate production optimization of the connecting unit simulation model.

[0281] A certain water drive oil reservoir (Embodiment 2) contains two special geological features and well types, i.e. horizontal wells and water bodies. When the connecting unit simulation model is applied, an equivalent simplified processing method is considered, i.e. the horizontal well is equivalent to three nodes, and the numerical water body is used to represent the edge and bottom water.

[0282] The initial connectivity model is established by using the coordinate information of the existing well points, considering the permeability, porosity and thickness information at the well points. The initial value of the connectivity parameter of each connection unit is obtained from the actual distance between nodes and geological data such as reservoir thickness, as shown in FIG. 19. The model combines geological information and can roughly reflect the interwell connectivity.

[0283] The model after history matching can better reflect the real production performance of the block and can be used as an auxiliary tool for subsequent reservoir development adjustment and optimization. FIG. 20 shows the fitting effect of some single wells.

[0284] The initial value of the connection unit conductivity and the initial value of the connected volume obtained after history matching are shown in FIG. 21, and the water injection splitting coefficient is shown in FIG. 22. These information, combined with the connectivity results of the CRM inversion and the field tracer monitoring data, can more reliably determine the changes of the injection-production correspondence, thereby guiding the subsequent injection-production adjustment.

[0285] Taking equation (21) as the objective function, setting 6 optimization steps and a step length of 6 months, and a total optimization period of 3 years. The liquid production of the production well is set to be 0.8-1.2 times the original liquid production, the water injection of the injection well is set to be 0.5-1.5 times the original injection, and the overall injection-production ratio of the reservoir is controlled within the range of 0.8-1.2. The injection-production system of a total of 226 wells (119 injection wells and 107 production wells) is adjusted. The optimal injection-production scheme obtained by the connection unit simulation model shows that 65 water injection wells need to be increased, 54 water injection wells need to be reduced, 17 oil wells need to be increased, and 90 oil wells need to be reduced. The injection-production system of some water injection wells is shown in FIG. 23, and the change of the liquid production of some production wells before and after optimization is shown in FIG. 24.

[0286] The third type is a flow network model (FNM)

[0287] The flow network model is another representative simplified physical model developed in recent years; it further introduces the step of grid division on the basis of CRMIP and INSIM, simplifying the multi-dimensional reservoir space into a series of one-dimensional flow networks between each well pair, coupling each one-dimensional model at the intersecting nodes, and the flow in each connection depends on the fluid volume and the corresponding reservoir permeability in the displacement region covered by the two wells. The two parameters are adjusted by using the output of the full-order reservoir model or the field observation data to calibrate the flow network model and then implement production optimization.

[0288] In an embodiment, the second type of prediction model is a pure data-driven proxy model or a non-proxy model in a reservoir numerical model, and the non-proxy model is a finite difference model and a streamline simulation model. Next, they are introduced in turn.

[0289] Fourth: Pure data-driven agent model

[0290] The pure data-driven agent model is a kind of machine learning agent model, wherein, for the pure data-driven agent model, the training and production optimization process of the agent is separated, and the constructed pure data-driven agent model needs to mine the internal law about the actual injection-production response in the production data to fully depict the historical dynamic, but cannot appear "overfitting" phenomenon to avoid reducing its ability to predict future dynamic. Fig. 25 shows the general process of pure data-driven agent model optimization in the field of reservoir production optimization, including:

[0291] 1. Define the production optimization problem;

[0292] 2. Prepare training data;

[0293] 3. Set agent options;

[0294] 4. Train the pure data-driven agent model;

[0295] 5. Test the pure data-driven agent model;

[0296] 6. Use the pure data-driven agent model for reservoir production optimization;

[0297] 7. Output the optimal decision variable.

[0298] Among them, steps 2-5 are the process of historical fitting of model parameters of various production optimization models; steps 6-7 are the process of solving the objective function corresponding to each production optimization model, obtaining the optimal decision variable corresponding to each production optimization model, and obtaining the development control scheme corresponding to each production optimization model.

[0299] Fifth: Streamline simulation model (SL)

[0300] Based on the production history data, the model parameters of the streamline simulation model (SL) are cyclically history-fitted by the following steps until the optimal model parameters are obtained, wherein the first type of prediction model provides an initial estimate of the interwell connectivity relationship for the streamline simulation model (SL), and the streamline simulation model (SL) calibrates the model parameters of the first type of prediction model; according to the updated reservoir dynamic understanding, the constraint conditions corresponding to the streamline simulation model are determined, the objective function corresponding to the streamline simulation model is solved, the optimal decision variable corresponding to the streamline simulation model is obtained, and the development control scheme corresponding to the streamline simulation model is obtained.

[0301] In addition to the above application in the comprehensive production optimization of oil reservoir, the streamline simulation model can also be used to complete the flow field diagnosis alone. The flow field diagnosis refers to using the streamline information output by the streamline simulator or the streamline information indirectly derived by the streamline tracking based on the simulation results of the finite difference simulator, diagnosing the production level of the injection-production well pair or single well based on the flow field diagnosis index (such as injection efficiency, production efficiency, movable residual oil, etc.), determining the streamline injection-production optimization strategy, and then optimizing the injection rate and liquid production rate of the single well according to the streamline injection-production optimization strategy.

[0302] In an embodiment, the production optimization based on the flow field diagnosis of the streamline simulation model is performed by the following steps:

[0303] Using the streamline information output by the streamline tracking in the streamline simulator or the streamline information obtained by the streamline tracking based on the simulation results of the finite difference simulator;

[0304] Diagnosing the production level of the injection-production well pair or single well based on the flow field diagnosis index, determining the streamline injection-production optimization strategy, and the flow field diagnosis index includes injection efficiency, production efficiency, movable residual oil displacement efficiency, and movable residual oil production efficiency, and the streamline injection-production optimization strategy includes the flow rate updating mode and the injection-production flow rate optimization criterion, and the injection-production flow rate optimization criterion is used to determine the flow rate updating weight.

[0305] Based on the streamline injection-production optimization strategy, the streamline injection-production optimization is performed to obtain the updated injection-production flow rate.

[0306] Adding constraint processing to the updated injection-production flow rate to obtain the optimal injection-production scheme, and the optimal injection-production scheme is the development control scheme.

[0307] In an embodiment, the flow rate updating mode includes: only updating the flow rate at the injection end and the flow rate at the production end of the injection-production well pair or single well according to the movable residual oil displacement efficiency; only updating the flow rate at the injection end and the flow rate at the production end of the injection-production well pair or single well according to the movable residual oil production efficiency; updating the flow rate at the injection end of the injection-production well pair or single well according to the movable residual oil displacement efficiency, and updating the flow rate at the production end of the injection-production well pair or single well according to the movable residual oil production efficiency; only updating the flow rate at the injection end and the flow rate at the production end of the injection-production well pair or single well according to the injection efficiency; only updating the flow rate at the injection end and the flow rate at the production end of the injection-production well pair or single well according to the production efficiency; updating the flow rate at the injection end of the injection-production well pair or single well according to the injection efficiency, and updating the flow rate at the production end of the injection-production well pair or single well according to the production efficiency.

[0308] In an embodiment, the flow rate optimization criterion includes the unbounded weight criterion and the bounded weight criterion.

[0309] The unbounded weight criterion does not consider the amplitude limit of the flow rate change, and only determines the flow rate updating weight according to the difference between the flow field diagnosis index and the average value.

[0310] The bounded weight criterion introduces a diagnostic indicator of difference and a linear scaling ratio of the diagnostic indicator of difference to limit the weight variation amplitude.

[0311] In an embodiment, the constraint processing comprises:

[0312] For the injection-production well pair, the updated injection-production flow rate is converted to single well by superposition;

[0313] For the injection-production well pair, if there is an upper bound of the total injection of all injection wells and an upper bound of the total production of all production wells, the optimized single well injection-production flow rate is scaled;

[0314] For the injection-production well pair, if there is a constraint condition of single well, the limit flow rate of the injection-production well that violates the pressure constraint is taken as the final optimized flow rate of the injection-production well that violates the pressure constraint, and the remaining injection-production is proportionally scaled to the flow rate of the injection-production well that does not violate the pressure constraint again; or only the flow rate of the injection-production well that violates the pressure constraint is reduced and the injection-production well that does not violate the pressure constraint is not processed.

[0315] For the single well, if there is a total flow rate constraint or a bottom hole pressure constraint, the optimized single well injection-production flow rate is scaled.

[0316] Based on the history-fitted reservoir numerical model, relevant streamline information (such as injection allocation factor, produced fluid allocation factor, flight time field along the streamline, saturation field, etc.) is generated by using the streamline post-processing tool of the streamline simulator or the finite difference simulator, the injection-production correspondence relationship of the target reservoir can be quantitatively identified, the injection direction and the produced fluid source are determined, the update of reservoir understanding is assisted, and thus the production optimization based on the streamline simulation flow field diagnosis is implemented.

[0317] Let the total number of injection wells be N inj , the total number of production wells be N pro , the injection rate of injection well i be w i , the produced fluid rate of production well j be q j , injection well i and production well j form an injection-production well pair (i, j); the injection fluid flow rate on the injection-production well pair (i, j) is q i,j (injection well metering), the produced fluid flow rate on the injection-production well pair (i, j) is q j,i (production well metering), the water cut of the injection-production well pair (i, j) at the production well j at the production end is f w,j,i , and the water cut of the production well j is f w,j . Let the injection allocation factor of injection well i to production well j be AFI i,j (0≤AFI i,j ≤1), and the produced fluid allocation factor of production well j to injection well i be AFP j,i (0≤AFP j,i ≤1), then

[0318] Two important physical quantities (flow field diagnostic indicators) for flow field diagnosis are defined from the perspective of injector-producer well pair in the following: mobile residual oil displacement ( / production) efficiency and injection ( / production) efficiency. Fig. 26 shows the pore volume of the kth streamline on the injector-producer well pair (i, j), which has no specific shape and is related to the flow rate carried by the streamline and the time-of-flight coordinate:

[0319] wherein: (V p ) i,j,k,s is the pore volume represented by the s th segment of the k th streamline on the injector-producer well pair (i, j), m 3 ; τ is the time-of-flight, s; is the time-of-flight coordinate of the starting position of the s th segment of the k th streamline on the injector-producer well pair (i, j), s; is the time-of-flight coordinate of the ending position of the s th segment of the k th streamline on the injector-producer well pair (i, j), s; q i,j,k is the flow rate carried by the k th streamline on the well pair (i, j), m 3 / s.

[0320] Based on the relevant streamline information generated by the streamline post-processing tool of the streamline simulator or the finite difference simulator, the residual oil saturation distribution along each streamline can be counted, and then the mobile residual oil volume on the streamline can be represented as:

[0321] wherein: V MRO,i,j,k is the mobile residual oil volume of the k th streamline on the injector-producer well pair (i, j), m 3 ; is the average value of the water saturation of the s th segment of the k th streamline on the well pair (i, j); Seg i,j,k is the total number of segments of the k th streamline on the well pair (i, j); S or is the residual oil saturation.

[0322] According to the mobile residual oil volume and the flow rate on the streamline, the mobile residual oil displacement efficiency of the injector-producer well pair can be defined as: the ratio of the mobile residual oil volume on the well pair to the flow rate at the injection end of the well pair, i.e.

[0323] wherein: EMRO i,j is the mobile residual oil displacement efficiency of the well pair (i, j), s; EMRO ave,i is the mobile residual oil displacement efficiency in the average sense of all injector-producer well pairs, s; N SL,i,j is the total number of streamlines connecting the injector-producer well pair (i, j).

[0324] The movable residual oil displacement efficiency of the well pair (i, j) reflects the matching relationship between the residual recoverable oil volume on the well pair and the current injection volume. The difference between the movable residual oil displacement efficiency and the average value thereof reflects the displacement balance degree of the injection fluid on the well pair. The higher the value is, the more unbalanced the displacement is. i,j Exceeding EMRO ave,i The more, the more insufficient the displacement on the injection-production well pair (i, j) is. This may be caused by a large residual oil volume and a small injection volume. Therefore, the injection flow rate q of the injection-production well pair should be increased i,j ; otherwise, if EMRO i,j is lower than EMRO ave,i , it indicates that the displacement on the injection-production well pair (i, j) is excessive. This may be caused by a small residual oil volume and a large injection volume. Therefore, the injection flow rate q of the injection-production well pair should be decreased i,j . According to the relative size of the movable residual oil displacement efficiency of the injection-production well pair and the average value thereof, the displacement balance degree of the injection-production well pair can be determined, and then the direction of flow rate increase or decrease can be determined. The increase or decrease amplitude, i.e., the flow rate update weight, is determined according to the injection-production flow rate optimization criterion.

[0325] Similarly, the movable residual oil recovery efficiency (the movable residual oil displacement efficiency corresponding to the production end) of the injection-production well pair (i, j) can also be defined at the production end as follows: the ratio of the movable residual oil volume to the flow rate of the well pair production end, i.e.,

[0326] In the formula, EMRO j,i is the movable residual oil recovery efficiency of the well pair (i, j), s; EMRO ave,j is the average movable residual oil recovery efficiency of all injection-production well pairs, s.

[0327] Similarly, according to the relative size of EMRO j,i and EMRO ave,j , the adjustment direction of the production end flow rate q j,i of the injection-production well pair (i, j) can be determined. If EMRO j,i is higher than EMRO ave,j , it indicates that the residual oil volume on the injection-production well pair (i, j) is relatively large and the liquid production volume is relatively small. The production end flow rate does not match the relatively large residual oil volume. Therefore, the production flow rate q of the injection-production well pair should be increased j,i ; otherwise, if EMRO j,i is lower than EMRO ave,j , it indicates that the residual oil volume on the injection-production well pair (i, j) is relatively small and the liquid production volume is relatively high. The production end flow rate does not match the relatively small residual oil volume. Therefore, the production flow rate q of the injection-production well pair should be decreased j,i . The specific value of the flow rate increase or decrease can be calculated according to the injection-production flow rate optimization criterion.

[0328] In addition to the mobile residual oil displacement ( / production) efficiency, the injection ( / production) efficiency of the well pair can also be defined according to the flow rate data:

[0329] wherein: EI i,j is the injection efficiency of the well pair (i, j), i.e. the ratio of the oil production at the production end of the well pair to the injection at the injection end; EI ave,i is the average injection efficiency of all well pairs, i.e. the injection efficiency of the entire study area; EP j,i is the production efficiency of the well pair (i, j), i.e. the ratio of the oil production at the production end of the well pair to the liquid production at the production end; EP ave,j is the average production efficiency of all well pairs, i.e. the production efficiency of the entire study area.

[0330] According to the relative size of EI i,j and EI ave,i , the increase or decrease of the injection flow rate q i,j can be determined, and according to the relative size of EP j,i and EP ave,j , the increase or decrease of the production flow rate q j,i can be determined, and the adjustment method is similar to the method determined according to the mobile residual oil displacement ( / production) efficiency, which will not be repeated here.

[0331] For a well pair, the following six flow rate updating methods can be determined according to the above flow field diagnostic indicators (mobile residual oil displacement efficiency, mobile residual oil production efficiency, injection efficiency, production efficiency):

[0332] 1) Update the flow rate q i,j at the injection end and the flow rate q j,i at the production end of the well pair only according to the mobile residual oil displacement efficiency;

[0333] 2) Update the flow rate q i,j at the injection end and the flow rate q j,i at the production end of the well pair only according to the mobile residual oil production efficiency;

[0334] 3) Update the flow rate q i,j at the injection end according to the mobile residual oil displacement efficiency, and update the flow rate q j,i at the production end according to the mobile residual oil production efficiency;

[0335] 4) Update the flow rate q i,j at the injection end and the flow rate q j,i at the production end of the well pair only according to the injection efficiency;

[0336] 5) Update the flow rate q i,jand the flow rate q at the production end j,i ;

[0337] 6) Update the flow rate q at the injection end according to the injection efficiency of the well pair i,j and the flow rate q at the production end according to the production efficiency of the well pair j,i .

[0338] When the injection-production flow rate optimization of the well pair is completed, constraints are added to the streamline injection-production optimization strategy, one of which is to convert the updated injection-production flow rate to a single well through superposition, i.e.

[0339] wherein: is the injection rate of the injection well i after streamline simulation optimization, m 3 / s; is the updated value of the loss injection rate of the injection well i after optimization; is the updated value of the flow rate at the injection well i of the well pair (i, j) after optimization, m 3 / s; is the liquid production rate of the production well j after streamline simulation optimization, m 3 / s; is the updated value of the injection rate of other energy sources (non-injection wells) to the production well j after optimization, m 3 / s; is the updated value of the flow rate at the production well j of the well pair (i, j) after optimization, m 3 / s; AQ i is the injection rate of the injection well i lost due to injection into non-reservoir, m 3 / s; r AQ is the proportion of the reduction in the loss injection rate.

[0340] The injection rate w of the injection well i i may not all be injected into the reservoir to directly provide energy support to the production well, but also may be lost to the water body (aquifer) or used to balance the compressibility of the reservoir system; similarly, the liquid production of the production well may not all come from the energy supply of the related injection well, but also may be partially affected by the water body (aquifer) or the compressibility of the reservoir system; therefore, the flow rate terms of the non-injection-production well pair are added to equations (89) and (90).

[0341] The second constraint is that if there is an upper bound for the total injection rate of all injection wells and an upper bound for the total production rate of all production wells, the optimized single-well injection-production flow rate also needs to be scaled, i.e.

[0342] wherein: is the injection rate of the injection well i after flow rate optimization and scaling, m 3 / s; w sum,ubTotal injection rate target or upper bound for all injectors, m 3 / s; Liquid production rate of production well j after flow rate optimization and scaling, m 3 / s; q sum,ub Total liquid production rate target or upper bound for all producers, m 3 / s.

[0343] If there are still single well pressure constraints, the limit flow rate of the injection-production well that violates the pressure constraint is taken as the final optimized flow rate of the injection-production well that violates the pressure constraint, and the excess injection-production rate is proportionally scaled again to the injection-production well that does not violate the pressure constraint; or only the flow rate of the injection-production well that violates the pressure constraint is reduced and the injection-production well that does not violate the pressure constraint is not processed. Here, a simple example is given, Table 1 gives an example of single well level flow rate scaling, a 2 injection 3 production model, the total injection rate and the liquid production rate upper bound are both 180 m 3 / d, and after flow line simulation optimization, the total injection-production flow rate is 200 m 3 / d, the single well injection-production rate needs to be scaled to meet the total flow rate upper bound constraint; however, since the injection well I3 and the production well P1 have bottom hole pressure limits, their flow rate upper bounds are 85 m 3 / d and 70 m 3 / d, respectively, so after scaling the single well flow rate again, the final injection-production optimization result is obtained.

[0344] Table 1 Single well level flow rate scaling example

[0345] The foregoing discusses the optimization idea of first optimizing the injection-production well pair flow rate according to the four flow field diagnostic indicators at the well pair level and then superimposing it to the single well level optimization, these flow field diagnostic indicators can also be defined at the single well level, and the flow rate addition step is omitted. The movable residual oil displacement efficiency, movable residual oil recovery efficiency, injection efficiency, and recovery efficiency at the single well level are defined as follows:

[0346] In the formula: EMRO i is the movable residual oil displacement efficiency of injection well i, s; EMRO j is the movable residual oil recovery efficiency of production well j, s; EI i is the injection efficiency of injection well i; EP j is the recovery efficiency of production well j.

[0347] Consistent with the foregoing flow field diagnosis at the well pair level, the movable residual oil displacement efficiency EMRO i , the movable residual oil recovery efficiency EMRO j , the injection efficiency EI i , and the recovery efficiency EP jCompare with the corresponding average reference value EMRO ave,i EMRO ave,j EI ave,i EP ave,j By comparing the data, the direction of increase or decrease in the horizontal flow rate of a single well can be determined, and there are also 6 flow rate update methods (similar to the flow rate update methods of the aforementioned well pairs):

[0348] 1) Update the flow rate w at the injection end of a single well based only on the movable residual oil displacement efficiency. i and the flow rate q at the extraction end j ;

[0349] 2) Update the flow rate w at the injection end of a single well only based on the movable remaining oil production efficiency. i and the flow rate q at the extraction end j ;

[0350] 3) Update the flow rate w at the injection end of a single well based on the movable residual oil displacement efficiency. i The flow rate q at the production end is updated based on the recoverable remaining oil production efficiency. j ;

[0351] 4) Update the flow rate w at the injection end of a single well only based on the injection efficiency. i and the flow rate q at the extraction end j ;

[0352] 5) Update the flow rate w at the injection end of a single well based only on the production efficiency. i and the flow rate q at the extraction end j ;

[0353] 6) Update the flow rate w at the injection end of a single well based on the injection efficiency. i The flow rate q at the extraction end is updated based on the extraction efficiency. j .

[0354] For a single well, the constraint treatment is as follows: if there is a total flow constraint or a bottom hole pressure constraint, the optimized single well injection and production flow rate is scaled in a manner similar to that in equations (92) and (93).

[0355] Thus far, 12 flow rate update methods based on four flow field diagnostic indicators have been proposed for both injection-production well pairs and single-well levels. The initial optimized flow rate value before scaling should be determined according to the "Injection-Production Flow Rate Optimization Criteria". For ease of writing and discussion, the injection-production flow rate to be optimized (q) will be used as the reference value. i,j q j,i or w i q j )Unified by q original This indicates that the corresponding updated injection / collection flow rate is q. opt The aforementioned eight flow field diagnostic indicators (well-to-level: EMRO)i,j , EMRO j,i , EI i,j , EP j,i ; single well horizontal: EMRO i , EMRO j , EI i , EP j ) are denoted by the efficiency symbol E, and its normalized value is denoted by β, which is illustrated by taking the linear normalization of the maximum value as an example:

[0356] wherein: E is the flow field diagnostic index (EMRO i,j , EMRO j,i , EI i,j , EP j,i or EMRO i , EMRO j , EI i , EP j ); E ave is the flow field diagnostic index in the average sense (EMRO ave,i , EMRO ave,j , EI ave,i , EP ave,j ); β is the normalized flow field diagnostic index, which is between 0 and 1; β ave is the normalized value of the flow field diagnostic index in the average sense.

[0357] As mentioned before, the increase or decrease amplitude of injection-production rate will be determined according to the difference degree between β and β ave , and the flow rate update weight W is defined as: q opt = W(β, β ave )·q original (100)

[0358] The injection-production rate optimization criterion specifies how to calculate the flow rate update weight W; wherein, the unbounded weight criterion does not consider the amplitude limit of flow rate change, and only determines the weight of flow rate update according to the difference degree between the flow field diagnostic index and the average value, i.e.

[0359] wherein: α is the weight index, and α > 0.

[0360] Figure 27 shows the change relationship of flow rate update weight W with normalized flow field diagnostic index β under different weight indexes in the unbounded weight criterion (taking β ave = 0.5 as an example); when 0 < α < 1, with the increase of α, W decreases in the range of β < β ave , and increases in the range of β > β ave; when a = 1, W varies linearly with β; when a > 1, W increases with β in the range of β ave ; when a > 1, W decreases with β in the range of β ave ; the selection of weight index a affects the result of flow optimization, so the weight index can be further optimized when using unbounded weight criterion for injection-production flow optimization.

[0361] The bounded weight criterion introduces diagnostic index difference Δβ and linear scaling ratio W of diagnostic index difference lim , which limits the variation range of weight, and its expression can be written as: W = 1 + W lim ·(Δβ) α (102)

[0362] wherein,

[0363] In the formula: Δβ is the difference degree of β and β ave ; W lim is the linear scaling ratio of Δβ; W lim,min is the linear scaling ratio of Δβ when β is less than β ave , -1≤W lim,min <0; W lim,max is the linear scaling ratio of Δβ when β is greater than β ave , 0<W lim,min ≤1; R is a constant for setting the calculation range [β min , β max ] of the normalized diagnostic index β, 0<R≤1; β min is the lower limit of the calculation range of β; β max is the upper limit of the calculation range of β.

[0364] By substituting formula (103)-(105) into formula (102), the flow update weight W is a bounded segmented function of β:

[0365] According to formula (106), the minimum value of weight W is 1+W lim,min , and the maximum value is 1+W lim,max , so when W lim,min takes -1, the flow of the well or well pair with poor performance (β<β min ) is reduced to 0; when W lim,max takes 1, the flow of the well or well pair with good performance (β>β max ) is increased to at most twice the original value.

[0366] Taking β ave =0.5 as an example, FIG. 28A)-FIG. 28D) show different a, R, W lim,min, W lim,max The relationship between the flow rate update weight W and the normalized flow field diagnostic indicator β. As shown in Fig. 28A), the weight index a influences the behavior of the W~β relationship curve in the interval [β min ,β max ]. With the increase of a, the curve moves up in [β min ,β ave ] and moves down in [β ave ,β max ]. As shown in Fig. 28B), the coefficient R influences the range of the non-limiting weight of the W~β relationship curve, and the length of the interval [β min ,β max ] increases with R, and at most increases to 2β ave . As shown in Fig. 28C), the scaling coefficient W lim,min influences the minimum weight of the W~β relationship curve, and the lower bound of the flow rate update weight W gradually tends to 1 in the process of W lim,min from -1 to 0. As shown in Fig. 28D), the scaling coefficient W lim,max influences the maximum weight of the W~β relationship curve, and the upper bound of the flow rate update weight W gradually tends to 2 in the process of W lim,max from 0 to 1. From Figs. 28A)-28D), it can be seen that the parameters a, R, W lim,min , W lim,max in the bounded weight criterion influence the calculation of the flow rate update weight, and further influence the optimization results of the injection and production flow rates. Generally speaking, the optimization effects of the streamline simulation optimization schemes obtained by setting different a, R, W lim,min , W lim,max values are different, but due to the streamline optimization process, the flow field diagnostic indicators (injection and production efficiency, remaining movable oil displacement efficiency, etc.) are balanced, and these optimal injection and production schemes can all obtain good oil increase and water control effects; if necessary, the variable parameters in the injection and production flow rate optimization criterion can be further optimized.

[0367] In summary, after the streamline simulation or the streamline tracking post-processing process of the finite difference simulation using the history-fitted reservoir numerical model, the four well-pair level flow field diagnostic indicators (movable remaining oil displacement efficiency, movable remaining oil recovery efficiency, injection efficiency, and recovery efficiency) and the four single-well level flow field diagnostic indicators can be used to evaluate the current water flooding effect and displacement balance of the injection and production well pair or single well, and further determine the direction of flow rate adjustment. For the injection and production well pair or single well, there are six flow rate update modes for each of them. Then, the injection and production flow rates of the well pair or single well can be optimized in combination with the injection and production flow rate optimization criterion (such as the unbounded weight criterion and the bounded weight criterion); and then the constraint processing is added.

[0368] Figure 29 outlines the basic flow of injection-production optimization of streamline simulation model at a single time step; injection-production optimization at current time step is based on the streamline information obtained from streamline simulation at last time step, while injection-production optimization at next time step is based on the updated injection-production flow rate obtained from injection-production optimization at current time step, thus, the whole streamline simulation optimization process is sequentially performed during the whole optimization period, only one complete forward reservoir simulation is needed to obtain the optimal injection-production scheme of all injection-production wells in the target area, and the number of wells involved in production optimization problem is not sensitive, which can effectively overcome the shortcomings of traditional production optimization methods based on simulator and optimization algorithm, such as too many simulation times and too high calculation cost, and the calculation efficiency is also higher than that of machine learning agent optimization method (agent optimization method based on simulator).

[0369] A specific embodiment (Example 3) is given below to illustrate the application of streamline simulation model alone for flow field diagnosis.

[0370] The reservoir porosity is 0.25; the number of grids is 23x23x3, the grid plane size is 30x30m 2 , and the grid longitudinal size is 15+10+10m; the production wells (PRO1, PRO2, PRO3, PRO4) are located at the center of the rectangular reservoir boundary, the four injection wells (INJ1, INJ2, INJ4, INJ5) are located at the four corners of the reservoir, and the one injection well (INJ3) is located at the center of the reservoir. The reservoir has a production history of 6902 days, the oil well is produced at a target oil production of 70m 3 / d, and the bottom hole pressure lower limit is set to 165 bar, and the water well is injected at a constant flow rate, and the bottom hole pressure upper limit is set to 600 bar.

[0371] The injection-production flow rate of the five injection wells and four production wells for the next 12 time steps (step length of 30d) is optimized by using the production optimization method based on streamline simulation model, and the bounded weight criterion parameters are set as: a = 1.5, R = 0.85, W lim,min = -1, W lim,max = 1, and the optimal injection-production scheme obtained by using injection efficiency and movable residual oil displacement efficiency for flow field diagnosis is shown in Figures 30 and 31, respectively.

[0372] Table 2 shows the comparison of the production indexes of the optimal injection-production scheme (i.e., development regulation scheme) obtained based on the streamline simulation model flow field diagnosis optimization and the benchmark scheme, from which it can be seen that the optimization schemes based on the streamline simulation model of the injection efficiency of well pairs and the movable residual oil displacement efficiency of well pairs both improve the water drive performance of the original scheme and the benchmark scheme to some extent, and the optimization scheme based on the streamline simulation of the movable residual oil displacement efficiency of well pairs obtains better oil increase and water control effect and higher net present value. This shows that increasing the injection rate of injection well INJ3 can displace more movable residual oil around it and improve the development effect.

[0373] Table 2 shows the comparison of the production indexes of the optimal injection-production scheme (i.e., development regulation scheme) obtained based on the streamline simulation model flow field diagnosis optimization and the benchmark scheme, from which it can be seen that the optimization schemes based on the streamline simulation model of the injection efficiency of well pairs and the movable residual oil displacement efficiency of well pairs both improve the water drive performance of the original scheme and the benchmark scheme to some extent, and the optimization scheme based on the streamline simulation of the movable residual oil displacement efficiency of well pairs obtains better oil increase and water control effect and higher net present value. This shows that increasing the injection rate of injection well INJ3 can displace more movable residual oil around it and improve the development effect. (Note: NPV represents net present value, N p represents cumulative oil production, W p represents cumulative water production, f w represents water cut.)

[0374] Sixth: numerical simulator surrogate model

[0375] For the numerical simulator surrogate model, the training and production optimization of the surrogate are coupled. The numerical simulator surrogate model does not need to accurately approximate all the real responses of the reservoir simulator at the beginning, but is updated by continuously "sampling" and constantly seeking valuable target function evaluation points (potential global optimal points) on the basis of the initial surrogate constructed by the initial experimental design, and then searches for the global optimal solution as quickly as possible under the limited number of function evaluations. Therefore, the numerical simulator surrogate model can be used for production optimization alone. The general process of production optimization based on the numerical simulator surrogate model in the field of reservoir production optimization is shown in FIG. 32.

[0376] Production optimization based on the numerical simulator surrogate model can actually be regarded as a special case of applying a general surrogate optimization algorithm to solve the problem of reservoir production optimization. At this time, the value of the objective function of the optimization problem needs to be obtained by running the numerical simulator surrogate model to predict the oil and water production performance under a certain injection-production setting, so as to complete the initial sampling, surrogate construction and updating. The surrogate optimization algorithm based on the numerical simulator surrogate model includes three core parts: initial experimental design, surrogate model selection and construction, and sampling strategy (sampling method); FIG. 33 shows the main process of the execution of the surrogate optimization algorithm through a simple one-dimensional optimization problem example. The general process of FIG. 32 is improved in the embodiments of the present application to obtain the process of production optimization based on the numerical simulator surrogate model.

[0377] In an embodiment, the following steps are used to perform production optimization based on the numerical simulator surrogate model:

[0378] determining a production optimization problem, the production optimization problem including decision variables, constraint conditions and an objective function;

[0379] based on the production optimization problem, determining a compromise sampling strategy, a type of the numerical simulator surrogate model and a sample set, the compromise sampling strategy taking a minimum of a response of the numerical simulator surrogate model and a maximum of a future search ability as dual objectives;

[0380] based on the sample set, constructing an initial numerical simulator surrogate model and taking the initial numerical simulator surrogate model as a current numerical simulator surrogate model, repeatedly performing the following steps until a termination condition is met, outputting an optimal decision variable when the termination condition is met, and obtaining a development control scheme corresponding to the optimal decision variable:

[0381] constructing a comprehensive evaluation function according to an evaluation index of the future search ability and an evaluation index of the response of the numerical simulator surrogate model and a compromise coefficient;

[0382] determining a final potential candidate point according to the comprehensive evaluation function;

[0383] in the current numerical simulator surrogate model, using an oil reservoir production optimization objective function evaluation subprogram to evaluate the final potential candidate point, wherein the oil reservoir production optimization objective function evaluation subprogram takes the decision variable as input;

[0384] adding the final potential candidate point to the sample set that has been evaluated by the function;

[0385] based on the sample set that has been evaluated by the function, updating the numerical simulator surrogate model;

[0386] taking the updated numerical simulator surrogate model as the current numerical simulator surrogate model.

[0387] In an embodiment, based on the production optimization problem, determining the sampling strategy, the type of the numerical simulator surrogate model and the sample set includes:

[0388] using an initial experimental design method to select a preset number of samples in a solution space of the production optimization problem for objective function evaluation, and generating a sample set required for constructing an initial numerical simulator surrogate model.

[0389] In an embodiment, after obtaining the initial potential candidate point set, further includes:

[0390] in the initial potential candidate point set, adding randomly sampled points near the current best feasible potential candidate point and Latin hypercube sampled points within a bounded constraint range.

[0391] In an embodiment, determining the final potential candidate point according to the comprehensive evaluation function includes:

[0392] constructing a constraint violation evaluation function according to the constraint conditions;

[0393] sampling the decision variables randomly to obtain an initial potential candidate point set under bounded constraints of the production optimization problem, wherein the bounded constraints are one of the constraint conditions;

[0394] screening feasible potential candidate points from the initial potential candidate point set based on the constraint violation evaluation function to obtain a feasible potential candidate point set;

[0395] selecting a feasible potential candidate point with a minimum comprehensive evaluation function as a final potential candidate point, wherein an evaluation index of future search capability in the comprehensive evaluation function adopts a minimum distance from a sample to a sample that has been evaluated by the objective function.

[0396] In an embodiment, screening feasible potential candidate points from the initial potential candidate point set based on the constraint violation evaluation function to obtain a feasible potential candidate point set comprises:

[0397] calculating a constraint violation evaluation function value;

[0398] excluding samples with a constraint violation evaluation function value greater than a constraint tolerance error in the initial potential candidate point set;

[0399] regarding the remaining set as a feasible potential candidate point set;

[0400] removing samples with a distance to a sample that has been evaluated by the objective function exceeding a tolerance distance from the feasible potential candidate point set.

[0401] Determining a sample set refers to using an initial experimental design method to select a preset number of samples in a solution space of the production optimization problem to evaluate the objective function, and generating a sample set required for constructing an initial numerical simulator proxy model to guide a subsequent sampling process to quickly cover a region where a global optimal solution is located, thereby completing fast optimization. Generally speaking, the final optimization effect of the proxy optimization algorithm should be less sensitive to the initial required sample set. In the present application, Latin hypercube sampling (LHS) and symmetric Latin hypercube sampling (SLHS) are mainly used for initial experimental design to generate a sample set required for constructing an initial numerical simulator proxy model.

[0402] For the type of numerical simulator agent model in the agent optimization algorithm, since it needs to be iteratively constructed and dynamically updated, its construction process cannot be too cumbersome and time-consuming. The present application selects the radial basis function interpolation model as the preset item, and can also use multivariate adaptive regression splines (MARS), polynomial regression, Kriging model, Gaussian process regression, and neural network models such as precise radial basis network (RBE), generalized regression neural network (GRNN), and radial basis network (RB). Most of these proxy models can be implemented through built-in functions in Matlab, open source libraries in Python, etc. Here, only the radial basis function interpolation model that is more suitable for the agent optimization algorithm is introduced in detail.

[0403] For samples that have been evaluated by functions, the estimated value of the radial basis function (RBF) interpolation model is the true value, which is an important advantage of the RBF interpolation model, and its expression can be written as: p(x)=p·c (108) p=[x T ,1] (109) x=[x1,x2,...,x n ] T (110) c=[c1,c2,...,c n ,c0] T (111) λ=[λ1,λ2,...,λ m ] T (112)

[0404] In the formula: x is a decision variable column vector; n is the dimension of the decision variable column vector (problem dimension); S surrogate (x) is the response value of the proxy model at x; m is the number of samples that have been evaluated by functions; x i is the i-th sample that has been evaluated by functions; RBF(·) is the radial basis function model; ||…|| is the Euclidean norm / distance (or vector modulus); λ i is the radial basis function coefficient of the i-th sample that has been evaluated by functions; p(x) is the linear tail of the radial basis function interpolation model; p is a row vector constructed by adding an element 1 to the decision variable vector; c is the coefficient vector of the linear tail of the radial basis function interpolation model.

[0405] If there are m samples that have been evaluated by functions, i.e., the objective function values of these m samples are known, then the relationship between the decision variable x and the objective function f un (x) needs to be fitted using formula (107). To solve the radial basis function coefficient vector λ and the coefficient vector c of the linear tail, the following linear equation system can be constructed:

[0406] wherein: R ij = R ji BF(||x i -x j ||), (1≤i,j≤m) (115) Fun= [f un (x1), f un (x2), f un (x3), …, f un (x m )] T (118)

[0407] wherein: R ij is the value of the radial basis function; Phi is the radial basis function matrix; x i,j is the jth element of the ith sample on which the function has been evaluated; P is the constructed sample matrix; λ is the radial basis function coefficient vector; c is the coefficient vector of the linear tail of the radial basis function interpolation model; Fun is the objective function vector corresponding to the samples on which the objective function has been evaluated; f un (x) is the objective function corresponding to the decision variable column vector x.

[0408] When the rank of the matrix P is n+1 (m≥n+1), the coefficient matrix in the linear equation system (113) is invertible, denoted as A rbf , then the to-be-solved model parameters of the radial basis function interpolation model can be expressed as: C=A rbf \B (119)

[0409] wherein,

[0410] According to equations (119) and (120), the to-be-solved matrix C of the RBF interpolation model as a surrogate model can be easily solved, and then the radial basis function coefficient vector λ and the linear tail coefficient vector c can also be obtained. It should be noted that there are 6 commonly used forms of radial basis functions RBF(·) to choose from, as shown in Table 3. For the surface spline function, when k=1, it becomes a linear function; when k=2 and 3, it is converted into a thin plate spline function and a cubic spline function, respectively.

[0411] The sampling strategy (or sampling method) is the most important component of the surrogate optimization, which directly determines the nature of the subsequent selected function evaluation points (potential candidate points) and the update of the surrogate model, and affects the search speed of the global optimal solution. In addition, for optimization problems containing complex linear and nonlinear constraints, the constraint handling method can be coupled in the sampling strategy to ensure that the subsequent function evaluation points are all feasible points, thereby accelerating the search for the global best feasible solution. The sampling strategy of balancing the response of the surrogate model and the future search ability is adopted to select the potential candidate points (the next batch of valuable function evaluation points), and the processing method of the complex constraint conditions is considered.

[0412] On the one hand, the optimization problem requires to find the optimal decision variable as much as possible to minimize the objective function, so the subsequent sampling points should make the response value of the surrogate model as small as possible; on the other hand, in order to avoid sampling into a local optimal solution, as many unknown solution spaces as possible need to be sampled; therefore, a compromise sampling strategy with the dual objectives of minimizing the response of the surrogate model and maximizing the future search ability is proposed, and a comprehensive evaluation function is constructed to reconcile the contradiction between the two objectives, the key of which is:

[0413] Table 3 Main types of radial basis functions

[0414] Key 1: Define the constraint violation evaluation function.

[0415] In the description of the optimization problem, the constraint violation evaluation function can be constructed according to the constraint conditions, which lays the foundation for the selection of subsequent samples. As shown in equation (123), a general optimization problem (written in the form of minimization) may contain bounded constraints, linear inequalities, linear equations, nonlinear inequalities, nonlinear equations and other constraint conditions:

[0416] In the description of the optimization problem, the constraint violation evaluation function can be constructed according to the constraint conditions, which lays the foundation for the selection of subsequent samples. As shown in equation (123), a general optimization problem (written in the form of minimization) may contain bounded constraints, linear inequalities, linear equations, nonlinear inequalities, nonlinear equations and other constraint conditions: un (x) is the objective function; x lb is the lower bound of the decision variable composed of the column vector; x ub is the upper bound of the decision variable composed of the column vector; A is the coefficient matrix of the linear inequality constraint; b is the column vector composed of the linear inequality constraint constant; A eq is the coefficient matrix of the linear equation constraint; b eq is the column vector composed of the linear equation constraint constant; non ineq is the column vector composed of the nonlinear inequality constraint; non eq is the column vector composed of the nonlinear equation constraint.

[0417] For any decision variable column vector x (sample), the constraint violation evaluation function is defined as:

[0418] where f penalty (x) is the constraint violation evaluation function (or penalty function) corresponding to the decision variable column vector x.

[0419] Since the bounded constraints can be directly satisfied when selecting initial samples, formula (124) only defines the penalty terms corresponding to linear inequality, linear equality, nonlinear inequality, and nonlinear equality constraints. It can be found that f penalty is positive when x violates any of the constraints, and is 0 when x satisfies all the constraints. Therefore, whether x is a feasible point can be determined according to whether f penalty (x) is less than a set constraint error.

[0420] Key 2: Determine initial potential candidate points.

[0421] Potential candidate points refer to samples to be selected as the next batch of samples for target function evaluation. Random sampling of decision variables under the bounded constraints of the production optimization problem can obtain an initial potential candidate point set. On this basis, the present embodiment further adds random sampling points near the current best feasible point, Latin hypercube sampling points within the bounded constraint range, and the like as supplements to the initial potential candidate point set. Since the speed of the proxy model for evaluating a sample is very fast, the number of samples obtained by each sampling method described above can be set to a relatively large number (for example, set to 50, 100, or 200 times the dimension of the optimization problem).

[0422] Key 3: Screen feasible potential candidate points.

[0423] The samples in the initial potential candidate point set obtained in the previous step only satisfy the bounded constraints of the production optimization problem, and may violate other constraints. Therefore, samples that do not satisfy the constraints need to be removed according to the defined constraint violation evaluation function (penalty function). The set remaining after excluding samples with a penalty function greater than a constraint tolerance error from the initial potential candidate point set is referred to as a feasible potential candidate point set. It should be noted that some points in the feasible potential candidate point set may be too close to known points (samples that have been evaluated for the target function), and it is not meaningful to evaluate the target function for these points as final potential candidate points. Therefore, samples in the feasible potential candidate point set that are more than a tolerance distance away from known points also need to be removed.

[0424] Key 4: Determine final potential candidate points.

[0425] The feasible potential candidate points obtained in the previous step satisfy all the constraints, but which points have the greatest value for function evaluation still need to be comprehensively evaluated according to their corresponding proxy model responses and future search capabilities. The minimum distance Ddistance As an evaluation index of "future search ability", the smaller the value is, the closer the distance to the sample that has been evaluated by the objective function is, and the weaker the search ability for the region that has not been evaluated by the objective function is; for the response of the surrogate model, the S surrogate is directly used. distance In order to avoid D surrogate , the two physical quantities are normalized, that is, the S distance is normalized.

[0426] Or

[0427] In the formula, S is the normalized value of the response of the surrogate model at the feasible potential candidate point x; is the normalized value of the minimum distance D surrogate from the feasible potential candidate point x to the sample that has been evaluated by the objective function; S distance is the set of surrogate response values corresponding to the set of feasible potential candidate points; D balance is the set of minimum distances from the feasible potential candidate points to the sample that has been evaluated by the objective function; max(·) is the maximum value operator; and min(·) is the minimum value operator.

[0428] A compromise coefficient C cv is introduced to control the balance between the minimum surrogate model response and the maximum future search ability, and a comprehensive evaluation function f cv is constructed to evaluate the value of the potential subsequent point:

[0429] In the formula, f balance is the comprehensive evaluation function; and C balance is the compromise coefficient for balancing the surrogate model response and the future search ability, and takes a value in the interval [0, 1].

[0430] The formula (127) comprehensively considers the dual objectives of the minimum surrogate model response and the maximum future search ability, and the value of C balance can reflect different search strategies; when C cv is close to 0, the distance term plays a dominant role, so that the feasible potential candidate point with the minimum comprehensive evaluation function f balance has strong future search (global optimization) ability; when C cv is close to 1, the surrogate term plays a dominant role, so that the feasible potential candidate point with the minimum comprehensive evaluation function f balance has strong local optimization ability. The "compromise" in the compromise sampling strategy proposed in the present application is reflected in the selection of the compromise coefficient C balance , which can be determined in the following three ways:

[0431] 1) Set to a fixed value, such as 0.4, 0.5, 0.6, etc.; if C balance = 1, the compromise sampling strategy degenerates into the agent minimization sampling strategy;

[0432] 2) Sequential iteration in a fixed array; for example, cyclically taking values at nine-equal points in [0, 1], ten-equal points in [0.1, 0.95], or using other fixed sequences;

[0433] 3) In addition to sequential iteration in a fixed array, the order of taking values is also adjusted in real time according to the sampling effect; for example, when the current optimal function value obtained after 3 consecutive samplings is lower than the previous optimal value, it indicates that the sampling is successful, and the sequence can be jumped to the later sequence, increasing C balance to increase the weight of the agent, thereby accelerating the search for a local optimal solution; if the function value of 3 consecutive samplings is higher than the previous optimal value, it indicates that the sampling fails, and the sequence can be jumped to the earlier sequence, decreasing C balance to increase the weight of the distance term, thereby accelerating the search for a global optimal solution.

[0434] After determining the compromise coefficient C balance , a batch (or one) of feasible potential candidate points with the minimum comprehensive evaluation function f cv is selected as the final potential candidate point, and then the final potential candidate point is subjected to objective function evaluation operation (i.e., evaluating the value of the final potential candidate point), and the final potential candidate point is added to the sample set that has been subjected to function evaluation for updating the numerical simulator agent model, and then the updated numerical simulator agent model is combined with the compromise sampling strategy to sample the next final potential candidate point, and the process is repeated until the maximum number of function evaluations or other termination conditions are reached.

[0435] In summary, the numerical simulator agent model optimization algorithm proposed in the embodiments of the present application uses Latin hypercube design and other initial test design methods to determine the initial sample set, uses a cubic spline function interpolation model as the preset type of the numerical simulator agent model, and selects the compromise sampling strategy (C balance sequential dynamic updating) as the preset sampling strategy, which can be expanded according to needs; this algorithm can solve optimization problems containing various complex constraints, and is particularly suitable for time-consuming constraint optimization problems such as oil reservoir production optimization.

[0436] A specific embodiment is given below to illustrate the specific application of production optimization based on a numerical simulator agent model.

[0437] Now consider a 5 injection 4 production heterogeneous reservoir model, keep all settings in example 3 unchanged, the permeability distribution in x, y, z directions is shown in figure 34. The production history of this reservoir for 6902 days is shown in figure 35, and the injection and production schedule of all wells is shown in table 4.

[0438] The injection strategy of injection wells in the production history of example 3 is to divide the total injection amount equally, and the production wells are produced at a constant oil rate in the early stage, and then converted to constant flow pressure production due to the low oil production rate, which needs to be improved in the adaptability of the permeability heterogeneity of the reservoir; therefore, the production schedule for the next 12 time steps (30d step length) is optimized, the upper limit of the bottom hole pressure of the injection well and the lower limit of the bottom hole pressure of the production well are kept unchanged, the lower limit of the single well injection rate is set to 5m 3 / d (if the optimized injection rate is lower than this value, the well is shut down), the upper limit of the single well injection rate is 297m 3 / d, the lower limit of the single well liquid production rate is 10m 3 / d, the upper limit of the single well liquid production rate is 287.1m 3 / d, the upper limit of the water cut of the production well is 0.98, the upper limit of the injection rate of the reservoir is 300m 3 / d, the upper limit of the liquid production rate of the reservoir is 290m 3 / d.

[0439] Table 4 injection and production schedule of example 3

[0440] Under the condition of meeting all the above constraints, the injection and production scheme obtained by dividing the total injection amount of the water well and the total liquid production rate of the well group to control the oil well is taken as the original scheme, and the injection and production scheme obtained by controlling the injection and production of the well group (development control scheme) is taken as the benchmark scheme. The injection and production flow rate is optimized by using the agent optimization algorithm proposed in the foregoing according to the flow shown in figure 32, and the injection and production flow rate setting of the optimized scheme and the benchmark scheme is shown in figure 36; figure 37 shows the change of the production optimization objective function (NPV, the income converted to the initial stage of production) with the number of simulator runs.

[0441] As can be seen from figure 37, after about 600 times of objective function evaluation, the increase amplitude of net present value (NPV) is not obvious, at which time the running can be terminated. In this example, the dimension of the production optimization problem is 12x(5+4)=108, and it is recommended that the maximum number of function evaluations be set to more than 6 times the dimension of the optimization problem. Table 5 gives the comparison of the production indicators of the injection and production scheme obtained by agent optimization, the original scheme and the benchmark scheme at the end of the optimization period, from which it can be seen that the optimized scheme makes the reservoir dynamic indicators better than the original scheme and the benchmark scheme through 12 months of injection and production adjustment, the oil production is increased by 866m 3 , and the water production is reduced by 7891m 3, the end-of-period water cut is reduced by 3.4 percent, and the net present value is increased by 7.5641 million dollars. This shows that the injection strategy of the original equal injection amount is not suitable for the injection-production control of the reservoir in Example 3, and the injection-production scheme under the control of the well group is based on the injection capacity or production capacity of the single well to allocate the flow, and the effect of increasing oil and controlling water is not good.

[0442] Table 5 Comparison of optimization results of the optimal injection-production scheme obtained by the numerical simulator agent model and the benchmark scheme

[0443] The seventh kind: streamline simulation agent model

[0444] Let the number of injection wells to be optimized in the reservoir production optimization problem be N inj , the number of production wells to be optimized be N pro , and the optimization time step be N t , then the dimension of the optimization variable is (N inj +N pro )×N t ; when the agent optimization method is used to solve such time-consuming production optimization problems, as the number of wells increases, the maximum number of function evaluations that need to be set to ensure optimization accuracy should also increase. Since numerical simulation generally takes a long time, the time for calling the numerical simulator for operation will increase rapidly, so the agent optimization is generally more suitable for low-dimensional production optimization problems. As mentioned earlier, the production optimization based on the streamline simulation flow field diagnosis is a sequential optimization process, which is not sensitive to the number of injection-production wells, and is more suitable for high-dimensional optimization problems; however, there is still room for further optimization in the variable parameters in the injection-production flow optimization criteria (such as α in the unbounded weight criterion and α, R, W lim,min , W lim,max , etc. in the bounded weight criterion), so the streamline simulation agent model can be used to optimize these parameters, that is, to run the optimization process of the streamline simulation under the framework of agent optimization; at this time, the dimension of the optimization problem is less than 5, and the production optimization problem involving a large number of injection-production wells and the time-consuming problem of calculation are effectively solved. The main steps of production optimization based on the streamline simulation agent model are summarized in FIG. 38. FIG. 39 details the basic elements of production optimization based on the streamline simulation agent model.

[0445] In an embodiment, the following steps are used to perform production optimization based on the streamline simulation agent model:

[0446] Determine the production optimization problem, which includes decision variables, constraint conditions and objective functions;

[0447] Based on the production optimization problem, determine the compromise sampling strategy, the type of streamline simulation agent model and the sample set, and the compromise sampling strategy takes the minimum response of the streamline simulation agent model and the maximum future search ability as the dual objectives;

[0448] Based on the sample set, an initial streamline simulation agent model is constructed, and as a current streamline simulation agent model, the following steps are repeatedly executed until a termination condition is met, and the optimal decision variable when the condition is met is output, and the development control scheme corresponding to the optimal decision variable is obtained:

[0449] According to the evaluation index of future search capability and the evaluation index of the response of the streamline simulation agent model, and the compromise coefficient, a comprehensive evaluation function is constructed;

[0450] According to the comprehensive evaluation function, a final potential candidate point is determined;

[0451] In the current streamline simulation agent model, a streamline optimization sub-function is used to evaluate the objective function of the final potential candidate point, wherein the streamline optimization sub-function takes the variable parameter of the flow optimization criterion when production optimization is performed based on the flow field diagnosis of the streamline simulation model as input, and calculates the objective function value based on the simulation result of the development control scheme obtained when production optimization is performed based on the flow field diagnosis of the streamline simulation model;

[0452] The final potential candidate point is added to the sample set that has been subjected to function evaluation;

[0453] Based on the current sample set that has been subjected to function evaluation, the streamline simulation agent model is updated;

[0454] The updated streamline simulation agent model is taken as the current streamline simulation agent model.

[0455] A specific embodiment is given below to illustrate the specific application of production optimization based on the streamline simulation agent model.

[0456] A production optimization method based on the streamline simulation agent model is used to compile a corresponding program to optimize the injection and production flow rates of the 9 wells in Example 3 for 12 time steps in the future. Take W lim,min =-1, W lim,max =1, take α and R as optimization variables in the bounded weight criterion, and Fig. 40 shows the calculation process of the streamline simulation agent optimization algorithm when the flow field diagnosis index is the injection efficiency of the well, and the optimal decision variables obtained are α=0.262341, R=0.822361, and the optimal injection and production scheme is shown in Fig. 41.

[0457] Take α, R, W lim,min , W lim,max as optimization variables in the bounded weight criterion, and the calculation process of the streamline simulation agent model optimization method when the flow field diagnosis index is the displacement efficiency of the movable residual oil of the well is shown in Fig. 42, and the optimal decision variables obtained are: α=0.304290, R=0.854140, W lim,min= -0.765057, W lim,max = 0.854954, and the optimal injection-production scheme is shown in FIG. 43.

[0458] As can be seen from FIGS. 41 and 43, the optimization mode of the two flow line simulation proxy models taking injection efficiency and movable residual oil displacement efficiency as flow field diagnosis indicators both obtain better net present value by increasing the injection rate of injection well INJ3; in FIG. 41, the flow rates allocated to injection wells INJ2 and INJ5 are lower than 5 m 3 / d, violating the minimum injection rate constraint, so these two injection wells are shut down in the subsequent optimization time steps, and the liquid production rate of production well POR4 is gradually reduced to a lower level, and the liquid production rates of PRO2 and PRO3 are increased. In FIG. 43, the injection rate of INJ3 is obviously higher than that of other injection wells from the second optimization time step, but the liquid production rates of the production wells are relatively balanced, and the liquid production rate of POR4 is slightly reduced. Therefore, attention should be paid to the consistency of the optimization schemes based on the two flow line simulation proxy models of well pair injection efficiency and well pair movable residual oil displacement efficiency in the flow rate adjustment of some injection-production wells, and these development regulation schemes can be referred to in actual reservoir injection-production adjustment.

[0459] Table 6 shows the comparison of production indicators of the development regulation schemes obtained by the flow line simulation proxy models with the base scheme, from which it can be seen that, compared with the original flow line simulation model, both the flow line simulation proxy models based on well pair injection efficiency and well pair movable residual oil displacement efficiency further improve the economic benefit (NPV), and the latter has a smaller improvement range; for example 3, the flow line simulation proxy optimization based on well pair movable residual oil displacement efficiency obtains the highest net present value and cumulative oil production, the flow line simulation proxy optimization based on well pair injection efficiency has the lowest water cut at the end of the optimization period, and the cumulative water production of the proxy optimization scheme is the smallest. All development regulation schemes are superior to the injection-production settings of the original scheme and the base scheme. Considering the efficiency, accuracy of optimization calculation and the pros and cons of optimization results, the injection-production scheme corresponding to the flow line simulation proxy model based on well pair movable residual oil displacement efficiency is the optimal development regulation scheme.

[0460] Table 6 Comparison of production indicators of development regulation schemes obtained by flow line simulation proxy models with base scheme

[0461] In the embodiments of the present application, the well interconnection relationship needs to be calibrated for multiple production optimization models, so that the first type of prediction model provides an initial estimate of the well interconnection relationship for the second type of prediction model, and the second type of prediction model calibrates the model parameters of the first type of prediction model.

[0462] After the history matching process, the streamline simulation (SL) model can directly output the stream line allocation factor, and the finite difference (FD) model can calculate the injection-production flow allocation relationship through stream line tracking. These SL-based connectivity relationships can in turn guide the adjustment of the connectivity coefficient of the time-varying resistance-capacitance model and the initial conductivity of the connected element simulation model, so that the flow allocation relationship obtained by these simplified physical models is basically consistent with the stream line results derived from the SL model or the FD model, thereby reducing the multi-solution of these simplified physical models. Therefore, the simplified physical model uses its speed advantage to provide an initial estimate of the interwell connectivity relationship for the reservoir numerical (SL or FD) model (although this estimate may not be decisive for the fitting process of the SL or FD model), and the reservoir numerical model uses its credibility advantage to calibrate the model parameters affecting the flow allocation in the simplified physical model after fitting. In this way, the simplified physical model and the reservoir numerical (SL or FD) model can be considered consistent.

[0463] Next, intelligent production optimization is implemented. The simplified physical model after the interwell connectivity calibration improves the shortcomings of strong history matching multi-solution and weak long-term prediction ability inherent in the simplified physical model, and the injection rate, liquid production rate (or bottom hole pressure) optimization settings obtained from the production optimization process can provide development control schemes for short-term real-time optimization needs; while the numerical simulator proxy model, the flow field diagnosis optimization method based on the streamline simulation model, and the streamline simulation proxy model can provide development control schemes for medium and long-term optimization needs; these production optimization processes should all impose constraints on individual wells and well groups set by the updated reservoir dynamic understanding. It should be noted that when multiple production optimization methods are available, the optimal injection-production schemes obtained by using CRM, LUSM, numerical simulator proxy model, various optimization methods based on different flow field diagnosis indicators, and various streamline simulation proxy models can be compared, and actual operation should focus on injection-production wells with consistent optimization adjustment directions (such as injection well INJ3 in Example 3). These injection-production optimization results obtained from the prediction model can be used as theoretical model understanding to assist the update of reservoir dynamic understanding. At this point, the reservoir dynamic understanding, the first type of prediction model (simplified physical model), and the second type of prediction model (non-proxy model in the reservoir numerical model) realize their coupled interaction (Figure 44).

[0464] In step 106, the adjustment effect feedback and strategy update are performed.

[0465] After determining the development control scheme and implementing the development control, the adjusted wells are monitored in real time, the control effect is evaluated, and the injection-production adjustment strategy is updated and summarized, and the injection-production optimization method based on the reservoir dynamic understanding is improved. Then, return to step 101 to repeat a new round of comprehensive production optimization steps.

[0466] The above comprehensive production optimization workflow reveals a comprehensive optimization concept of benign coupling interaction between each optimization model (reservoir dynamic understanding, simplified physical model, reservoir numerical model). The "comprehensive optimization" is reflected in the use of production optimization based on reservoir dynamic understanding, the first type of prediction model, and the second type of prediction model to meet the development regulation and control requirements of different frequencies such as very short term, short term, and medium and long term, respectively. The "coupling interaction" is reflected in the constraints imposed by reservoir dynamic understanding on the first type of prediction model and the second type of prediction model when performing production optimization, the supplement of the optimization results of the first type of prediction model and the second type of prediction model to the reservoir dynamic understanding, the auxiliary role of the connected relationship obtained by the history matching of the first type of prediction model as an initial estimate to the history matching process relied on by the second type of prediction model, and the calibration role of the flow line distribution coefficient obtained by the second type of prediction model to the history matching parameters of the first type of prediction model. According to the comprehensive production optimization workflow, a series of intelligent reservoir closed-loop production optimization technologies can be gradually developed and formed, which have the advantages of complementary advantages, real-time and long-term consideration, and physical meaning and data-driven, and are suitable for different optimization requirements.

[0467] The following gives a specific embodiment again to illustrate the specific application of the intelligent comprehensive production optimization of the reservoir.

[0468] A certain actual reservoir model (Example 4) has a production history of nearly 12 years. The reservoir model currently has 5 water injection wells and 9 oil production wells. The injection wells are all in normal operation, but only 3 oil production wells (P1, P2, and P3) are currently in production, and the remaining production wells (P4-P9) are shut down due to insufficient pressure or pump failure. The history matching results of TV-CRM show that only injection wells I1 and I3 have a direct effect on the surrounding production wells, and the remaining injection wells mainly act on the formation. Referring to the connectivity evaluation results of TV-CRM, after history matching of the non-proxy model in the reservoir numerical model, the simplified physical model is calibrated, and finally the injection-production correspondence of each well of Example 4 is obtained as shown in FIG. 45A) and FIG. 45B). FIG. 45A) is the interwell connectivity coefficient inverted by TV-CRM in the embodiment of the application; FIG. 45B) is the injection flow distribution coefficient derived from the flow line simulation at the end of the history matching period, wherein the length of the high of the isosceles triangle between the water injection well (triangle) and the oil production well (circle) represents the strength of the interwell connectivity or the size of the injection distribution coefficient.

[0469] According to the reservoir understanding obtained from the main control factor analysis, periodic dynamic analysis, and node system analysis, the upper limit of the single-well injection rate of the reservoir is 560 m 3 / d, and the lower limit is 0 m 3 / d, the upper limit of the single-well liquid production rate is 155 m 3 / d, and the lower limit is 0 m 3 / d. Given the optimization example effect of the foregoing, an optimization method based on the streamline simulation proxy model is used to optimize the injection-production well flow of Example 4; the optimization period is also set to 12 months, each time step is 30d, and the upper bound constraint of the total injection volume of the reservoir is set to 594m 3 / d, and the upper bound constraint of the total liquid production volume of the reservoir is set to 162m 3 / d, and the injection-production system controlled by the well group is taken as the benchmark scheme.

[0470] Since the grid number (184.3344x10 4 ) of the model is large, the simulation takes a long time, and here only a and R in the bounded weight criterion are optimized, and the maximum number of numerical simulation runs (the number of times of evaluating the objective function) is not more than 50. FIG. 46 and FIG. 47 respectively show the optimization calculation process of the streamline simulation proxy model based on the injection efficiency of the well pair and the displacement efficiency of the movable residual oil of the well pair.

[0471] FIG. 48, FIG. 49 and FIG. 50 respectively show the injection-production flow allocation plan of the benchmark scheme of Example 4 under the flow control of the well group, the development regulation scheme of the streamline simulation proxy model based on the injection efficiency of the well pair (a = 0.837076, R = 0.091365), and the development regulation scheme of the streamline simulation proxy model based on the displacement efficiency of the movable residual oil of the well pair (a = 1.037813, R = 0.597812), in which the first two time steps are the last two time steps of the history matching, and time steps 3-14 represent the optimization period; and the comparison of the production indexes of each development regulation scheme is shown in Table 7. As can be seen from the table, since the injection capacity of I3 is strong, the injection volume of I3 is set to the maximum in the benchmark scheme, and the flow allocated to the remaining injection wells is small, and for the same reason, the liquid volume of production well P3 is maintained at a high value in the benchmark scheme, and the liquid volume of other oil wells is low, and this injection-production setting does not change significantly in the entire optimization period. The development regulation scheme of the streamline proxy model based on the injection efficiency is to first decrease and then increase the injection volume of I3, and first increase and then decrease the injection volume of I1, and the liquid volume is gradually dominated by P2, and the liquid volume of P1 and P3 gradually decreases. For the development regulation scheme of the streamline proxy model based on the displacement efficiency of the movable residual oil, the regulation direction of the liquid volume is consistent with that of the streamline proxy optimization scheme based on the injection efficiency, the liquid volume is dominated by P2 in the later period, and the liquid volume of P1 and P3 is very small, and the difference is that this scheme gradually reduces the injection volume of I3, and gradually increases the injection volume of I2, and first increases and then decreases the injection volume of I4, and first increases and then decreases the injection volume of I5, and finally I2, I5 and I4 jointly dominate the water drive process. The change of the above injection-production system with time reflects the difference in the optimization logic of different flow field diagnosis indexes, and it needs to be noted that the development regulation schemes of the two streamline proxy models show a consistent regulation idea for the liquid volume, and therefore this optimization result can be referred to in the injection-production optimization method based on dynamic understanding, and the liquid volume of oil well P2 can be appropriately increased.

[0472] Table 7 Comparison of production indexes of each development regulation scheme of Example 4

[0473] The intelligent comprehensive production optimization method and workflow (and various production optimization models or methods attached thereto) first proposed in the present application are applied to injection-production adjustment of various types of water drive reservoirs, and good field application effects are preliminarily obtained, and the effects of oil increase and water control are relatively obvious. The program corresponding to the intelligent comprehensive production optimization method can quickly obtain development regulation schemes with different adjustment frequencies (daily, weekly, monthly) within a few minutes to a few hours. The production optimization based on the reservoir numerical model also needs to develop an automatic history matching technology to further speed up the iteration speed of the optimization closed loop. The comprehensive production optimization method deeply integrates reservoir dynamic understanding, the first type of prediction model (simplified physical model (CRM, LUSM, etc.)) and the second type of prediction model under the closed-loop reservoir management framework, realizes the interaction and linkage between reservoir understanding optimization and intelligent optimization methods and between intelligent optimization methods, can assist real-time injection-production regulation and control and the formulation of medium and long-term development adjustment decisions in water drive oilfield sites, can greatly save the economic and time costs consumed by traditional production optimization, and can improve the efficiency, accuracy and intelligent level of reservoir development regulation.

[0474] The present application also proposes a reservoir comprehensive production optimization device, which has a similar principle to the reservoir comprehensive production optimization method, and will not be described here.

[0475] FIG. 51 is a schematic diagram of a reservoir comprehensive production optimization device in an embodiment of the present application, which includes:

[0476] The production history data obtaining module 5101 is configured to obtain production history data of a target reservoir.

[0477] The reservoir dynamic understanding updating module 5102 is configured to update reservoir dynamic understanding according to the production history data, including regional level understanding, single well level understanding and production optimization model understanding.

[0478] The history matching module 5103 is configured to perform cyclic history matching on model parameters of a plurality of production optimization models based on the production history data until optimal model parameters are obtained, the plurality of production optimization models including the first type of prediction model and the second type of prediction model, by the following steps: performing history matching on the model parameters of each production optimization model, and calibrating interwell connectivity of the plurality of production optimization models, so that the first type of prediction model provides an initial estimation of the interwell connectivity for the second type of prediction model, and the second type of prediction model calibrates the model parameters of the first type of prediction model.

[0479] The comprehensive production optimization module 5104 is implemented to determine the constraint condition corresponding to each production optimization model according to the updated reservoir dynamic understanding, to solve the objective function corresponding to each production optimization model, to obtain the optimal decision variable corresponding to each production optimization model, and to obtain the development control scheme corresponding to each production optimization model; and to obtain the comprehensive development control scheme of the reservoir according to the development control scheme corresponding to each production optimization model, the comprehensive development control scheme of the reservoir being used to guide the development operation of the target reservoir.

[0480] The adjustment effect feedback and updating module 5105 is implemented to update the reservoir dynamic understanding according to the control effect obtained through real-time monitoring of the development operation.

[0481] In an embodiment, the reservoir dynamic understanding updating module is specifically used to:

[0482] According to the production history data, the development phase and regional development characteristics are divided in the reservoir development process, and a plurality of main control factors affecting the reservoir and regional development effect and producing the observed production dynamic are analyzed from a plurality of development aspects;

[0483] The reservoirs are classified according to different main control factors;

[0484] For the reservoir classification conforming to the regional overall law, the regional level understanding is obtained;

[0485] For the reservoir classification not conforming to the regional overall law, the single-well level understanding is obtained through the periodic reservoir dynamic analysis and the node system analysis;

[0486] The theoretical model understanding is obtained according to the production optimization model;

[0487] The regional level understanding is corrected through the single-well level understanding, and the single-well level understanding is verified through the regional level understanding;

[0488] The regional or single-well level understanding is updated through the theoretical model understanding, and the theoretical model understanding is constrained through the single-well level understanding;

[0489] The reservoir dynamic understanding is updated according to the regional level understanding, the single-well level understanding and the production optimization model understanding;

[0490] The reasonable development technical strategy is determined according to the updated reservoir dynamic understanding;

[0491] The constraint condition corresponding to each production optimization model is determined according to the updated reservoir dynamic understanding, including:

[0492] The constraint condition corresponding to each production optimization model is determined according to the reasonable development technical strategy.

[0493] In an embodiment, the reservoir dynamic understanding updating module is specifically used to:

[0494] According to the rational development technical strategy, the reservoir state is diagnosed, and the problem well is determined;

[0495] The development regulation scheme is generated for the problem well, and the development regulation scheme is used for adjusting operation of the problem well;

[0496] According to the real-time monitoring and effect feedback of the adjusting operation, the reservoir dynamic understanding and the development regulation scheme are updated.

[0497] In an embodiment, the first type of prediction model is a capacitance resistance model, a connection unit model or a flow network model.

[0498] In an embodiment, the history fitting module is further used for:

[0499] On the basis of the approximation of the oil-water two-phase material balance equation and the deliverability equation, the flow control equation and the saturation control equation are determined;

[0500] According to the flow control equation and the saturation control equation, the time-varying capacitance resistance model is established.

[0501] In an embodiment, the flow control equation of the time-varying capacitance resistance model is as follows:

[0502] Wherein, q t (t k+1 ) and q t (t k ) are the underground liquid production at t k+1 and t k , respectively, J(t k+1 ) and J(t k ) are the productivity index at t k+1 and t k , respectively, τ(t k+1 ) and τ(t k ) are the time constant at t k+1 and t k , respectively, w(t k+1 ) is the ground injection at t k+1 , B w is the water phase volume coefficient, p wf (t k+1 ) and p wf (t k ) are the bottom hole pressure at t k+1 and t k , respectively;

[0503] The saturation control equation of the time-varying capacitance resistance model is as follows:

[0504] Wherein, S w (tk+1 ) and S w (t k ) are the water saturation of the production well control volume at time t k+1 and t k , V p is the pore volume of the production well control volume, w is the surface injection rate, q w is the surface water production rate, C w is the water phase compressibility, C φ is the rock compressibility, C t is the overall compressibility, q t is the subsurface liquid production rate, S o (t k+1 ) and S o (t k ) are the oil saturation of the production well control volume at time t k+1 and t k , q o is the surface oil production rate, B o is the oil phase volume factor, C o is the oil compressibility.

[0505] In an embodiment, the history matching module is specifically configured to:

[0506] when the production optimization model is a time-varying capacitance resistance model, determining the model parameters to be the connectivity factor, the control volume pore volume, the deliverability index and the saturation;

[0507] determining the initial values of the model parameters;

[0508] determining the evolution process of the model parameters;

[0509] based on the initial values and the evolution process of the model parameters, iteratively solving a constraint optimization problem constituted by the model errors to determine the model parameters when the model errors are less than a preset threshold as the optimized model parameters.

[0510] In an embodiment, the model errors are as follows:

[0511] wherein, N pro is the number of production wells; N inj is the number of injection wells; error is the model estimation error of the subsurface liquid production rate of the production well at a certain time; N tim is the number of production history data points; E sum is the total sum of the subsurface liquid production rate estimation errors of all production wells at all time sequences, is the calculated value of the subsurface liquid production rate of the jth production well at k+1 time, is the true value of the subsurface liquid production rate of the jth production well at k+1 time;

[0512] The constrained optimization problem is as follows:

[0513] where S wc is the irreducible water saturation; S or is the residual oil saturation; V pt is the total pore volume of the production system; f ij is the connectivity coefficient of the injection well i to the production well j, is the average water saturation in the control volume of the production well j at the initial time; V p,j is the pore volume of the production well j.

[0514] In an embodiment, the comprehensive production optimization module is specifically used for:

[0515] When the production optimization model is a time-varying capacitance resistance model, the water injection amount of the injection well and the bottom hole pressure of the production well are taken as decision variables, and initial values of the decision variables are determined;

[0516] According to the initial values of the decision variables, the subsurface liquid production is calculated through the fractional flow model, and the surface water production and the water saturation are updated;

[0517] Based on the surface water production and the water saturation, the surface oil production and the subsurface liquid production are predicted through the fractional flow model;

[0518] According to a reasonable development technical strategy, the constraint conditions corresponding to the time-varying capacitance resistance model are determined, the production optimization objective function corresponding to the time-varying capacitance resistance model is iteratively solved, and the decision variables when the objective function reaches the optimal value are determined as the optimal decision variables;

[0519] According to the optimal decision variables, a development control scheme is obtained.

[0520] In an embodiment, the fractional flow model represents the relationship between the surface oil production and the subsurface liquid production, and the relationship between the surface water production and the subsurface liquid production through the saturation inverted by the historical fitting of the time-varying capacitance resistance model; the relationship between the surface oil production and the subsurface liquid production is as follows:

[0521] The relationship between the surface water production and the subsurface liquid production is as follows:

[0522] where F o,j is the ratio of the surface oil production to the subsurface liquid production of the production well j; q o,j is the surface oil production of the production well j; q t,j is the subsurface liquid production of the production well j; a i is a constant coefficient for predicting the oil production; S w,j is the water saturation in the control volume of the production well j; N is the polynomial degree; and F w,jis the ratio of surface water production and subterranean liquid production of production well j; q w,j is the surface water production of production well j; β i is a constant coefficient for predicting water production.

[0523] In an embodiment, the objective function corresponding to the time-varying capacitance resistance model is the net present value, and the net present value is as follows:

[0524] The constraint condition corresponding to the time-varying capacitance resistance model is as follows:

[0525] wherein w lim is the injection limit of a single well; w sum is the upper limit of total injection; q sum is the upper limit of total liquid production; NPV is the net present value; p wf,min is the lower limit of bottom hole pressure of a production well; p wf,max is the upper limit of bottom hole pressure of a production well; N t is the total number of optimization time steps; N pro is the number of production wells; N inj is the number of injection wells; r o (t k ) is the oil price at time t k ; b is the annual discount rate, q o,j (t k ) is the surface oil production of production well j at time t k ; r w (t k ) is the water treatment cost at time t k ; q w,j (t k ) is the surface water production of production well j at time t k ; r wi (t k ) is the injection cost at time t k ; w i (t k ) is the surface injection of injection well i at time t k ; t opt is the time when the optimization period starts.

[0526] In an embodiment, the first type of prediction model is a connected element simulation model, and the connected element simulation model is an improved connected element model.

[0527] The history fitting module is specifically configured to:

[0528] In the production optimization model is a connection unit simulation model, the model parameters are determined to be the conductivity, the pore volume and the productivity index, the nodes in the connection unit simulation model include well nodes and non-well nodes, the non-well nodes are arranged between well nodes with a distance exceeding a preset range or at positions requiring fine simulation in a region of interest;

[0529] The initial value of the model parameter is determined;

[0530] The node saturation is determined by using a front tracking algorithm in streamline simulation;

[0531] The water cut of each node is determined according to the node saturation;

[0532] The oil production of each node is calculated according to the water cut of each node;

[0533] A constraint optimization problem composed of the oil production of each node and the bottom hole flowing pressure is calculated as a constraint optimization problem of the connection unit simulation model;

[0534] The constraint optimization problem of the connection unit simulation model is iteratively solved based on the initial value of the model parameter, and the model parameter when the model error is less than a preset threshold is determined as the optimized model parameter.

[0535] In an embodiment, the history matching module is specifically configured to:

[0536] The pressure value at each node is calculated;

[0537] The flow rate on the connection unit is calculated according to the pressure value at each node, wherein two nodes constitute a connection unit;

[0538] The node saturation is determined by using a front tracking algorithm in streamline simulation.

[0539] In an embodiment, the history matching module is specifically configured to:

[0540] The saturation solving problem is determined;

[0541] The initial condition in the saturation solving problem is simplified to obtain an initial saturation;

[0542] The initial saturation is simplified as a multi-segment constant function to form a plurality of saturation solving sub-problems;

[0543] A new saturation solving sub-problem generated by possible collision of saturation fronts of adjacent two sub-problems in a moving process is determined:

[0544] After the front collision, a new saturation solving sub-problem is solved according to the collision point, a next possible collision is searched and a corresponding new front is inserted, until the collision exceeds the space and time boundary, a saturation profile at the end of the time step is determined, the saturation profile at the end of the time step is a piecewise constant function, and the saturation profile at the end of the time step is taken as an initial condition for front tracking in the next time step.

[0545] According to the saturation profile at the end of the time step, the node saturation is determined.

[0546] In an embodiment, the new saturation solving sub-problem is as follows:

[0547] Wherein, S w is the water saturation, x is a one-dimensional coordinate, f w is the water content, τ is the characteristic velocity, x collision is the position where the saturation collision occurs; t collision is the time when the saturation collision occurs, S w (x, t collision ) is the water saturation of x at t collision ;

[0548] In an embodiment, the constraint optimization problem of the connecting unit simulation model is as follows:

[0549] Wherein, N w is the number of nodes related to i, N tim is the number of history fitting time steps, is the calculated value of the bottom hole pressure of node i at t n , is the true value of the bottom hole pressure of node i at t n , is the calculated value of the formation oil production of node i at t n , is the true value of the formation oil production of node i at t n , is the initial value of the conductivity of the connecting unit composed of nodes i and j in the kth layer, is the initial value of the pore volume of the connecting unit composed of nodes i and j in the kth layer, S wc is the irreducible water saturation, is the initial value of the water saturation in the control range of node i in the kth layer, S or is the residual oil saturation, u is the decision variable of the history fitting of the connecting unit simulation model; V pt is the total pore volume; p wf,inj is the bottom hole pressure of the injection well node; p injis the average formation pressure of the injection well node; p wf,pro is the bottom hole pressure of the production well node; p pro is the average formation pressure of the production well node.

[0550] In an embodiment, the second type of prediction model is a non-proxy model in a pure data-driven proxy model or a reservoir numerical model, and the non-proxy model is a finite difference model and a streamline simulation model.

[0551] In an embodiment, implementing the comprehensive production optimization module is further used for:

[0552] Production optimization is performed based on the flow field diagnosis of the streamline simulation model by using the following steps:

[0553] The streamline information output by the streamline tracking in the streamline simulator, or the streamline information obtained by streamline tracking from the simulation results of the finite difference simulator, is used;

[0554] Based on the flow field diagnosis index, the production level of the injection-production well pair or single well is diagnosed, and the streamline injection-production optimization strategy is determined, the flow field diagnosis index includes injection efficiency, production efficiency, movable residual oil displacement efficiency, and movable residual oil production efficiency, and the streamline injection-production optimization strategy includes a flow rate updating method and an injection-production flow rate optimization criterion, and the injection-production flow rate optimization criterion is used to determine a flow rate updating weight.

[0555] Based on the streamline injection-production optimization strategy, streamline injection-production optimization is performed to obtain updated injection-production flow rates.

[0556] The updated injection-production flow rates are subjected to constraint processing to obtain an optimal injection-production scheme, and the optimal injection-production scheme is a development control scheme.

[0557] In an embodiment, the flow rate updating method includes: only updating the flow rate at the injection end and the flow rate at the production end of the injection-production well pair or single well according to the movable residual oil displacement efficiency; only updating the flow rate at the injection end and the flow rate at the production end of the injection-production well pair or single well according to the movable residual oil production efficiency; updating the flow rate at the injection end of the injection-production well pair or single well according to the movable residual oil displacement efficiency, and updating the flow rate at the production end of the injection-production well pair or single well according to the movable residual oil production efficiency; only updating the flow rate at the injection end and the flow rate at the production end of the injection-production well pair or single well according to the injection efficiency; only updating the flow rate at the injection end and the flow rate at the production end of the injection-production well pair or single well according to the production efficiency; updating the flow rate at the injection end of the injection-production well pair or single well according to the injection efficiency, and updating the flow rate at the production end of the injection-production well pair or single well according to the production efficiency.

[0558] In an embodiment, the flow rate optimization criterion includes an unbounded weight criterion and a bounded weight criterion.

[0559] The unbounded weight criterion does not consider the amplitude limit of the flow change, and only determines the flow update weight according to the difference between the flow field diagnostic index and the average value;

[0560] The bounded weight criterion introduces the diagnostic index difference and the linear scaling ratio of the diagnostic index difference to limit the weight change amplitude.

[0561] In an embodiment, the constraint processing includes:

[0562] For the injection-production well pair, the updated injection-production flow is converted to single well by superposition;

[0563] For the injection-production well pair, if there is an upper bound of the total injection of all injection wells and an upper bound of the total production of all production wells, the optimized single well injection-production flow is scaled;

[0564] For the injection-production well pair, if there is a constraint condition of the single well, the limit flow of the injection-production well that violates the pressure constraint is taken as the final optimized flow of the injection-production well that violates the pressure constraint, and the remaining injection-production is proportionally scaled to the flow of the injection-production well that does not violate the pressure constraint; or only the flow of the injection-production well that violates the pressure constraint is reduced, and the injection-production well that does not violate the pressure constraint is not processed.

[0565] For the single well, if there is a total flow constraint or a bottom hole pressure constraint, the optimized single well injection-production flow is scaled.

[0566] In an embodiment, the comprehensive production optimization module is further used for:

[0567] The production optimization is performed based on the numerical simulator surrogate model by using the following steps:

[0568] A production optimization problem is determined, the production optimization problem including decision variables, constraint conditions and an objective function;

[0569] Based on the production optimization problem, a compromise sampling strategy, a type of the numerical simulator surrogate model and a sample set are determined, the compromise sampling strategy taking the minimum response of the numerical simulator surrogate model and the maximum future search ability as dual objectives;

[0570] Based on the sample set, an initial numerical simulator surrogate model is constructed and taken as a current numerical simulator surrogate model, and the following steps are repeatedly executed until a termination condition is met, the optimal decision variable when the condition is met is output, and a development control scheme corresponding to the optimal decision variable is obtained:

[0571] According to the evaluation index of the future search ability and the evaluation index of the response of the numerical simulator surrogate model, a compromise coefficient, a comprehensive evaluation function is constructed;

[0572] According to the comprehensive evaluation function, a final potential candidate point is determined;

[0573] In the current numerical simulator proxy model, an oil reservoir production optimization objective function evaluation subprogram is used to evaluate the objective function of the final potential candidate point, wherein the oil reservoir production optimization objective function evaluation subprogram takes the decision variable as input;

[0574] The final potential candidate point is added to the sample set that has been evaluated by the function;

[0575] Based on the sample set that has been evaluated by the function, the numerical simulator proxy model is updated;

[0576] The updated numerical simulator proxy model is used as the current numerical simulator proxy model.

[0577] In an embodiment, the comprehensive production optimization module is further configured to:

[0578] An initial trial design method is used to select a preset number of samples in the solution space of the production optimization problem for objective function evaluation, and an initial sample set required for constructing the numerical simulator proxy model is generated.

[0579] In an embodiment, the comprehensive production optimization module is further configured to:

[0580] A constraint violation evaluation function is constructed according to the constraint condition;

[0581] Under the bounded constraint of the production optimization problem, the decision variable is randomly sampled to obtain an initial potential candidate point set, wherein the bounded constraint is one of the constraint conditions;

[0582] Based on the constraint violation evaluation function, the feasible potential candidate points are screened from the initial potential candidate point set to obtain a feasible potential candidate point set;

[0583] The feasible potential candidate point with the minimum comprehensive evaluation function is selected as the final potential candidate point, wherein the evaluation index of the future search capability in the comprehensive evaluation function uses the minimum distance from the sample to the sample that has been evaluated by the objective function.

[0584] In an embodiment, the comprehensive production optimization module is further configured to:

[0585] After obtaining the initial potential candidate point set, the randomly sampled points near the current best feasible point and the Latin hypercube sampled points within the bounded constraint range are added to the initial potential candidate point set.

[0586] In an embodiment, the comprehensive production optimization module is further configured to:

[0587] The constraint violation evaluation function value is calculated;

[0588] samples in the initial potential candidate point set whose constraint violation evaluation function value is greater than the constraint tolerance error are excluded;

[0589] the remaining set is referred to as a feasible potential candidate point set;

[0590] in the feasible potential candidate point set, samples whose distance to samples that have been evaluated by the objective function exceeds the allowable distance are removed.

[0591] In an embodiment, implementing the comprehensive production optimization module is further used for:

[0592] production optimization is performed based on the streamline simulation surrogate model using the following steps:

[0593] determining a production optimization problem, the production optimization problem including decision variables, constraint conditions and an objective function;

[0594] based on the production optimization problem, determining a compromise sampling strategy, a type of streamline simulation surrogate model and a sample set, the compromise sampling strategy taking the minimum response of the streamline simulation surrogate model and the maximum future search ability as dual objectives;

[0595] based on the sample set, constructing an initial streamline simulation surrogate model, and taking the initial streamline simulation surrogate model as a current streamline simulation surrogate model, repeating the following steps until a termination condition is met, outputting an optimal decision variable when the condition is met, and obtaining a development control scheme corresponding to the optimal decision variable:

[0596] constructing a comprehensive evaluation function according to an evaluation index of future search ability and an evaluation index of response of the streamline simulation surrogate model, and a compromise coefficient;

[0597] determining a final potential candidate point according to the comprehensive evaluation function;

[0598] in the current streamline simulation surrogate model, using a streamline optimization sub-function to evaluate the objective function for the final potential candidate point, wherein the streamline optimization sub-function takes as input variable parameters of a flow optimization criterion when production optimization is performed based on flow field diagnosis of the streamline simulation model, and calculates the objective function value based on simulation results of the development control scheme obtained when production optimization is performed based on the flow field diagnosis of the streamline simulation model;

[0599] adding the final potential candidate point to the sample set that has been evaluated by the function;

[0600] based on the sample set that has been evaluated by the function, updating the streamline simulation surrogate model;

[0601] taking the updated streamline simulation surrogate model as the current streamline simulation surrogate model.

[0602] In summary, the method and device provided by the embodiments of the present application are aimed at the many deficiencies of the existing oil reservoir production optimization methods (for example, the oil reservoir engineering method is strong in experience, poor in transplantability, and low in precision, the traditional numerical model method has high data requirements, high calculation cost, time-consuming and laborious, and it is difficult to obtain an optimal solution, the simplified physical model has strong multi-solution, weak long-term prediction ability, low reliability, and limited application range, the machine learning agent modeling optimization method still takes a long time in the agent construction stage, and the optimization efficiency is low for high-dimensional problems, and the streamline simulation optimization method needs to rely on the fitting results of the geological model, the flow field diagnosis method suitable for the whole life cycle of the oil reservoir has not been effectively defined, and the optimization method for the flow updating criterion has not been clearly proposed), starting from the definition of the optimization problem, the main principles and advantages and disadvantages of several production optimization methods (models) such as production optimization based on reservoir dynamic understanding, time-varying capacitance resistance model, connected unit simulation model, numerical simulator agent model, streamline simulation model, and streamline simulation agent optimization model are explained, and the implementation process is improved and explained.

[0603] Compared with the prior art, the method and device provided by the embodiments of the present application have at least the following beneficial effects:

[0604] The blank in the research on the coupling and interaction optimization of the reservoir dynamic understanding, the first type of prediction model and the second type of prediction model and the intelligent comprehensive production optimization method is filled, which is specifically embodied in the following aspects:

[0605] 1. For reservoir dynamic understanding, the regional level understanding, single well level understanding, and theoretical model understanding obtained by the main factor analysis, periodic dynamic analysis, node system analysis, and other production optimization methods are taken as the source of reservoir dynamic understanding update and are solidified as a reasonable development technical strategy to diagnose the reservoir state, and then an oil reservoir engineering optimization method for implementing oil and water well parameter control according to the dynamically updated reservoir understanding and development control scheme is proposed, which can be applied to high-frequency or irregular injection-production adjustment in oilfield field.

[0606] 2. For the capacitance resistance model, the changes of the volume coefficient, saturation, productivity index, time constant, and bottom hole pressure are considered, the saturation and flow control equations are derived on the basis of the oil-water two-phase material balance equation and the approximate method of the productivity equation, the time-varying capacitance resistance model (TV-CRM) is established, the connectivity evaluation and production optimization method based on the TV-CRM are proposed, and the calculation efficiency and precision of the CRM are improved by using the proposed agent optimization algorithm, and the efficiency and precision of the production optimization are improved.

[0607] 3. For the INSIM type simplified physical simulation, the multi-layer non-well node flow simulation model is extended and the saturation tracking algorithm is improved on the basis of the initial interwell numerical simulation model (INSIM), the connected unit simulation model (LUSM) is developed, and the efficiency and precision of the production optimization are improved.

[0608] 4、For numerical simulator proxy model and streamline simulation proxy model, radial basis function interpolation (or other machine learning regression model) is used to construct objective function proxy and define constraint violation evaluation function (or penalty function) to handle various linear or nonlinear equality / inequality constraints. Based on compromise sampling strategy, a proxy optimization algorithm for time-consuming constraint optimization problem is proposed, which can be applied to various production optimization scenarios to improve the efficiency and accuracy of production optimization.

[0609] 5、For streamline simulation optimization method, a flow field diagnosis method is proposed to evaluate injection-production control effect according to physical indicators such as injection-production efficiency and movable residual oil displacement efficiency from well pair and single well level. The constraint proxy optimization algorithm is combined with the injection-production flow optimization criterion in the flow field diagnosis optimization method of streamline simulation, and a proxy optimization method for streamline simulation model is proposed, which can further improve the optimization effect of the flow field diagnosis optimization method and overcome the deficiency of proxy optimization for solving high-dimensional optimization problems.

[0610] 6、Through comprehensive comparative analysis, the internal relationship between various production optimization methods is established by using connectivity evaluation (in the model calibration / history matching stage), and a comprehensive production optimization concept, workflow and implementation method are proposed under the framework of intelligent closed-loop reservoir management, which can adapt to different adjustment frequencies, consider real-time and long-term, balance efficiency and accuracy, and interact and cooperate with multiple optimization models.

[0611] 7、The comprehensive production optimization method integrates the advantages of reservoir dynamic understanding, the first type of prediction model (simplified physical model) and the second type of prediction model (pure data-driven proxy model and non-proxy model in reservoir numerical model), wherein the production optimization based on reservoir dynamic understanding (semi-empirical formula) is not dependent on physical or numerical prediction model, mainly uses the advantage of convenient operation to provide extremely short time or high frequency (such as daily) optimization suggestions; the production optimization based on the first type of prediction model mainly uses the speed efficiency advantage to provide short-term optimization suggestions, and the production optimization based on the second type of prediction model mainly uses the accuracy (reliability) advantage to provide long-term decision reference; the reservoir dynamic understanding (and development technical strategy, injection-production adjustment strategy, etc.) summarized from the main control factor analysis, periodic dynamic analysis, node system analysis, etc. is taken as a constraint to make the result of comprehensive production optimization more reasonable; in turn, the production optimization results of the first type of prediction model (simplified physical model) and the second type of prediction model (pure data-driven proxy model and reservoir numerical model) form new reservoir dynamic understanding together with the periodic dynamic analysis understanding; according to the comprehensive production optimization concept proposed in the application, various production optimization methods meeting different optimization requirements can be gradually developed and formed to form an intelligent reservoir closed-loop production optimization method and technology series with complementary advantages of each taking its strengths, real-time and long-term consideration, equal emphasis on physical meaning and data-driven, and significant improvement of optimization efficiency and accuracy.

[0612] The embodiment of the application further provides a computer device, and FIG. 52 is a schematic diagram of the computer device in the embodiment of the application. The computer device 5200 comprises a memory 5210, a processor 5220, and a computer program 5230 stored in the memory 5210 and capable of running on the processor 5220. The processor 5220 implements the above-mentioned reservoir comprehensive production optimization method when executing the computer program 5230.

[0613] The embodiment of the application further provides a computer readable storage medium, which stores a computer program. The computer program is executed by a processor to implement the above-mentioned reservoir comprehensive production optimization method.

[0614] The embodiment of the application further provides a computer program product, which comprises a computer program. The computer program is executed by a processor to implement the above-mentioned reservoir comprehensive production optimization method.

[0615] Those skilled in the art should understand that the embodiments of the application can be provided as a method, a system, or a computer program product. Therefore, the application can adopt a complete hardware embodiment, a complete software embodiment, or an embodiment combining software and hardware aspects. Moreover, the application can adopt the form of a computer program product implemented on one or more computer usable storage media (including but not limited to disk storage, CD-ROM, optical storage, etc.) containing computer usable program codes.

[0616] The computer program instructions can also be loaded onto a computer or other programmable data processing apparatus to cause a series of operational steps to be performed on the computer or other programmable apparatus to produce a computer implemented process such that the instructions which execute on the computer or other programmable apparatus provide steps for implementing the functions specified in the flowchart block or blocks or in the flowchart one or more blocks.

[0617] These computer program instructions can also be stored in a computer readable memory that can direct a computer or other programmable data processing apparatus to function in a particular manner, such that the instructions stored in the computer readable memory produce an article of manufacture including instructions which implement the function specified in the flowchart block or blocks or in the flowchart one or more blocks.

[0618] These computer program instructions can also be loaded onto a computer or other programmable data processing apparatus to cause a series of operational steps to be performed on the computer or other programmable apparatus to produce a computer implemented process such that the instructions which execute on the computer or other programmable apparatus provide steps for implementing the functions specified in the flowchart block or blocks or in the flowchart one or more blocks.

[0619] The specific embodiments described above are illustrative of specific embodiments of the application. It should be understood that the application can be practiced otherwise than as specifically written without departing from the spirit and scope of the present application. Accordingly, the present application is not intended to be limited to the specific form set forth herein, but on the contrary, is intended to cover such alternatives, modifications, and equivalents, as can be included within the spirit and scope of the present application.

Claims

1. A method of integrated reservoir production optimization, characterized in that, The method comprises the following steps: obtaining production history data of a target reservoir; updating reservoir dynamic knowledge according to regional level knowledge, single well level knowledge and production optimization model knowledge based on the production history data; performing cyclic history matching on model parameters of multiple production optimization models including a first type of prediction model and a second type of prediction model based on the production history data until optimal model parameters are obtained, wherein the model parameters of each production optimization model are history matched, and the interwell connectivity relationship of the multiple production optimization models is calibrated, so that the first type of prediction model provides an initial estimation of the interwell connectivity relationship for the second type of prediction model, and the model parameters of the first type of prediction model are calibrated by the second type of prediction model; determining the constraint conditions corresponding to each production optimization model according to the updated reservoir dynamic knowledge, solving the objective function corresponding to each production optimization model, obtaining the optimal decision variable corresponding to each production optimization model, and obtaining the development regulation scheme corresponding to each production optimization model; obtaining a comprehensive development regulation scheme of the reservoir according to the development regulation scheme corresponding to each production optimization model, wherein the comprehensive development regulation scheme is used to guide the development operation of the target reservoir; and updating the reservoir dynamic knowledge according to the regulation effect obtained by real-time monitoring of the development operation. The method for updating the reservoir dynamic knowledge according to the regional level knowledge, the single well level knowledge and the production optimization model knowledge based on the production history data comprises the following steps:

2. The method of claim 1, wherein, analyzing multiple main control factors affecting the reservoir and regional development effect and producing the observed production dynamic from multiple development aspects according to the production history data during the development of the reservoir; classifying the reservoirs according to different main control factors; obtaining the regional level knowledge for the reservoir classification conforming to the overall law of the region; obtaining the single well level knowledge from the periodic reservoir dynamic analysis and the node system analysis for the reservoir classification not conforming to the overall law of the region; obtaining the theoretical model knowledge according to the production optimization model; correcting the regional level knowledge by the single well level knowledge, and verifying the single well level knowledge by the regional level knowledge; updating the regional or single well level knowledge by the theoretical model knowledge, and constraining the theoretical model knowledge by the single well level knowledge; updating the reservoir dynamic knowledge according to the regional level knowledge, the single well level knowledge and the production optimization model knowledge; and determining a reasonable development technical strategy according to the updated reservoir dynamic knowledge. The method for determining the constraint conditions corresponding to each production optimization model according to the updated reservoir dynamic knowledge comprises the following steps: determining the constraint conditions corresponding to each production optimization model according to the reasonable development technical strategy. The method further comprises the following steps:

3. The method of claim 1, wherein, performing reservoir state diagnosis according to the reasonable development technical strategy to determine problem wells; generating a development regulation scheme for the problem wells, wherein the development regulation scheme is used to adjust the operation of the problem wells; and updating the reservoir dynamic knowledge and the development regulation scheme according to the real-time monitoring and effect feedback of the adjustment operation. The first type of prediction model is a capacitance resistance model, a connection unit model or a flow network model. ​ 4. The method of claim 1, wherein, ​ 5. The method of claim 4, wherein, The capacitance resistance model is a time-varying capacitance resistance model, the time-varying capacitance resistance model is a modified capacitance resistance model, and the method further comprises: On the basis of approximations of oil-water two-phase material balance equations and deliverability equations, a flow rate control equation and a saturation control equation are determined; and A time-varying capacitance resistance model is established according to the flow rate control equation and the saturation control equation.

6. The method of claim 5, wherein, The flow control equation of the time-varying capacitance resistance model is as follows: wherein q t (t k+1 ) and q t (t k ) are the subterranean liquid production rates at time t k+1 and time t k , respectively, J(t k+1 ) and J(t k ) are the productivity indices at time t k+1 and time t k , respectively, τ(t k+1 ) and τ(t k ) are the time constants at time t k+1 and time t k , respectively, w(t k+1 ) is the surface injection rate at time t k+1 , B w is the water phase volume factor, p wf (t k+1 ) and p wf (t k ) are the bottom hole pressures at time t k+1 and time t k , respectively. The saturation control equation of the time-varying capacitance resistance model is as follows: wherein S w (t k+1 ) and S w (t k ) are the water saturation of the production well control volume at time t k+1 and time t k , respectively, V p is the pore volume of the production well control volume, w is the surface injection rate, q w is the surface water production rate, C w is the water phase compressibility, C φ is the rock compressibility, C t is the overall compressibility, q t is the subsurface liquid production rate, S o (t k+1 ) and S o (t k ) are the oil saturation of the production well control volume at time t k+1 and time t k , respectively, q o is the surface oil production rate, B o is the oil volume factor, C o is the oil compressibility.

7. The method of claim 6, wherein, Based on production history data, model parameters of a plurality of production optimization models are history-fitted, including: When the production optimization model is a time-varying capacitance resistance model, the model parameters are determined to be a connectivity factor, a control volume pore volume, a deliverability index, and a saturation; Initial values of the model parameters are determined; An evolution process of the model parameters is determined; and Based on the initial values of the model parameters and the evolution process, a constraint optimization problem constituted by model errors is iteratively solved, and the model parameters when the model errors are less than a preset threshold are determined as optimized model parameters.

8. The method of claim 7, wherein, The model error is as follows: where N pro is the number of production wells; N inj is the number of injection wells; error is the model estimation error of the subsurface liquid production rate of a production well at any time; N tim is the number of production history data points; E sum is the sum of the estimation error of the subsurface liquid production rate of all production wells over all time series; and Qj(k+1) = Qj(k) + (Qj(k) - Qj(k-1)) * (1 - e^(-Dj(k) * Δt)) (1) A true value of a subsurface liquid production rate of the jth production well at the k+1 time point; The constrained optimization problem is as follows: where S wc is the irreducible water saturation; S or is the residual oil saturation; V pt is the total pore volume of the production system, f ij is the connectivity coefficient of the injection well i to the production well j, Sj(t) = Vj(t) / Vj p,j Vj = Vj(t) = Vj(0) + ∑∑∑∑∑∑∑∑∑∑∑∑∑∑∑∑∑∑∑∑∑∑∑ 9. The method of claim 6, wherein, According to the updated reservoir performance understanding, a constraint condition corresponding to each production optimization model is determined, a target function corresponding to each production optimization model is solved, optimal decision variables corresponding to each production optimization model are obtained, and a development control scheme corresponding to each production optimization model is obtained, including: When the production optimization model is a time-varying capacitance resistance model, a water injection rate of an injection well and a bottom hole pressure of a production well are taken as decision variables, and initial values of the decision variables are determined; According to the initial values of the decision variables, a subsurface liquid production rate is calculated through a flow rate distribution model, and a surface water production rate and a water saturation are updated; Based on the surface water production rate and the water saturation, a surface oil production rate and a subsurface liquid production rate are predicted through the flow rate distribution model; According to a reasonable development technical strategy, a constraint condition corresponding to the time-varying capacitance resistance model is determined, a production optimization target function corresponding to the time-varying capacitance resistance model is iteratively solved, and the decision variables when the target function reaches an optimal value are determined as optimal decision variables; and According to the optimal decision variables, the development control scheme is obtained.

10. The method of claim 9, wherein, The split flow model expresses the relationship between the surface oil production and the subsurface liquid production, the relationship between the surface water production and the subsurface liquid production by the time-varying capacitance resistance model history matching inversion saturation; the relationship between the surface oil production and the subsurface liquid production is as follows: The relationship between the surface water production and the subterranean liquid production is as follows: where F o,j is the ratio of surface oil production to subsurface fluid production for producer j; q o,j is the surface oil production for producer j; q t,j is the subsurface fluid production for producer j; a i is a constant factor for predicting oil production; S w,j is the water saturation in the control volume for producer j; N is the polynomial order; F w,j is the ratio of surface water production to subsurface fluid production for producer j; q w,j is the surface water production for producer j; β i is a constant factor for predicting water production.

11. The method of claim 9, wherein, The objective function corresponding to the time-varying capacitance resistance model is the net present value, which is as follows: The constraint conditions corresponding to the time-varying capacitance resistance model are as follows: where w lim is the injection limit for a single well; w sum is the upper bound for total injection; q sum is the upper bound for total liquid production; NPV is the net present value; p wf,min is the lower bound for the bottom hole pressure of a production well; p wf,max is the upper bound for the bottom hole pressure of a production well; N t is the total number of optimization time steps; N pro is the number of production wells; N inj is the number of injection wells; r o (t k ) is the oil price at time t k ; b is the annual discount rate; q o,j (t k ) is the surface oil production of production well j at time t k ; r w (t k ) is the water production cost at time t k ; q w,j (t k ) is the surface water production of production well j at time t k ; r wi (t k ) is the injection cost at time t k ; w i (t k ) is the surface injection of injection well i at time t k ; t opt is the initial time of the optimization period.

12. The method of claim 4, wherein, The first type of prediction model includes a connection unit simulation model, and the connection unit simulation model is a modified connection unit model; Based on production history data, model parameters of a plurality of production optimization models are history-fitted, including: When the production optimization model is a connection unit simulation model, the model parameters are determined to be a conductivity, a pore volume, and a deliverability index, nodes in the connection unit simulation model include well nodes and non-well nodes, and the non-well nodes are arranged between well nodes with a distance exceeding a preset range or at positions requiring fine simulation in a region of interest; Initial values of the model parameters are determined; A front tracking algorithm in streamline simulation is used to determine node saturations; According to the node saturations, water cut rates of the nodes are determined; According to the water cut rates of the nodes, oil production rates of the nodes are calculated; A constraint optimization problem constituted by the oil production rates and bottom hole flow pressures of the nodes is calculated as a constraint optimization problem of the connection unit simulation model; and Based on the initial values of the model parameters, the constraint optimization problem of the connection unit simulation model is iteratively solved, and the model parameters when the model errors are less than a preset threshold are determined as optimized model parameters.

13. The method of claim 12, wherein, Determining the node saturation degree by using the front tracking algorithm in streamline simulation, comprising: calculating the pressure value at each node; calculating the flow rate on the connecting unit according to the pressure value at each node, wherein two nodes constitute a connecting unit; and Determining the node saturation degree by using the front tracking algorithm in streamline simulation.

14. The method of claim 13, wherein, Determining the node saturation degree by using the front tracking algorithm in streamline simulation, comprising: determining the saturation degree solving problem; simplifying the initial condition in the saturation degree solving problem to obtain the initial saturation degree; simplifying the initial saturation degree into a multi-segment constant function to form a plurality of saturation degree solving sub-problems; determining the new saturation degree solving sub-problems possibly generated by the collision of the saturation degree fronts of two adjacent sub-problems during movement: after the collision of the fronts, solving the new saturation degree solving sub-problems according to the collision point, finding the next possible collision and inserting the corresponding new front, until the collision exceeds the space and time boundary, determining the saturation profile at the end of the time step, the saturation profile at the end of the time step being a piecewise constant function, taking the saturation profile at the end of the time step as the initial condition for the front tracking of the next time step; and determining the node saturation degree according to the saturation profile at the end of the time step.

15. The method of claim 14, wherein, The new saturation solving sub-problem is as follows: where S w is the water saturation, x is the one-dimensional coordinate, f w is the water fraction, τ is the characteristic velocity, x collision is the position where the saturation collision occurs; t collision is the time when the saturation collision occurs; S w (x, t collision ) is the water saturation at x at t collision ; S w,left is the initial saturation corresponding to the saturation initial value on the left side of the collision position; S w,right is the saturation initial value on the right side of the collision position; t is the production time.

16. The method of claim 12, wherein, The constraint optimization problem of the connection unit simulation model is as follows: where N w is the number of nodes related to i, N tim is the number of historical fitting time steps, the calculated value of the bottom-hole pressure for node i at time t, n the calculated value of the bottom-hole pressure for node i at time t, the true value of the bottom-hole pressure for node i at time t n ​ For node i in t n Calculated formation oil production at any given time. the true value of the formation oil production for node i at time t n ​ the initial value of the conductance of the connection unit consisting of nodes i and j in the kth layer, S is the initial value of the pore volume of the connection unit consisting of nodes i and j in the kth layer wc S is the initial value of the pore volume of the connection unit consisting of nodes i and j in the kth layer S is the initial value of water saturation for the kth layer node i control range or S is the residual oil saturation, u is the decision variable; V pt S is the total pore volume; p wf,inj S is the bottom hole pressure of the injection well node; p inj S is the average formation pressure of the injection well node; p wf,pro S is the bottom hole pressure of the production well node; p pro S is the average formation pressure of the production well node.

17. The method of claim 1, wherein, The second type of prediction model is a pure data-driven proxy model or a non-proxy model in a reservoir numerical model, and the non-proxy model is a finite difference model and a streamline simulation model.

18. The method of claim 17, wherein, Further comprising: performing production optimization based on the flow field diagnosis of the streamline simulation model by using the following steps: obtaining the streamline information output by the streamline tracking in the streamline simulator, or obtaining the streamline information by performing streamline tracking on the simulation results of the finite difference simulator; diagnosing the production level of the injection-production well pair or single well based on the flow field diagnosis index, determining the streamline injection-production optimization strategy, the flow field diagnosis index including injection efficiency, production efficiency, movable residual oil displacement efficiency, and movable residual oil production efficiency, the streamline injection-production optimization strategy including a flow rate updating mode and an injection-production flow rate optimization criterion, the injection-production flow rate optimization criterion being used to determine the flow rate updating weight; performing streamline injection-production optimization based on the streamline injection-production optimization strategy to obtain updated injection-production flow rates; and adding constraint processing to the updated injection-production flow rates to obtain an optimal injection-production scheme, the optimal injection-production scheme being a development control scheme.

19. The method of claim 18, wherein, The flow rate updating mode includes: updating only the flow rates at the injection end and the production end of the injection-production well pair or single well according to the movable residual oil displacement efficiency; updating only the flow rates at the injection end and the production end of the injection-production well pair or single well according to the movable residual oil production efficiency; updating the flow rate at the injection end of the injection-production well pair or single well according to the movable residual oil displacement efficiency, and updating the flow rate at the production end of the injection-production well pair or single well according to the movable residual oil production efficiency; updating only the flow rates at the injection end and the production end of the injection-production well pair or single well according to the injection efficiency; updating only the flow rates at the injection end and the production end of the injection-production well pair or single well according to the production efficiency; updating the flow rate at the injection end of the injection-production well pair or single well according to the injection efficiency, and updating the flow rate at the production end of the injection-production well pair or single well according to the production efficiency.

20. The method of claim 18, wherein, The flow rate optimization criterion includes an unbounded weight criterion and a bounded weight criterion. The flow rate optimization criterion includes an unbounded weight criterion and a bounded weight criterion. The unbounded weight criterion does not consider the amplitude limit of the flow variation, and only determines the flow update weight according to the difference between the flow field diagnostic index and the average value; And The bounded weight criterion introduces the diagnostic index difference and the linear scaling ratio of the diagnostic index difference to limit the weight variation amplitude.

21. The method of claim 18, wherein, The constraint processing includes: For the injection-production well pair, the updated injection-production flow is converted to a single well by superposition; For the injection-production well pair, if there is an upper bound of the total injection of all injection wells and an upper bound of the total production of all production wells, the optimized single well injection-production flow is scaled; For the injection-production well pair, if there is a constraint condition of a single well, the limit flow of the injection-production well that violates the pressure constraint is taken as the final optimized flow of the injection-production well that violates the pressure constraint, and the remaining injection-production is proportionally scaled to the flow of the injection-production well that does not violate the pressure constraint again; or only reduce the flow of the injection-production well that violates the pressure constraint without processing the injection-production well that does not violate the pressure constraint; and For a single well, if there is a total flow constraint or a bottom hole pressure constraint, the optimized single well injection-production flow is scaled.

22. The method of claim 1, wherein, Further comprising: performing production optimization based on a numerical simulator surrogate model by using the following steps: determining a production optimization problem, the production optimization problem including decision variables, constraint conditions and an objective function; based on the production optimization problem, determining a compromise sampling strategy, a type of numerical simulator surrogate model and a sample set, the compromise sampling strategy taking the minimum numerical simulator surrogate model response and the maximum future search ability as dual objectives; based on the sample set, constructing an initial numerical simulator surrogate model and taking it as a current numerical simulator surrogate model, repeatedly executing the following steps until a termination condition is met, outputting an optimal decision variable when the condition is met, and obtaining a development control scheme corresponding to the optimal decision variable: constructing a comprehensive evaluation function according to an evaluation index of future search ability and an evaluation index of numerical simulator surrogate model response, and a compromise coefficient; determining a final potential candidate point according to the comprehensive evaluation function; in the current numerical simulator surrogate model, using an oil reservoir production optimization objective function evaluation subroutine to evaluate the objective function of the final potential candidate point, wherein the oil reservoir production optimization objective function evaluation subroutine takes the decision variable as input; adding the final potential candidate point to the sample set that has been evaluated by the function; based on the sample set that has been evaluated by the function, updating the numerical simulator surrogate model; and taking the updated numerical simulator surrogate model as the current numerical simulator surrogate model for determining the potential candidate point in the next iteration step.

23. The method of claim 22, wherein, Based on the production optimization problem, the sampling strategy, the type of numerical simulator surrogate model and the sample set are determined, including: using an initial test design method to select a predetermined number of samples in the solution space of the production optimization problem for objective function evaluation, and generating a sample set required for constructing an initial numerical simulator surrogate model.

24. The method of claim 22, wherein, According to the comprehensive evaluation function, the final potential candidate point is determined, including: constructing a constraint violation evaluation function according to the constraint condition; In the bounded constraint of the production optimization problem, random sampling is performed on the decision variables to obtain an initial potential candidate point set, wherein the bounded constraint is one of the constraint conditions; Based on the constraint violation evaluation function, feasible potential candidate points are screened from the initial potential candidate point set to obtain a feasible potential candidate point set; and The feasible potential candidate point with the minimum comprehensive evaluation function is selected as the final potential candidate point, wherein the evaluation index of the future search capability in the comprehensive evaluation function adopts the minimum distance from the sample to the sample that has been subjected to the objective function evaluation.

25. The method of claim 24, wherein, After obtaining the initial potential candidate point set, the following steps are further included: Random sampling points near the current best feasible point and Latin hypercube sampling points within the bounded constraint range are added to the initial potential candidate point set.

26. The method of claim 25, wherein, Based on the constraint violation evaluation function, feasible potential candidate points are screened from the initial potential candidate point set to obtain a feasible potential candidate point set, including: The constraint violation evaluation function value is calculated; Samples in the initial potential candidate point set with a constraint violation evaluation function value greater than the constraint tolerance error are excluded; The remaining set is referred to as the feasible potential candidate point set; and In the feasible potential candidate point set, samples with a distance to the sample that has been subjected to the objective function evaluation exceeding the allowable distance are removed.

27. The method of claim 1, wherein, Further included are the following steps for production optimization based on the streamline simulation proxy model: A production optimization problem is determined, the production optimization problem including decision variables, constraint conditions and an objective function; Based on the production optimization problem, a compromise sampling strategy, types of the streamline simulation proxy model and a sample set are determined, the compromise sampling strategy taking the minimum response of the streamline simulation proxy model and the maximum future search capability as dual objectives; Based on the sample set, an initial streamline simulation proxy model is constructed and used as a current streamline simulation proxy model, and the following steps are repeatedly executed until a termination condition is met, the optimal decision variable when the condition is met is output, and a development control scheme corresponding to the optimal decision variable is obtained: According to the evaluation index of the future search capability and the evaluation index of the response of the streamline simulation proxy model, a compromise coefficient is used to construct a comprehensive evaluation function; The final potential candidate point is determined according to the comprehensive evaluation function; In the current streamline simulation proxy model, the objective function evaluation of the final potential candidate point is performed using a streamline optimization sub-function, wherein the streamline optimization sub-function takes the variable parameter of the flow optimization criterion when the flow field diagnosis based on the streamline simulation model is used for production optimization as input, and calculates the objective function value using the simulation result of the development control scheme obtained when the flow field diagnosis based on the streamline simulation model is used for production optimization; The final potential candidate point is added to the sample set that has been subjected to the function evaluation; and The streamline simulation proxy model is updated based on the current sample set that has been subjected to the function evaluation; and The updated streamline simulation proxy model is used as the current streamline simulation proxy model.

28. An integrated reservoir production optimization apparatus, comprising: Included are: A production history data obtaining module for obtaining production history data of a target reservoir; A reservoir dynamic understanding updating module for updating reservoir dynamic understanding according to the production history data, performing regional horizontal understanding, single well horizontal understanding and production optimization model understanding; The history matching module is configured to perform cyclic history matching on model parameters of a plurality of production optimization models, including a first type of prediction model and a second type of prediction model, based on production history data, until optimal model parameters are obtained, by: performing history matching on the model parameters of each production optimization model, and calibrating interwell connectivity for the plurality of production optimization models, such that the first type of prediction model provides an initial estimation of interwell connectivity for the second type of prediction model, and the second type of prediction model calibrates the model parameters of the first type of prediction model; The integrated production optimization module is configured to determine constraint conditions corresponding to each production optimization model according to the updated reservoir dynamic understanding, solve an objective function corresponding to each production optimization model, obtain optimal decision variables corresponding to each production optimization model, and obtain a development control scheme corresponding to each production optimization model; and obtain a comprehensive development control scheme of the reservoir according to the development control schemes corresponding to each production optimization model, the comprehensive development control scheme of the reservoir being used to guide development operations of the target reservoir. The adjustment effect feedback and update module is configured to update the reservoir dynamic understanding according to a control effect obtained by real-time monitoring of the development operations. The processor implements the method of any one of claims 1 to 27 when executing the computer program.

29. A computer device comprising a memory, a processor, and a computer program stored on the memory and executable on the processor, wherein the computer program comprises the steps of: The computer readable storage medium stores the computer program, and the computer program is executed by the processor to implement the method of any one of claims 1 to 27.

30. A computer-readable storage medium, characterized in that, The computer program product comprises the computer program, and the computer program is executed by the processor to implement the method of any one of claims 1 to 27.

31. A computer program product, characterised in that, ​

Images