Full-grid scale automatic history matching method for water drive reservoirs based on perturbation pulse inversion

Abstract

The application discloses a water drive reservoir full-grid scale automatic history matching method based on disturbance pulse inversion and belongs to the technical field of reservoir numerical simulation. The method is characterized in that: a time-flow line set is obtained by discretizing a seepage field by using streamline simulation; a saturation propagation control equation along a flow line is obtained by using a characteristic line method to solve a water phase seepage flow equation; a water saturation disturbance transmission mechanism is constructed; a water cut error mapping with a production well is established and a linear equation set is formed; a saturation field optimal correction amount is obtained by solving a least square problem with a water cut fitting error as a constraint; a saturation field is updated and new flow lines are generated according to the saturation field; and iteration is performed until the fitting reaches a standard; and a permeability field is corrected based on a time sequence correlation coefficient of injection and production wells and a production well pressure prediction error. The application has the characteristics of being insensitive to an initial saturation, high inversion efficiency and high fitting precision, and is suitable for large-scale complex water drive reservoir history matching and prediction.

Description

Technical Field

[0001] This invention belongs to the field of reservoir numerical simulation technology, specifically relating to an automatic history fitting method for water-drive reservoirs at the full grid scale based on perturbation pulse inversion. Background Technology

[0002] Reservoir history fitting is a crucial step in reservoir numerical simulation. Its goal is to achieve an optimal match between the numerical simulation results and actual oilfield production dynamics by adjusting geological model parameters. An accurate model is fundamental for production prediction, development method optimization, and enhanced oil recovery. Traditional reservoir history fitting methods mainly include manual parameter tuning, stochastic optimization algorithms, gradient-based optimization methods, streamline-assisted methods, and data-driven methods. Manual parameter tuning methods rely heavily on engineer experience, are highly subjective, and inefficient. Stochastic optimization algorithms, such as genetic algorithms and particle swarm optimization, are computationally expensive, prone to getting trapped in local optima, and struggle to guarantee solution stability. Gradient-based optimization methods, such as the adjoint state method, have faster convergence speeds, but are highly sensitive to the selection of the initial model, and these methods often focus on static parameter inversion, neglecting the direct impact of dynamic changes in the water saturation field on production data. Streamline-assisted methods reduce computational complexity, but typically fail to achieve dynamic inversion of the saturation field, thus failing to fully utilize the saturation information in the water cut data. Data-driven methods, such as machine learning-based alternative models, heavily rely on large-scale, high-quality training datasets, lack physical interpretability, and have insufficient generalization ability when data is scarce. In summary, when facing large-scale, highly heterogeneous water-drive reservoirs, existing technologies struggle to balance computational efficiency and inversion accuracy, are highly dependent on the initial model, and generally fail to effectively utilize dynamic data such as water cut to achieve coordinated and efficient inversion of the permeability and saturation fields. Therefore, there is an urgent need for a robust automatic history fitting method that does not rely on an accurate initial model, can simultaneously invert key parameters, and is computationally efficient, in order to improve the predictive reliability of reservoir models and provide solid support for development decisions. Summary of the Invention

[0003] The purpose of this invention is to overcome the shortcomings of existing automatic reservoir history fitting methods in terms of efficiency, accuracy, dependence on initial models, and comprehensiveness of parameter inversion, and to provide an automatic history fitting method for water-driven reservoirs at the full grid scale based on perturbation pulse inversion. This method does not rely on an accurate initial water saturation model and can simultaneously, efficiently, and stably invert the permeability field and water saturation field that highly match the production history data.

[0004] To achieve the above objectives, the technical solution of the present invention is as follows:

[0005] An automatic history fitting method for water-drive reservoirs at the full grid scale based on perturbation pulse inversion includes the following steps:

[0006] Step S1: Based on the streamline simulation framework, the seepage field is discretized into multiple continuous time intervals and corresponding streamline sets;

[0007] Step S2: Solve the aqueous phase seepage equation using the method of characteristics to obtain the governing equation for the propagation of saturation along streamlines;

[0008] Step S3: Based on the governing equations, construct the perturbation transmission mechanism of water saturation micro-perturbation, establish its linear mapping relationship with the water cut error of the production well, and form a system of linear equations;

[0009] Step S4: Using the historical fitting error of water content as a constraint, construct and solve the least squares optimization problem to obtain the optimal correction amount of the water saturation field;

[0010] Step S5: Update the water saturation field using the optimal correction amount and generate a new streamline distribution;

[0011] Step S6: Repeat steps S2 to S5 for iteration until the moisture content fitting meets the accuracy requirements and the saturation field history fitting is achieved.

[0012] Step S7: Correct the permeability field based on the time series correlation coefficient between injection and production wells and the pressure prediction error of production wells.

[0013] Furthermore, the specific process of step S1 is as follows:

[0014] Define a continuous time discrete point sequence The entire simulation time is divided into a continuous time interval ;

[0015] In each time interval Internally, using streamline simulation, the transient flow field is represented as a set of discrete streamlines starting from the injection well and ending at the production well. ;

[0016] Along each streamline Using arc length coordinates The space is parameterized, where This is the total length of the streamline;

[0017] Using the streamlines as a framework, the reservoir seepage field is discretized into a series of dynamic flow tube units, and the fluid flow behavior within each flow tube unit is dominated by the streamlines it is attached to.

[0018] In subsequent iterations, the streamline distribution and the flow tube unit division are dynamically updated based on the updated reservoir parameter field.

[0019] Furthermore, the specific process of step S2 is as follows:

[0020] Taking oil-water two-phase flow as an example, traditional numerical simulation techniques simulate the evolution of the seepage field by solving the mass conservation equation that controls the underground flow. The conservation equation is expressed as follows:

[0021] ;

[0022] ;

[0023] in, Porosity For crude oil density, Oil saturation For absolute penetration rate, The relative permeability of the oil phase. For crude oil viscosity, For the flow rate of oil, The density of water, Water saturation The relative permeability of the aqueous phase. The viscosity of water, For water flow rate, For time, Let the coordinates be the arc length along the streamline. For pressure;

[0024] Assuming incompressible flow, and defining the Darcy velocity of water as... Simplified to:

[0025] ;

[0026] in, Darcy's velocity for water;

[0027] Solve using the method of characteristics, let The saturation propagation equation along the streamline is obtained as follows:

[0028] ;

[0029] Define the flight time along the streamline:

[0030] ;

[0031] in, Let these be the initial arc length coordinates along the streamline. Let the coordinates be the arc length along the streamline. For position arrive Flight time, For the total Darcy speed, This is the flow rate function of the water.

[0032] Integrating the saturation equation along the streamlines, we obtain the result at position and time The expression for water saturation:

[0033] ;

[0034] in, For position and time water saturation For position And the initial water saturation at time 0, From position Arrive at the location The time.

[0035] Furthermore, the specific process of step S3 is as follows:

[0036] Based on the aforementioned saturation propagation equation, when for mesh cells Apply perturbation to the initial water saturation hour:

[0037] ;

[0038] in, To disturb the water saturation;

[0039] By utilizing simplified discretization It can identify grid cells that affect the water cut of a production well within a specific time interval; thus, the impact of saturation disturbance in any grid cell on the dynamic response of the production well can be quantitatively tracked.

[0040] Consider the location starting from the injection well and ending at the production well. streamline On the grid cells ;set up , making ; applied to the unit Disturbance on initial water saturation It will propagate along the streamline and in time Then arrived at the production well :

[0041] ;

[0042] in, Indicates from position Arrive at the location Flight time;

[0043] Through calculation It can be determined that the saturation disturbance is in the production well. The arrival time; by using the error between the predicted value and the observed water content, the grid cell can be derived. Optimal correction for water saturation Flight time is:

[0044] ;

[0045] In time interval Inside, for those passing through grid cells And eventually reach the production well streamline , applied to disturbances on Will have an impact on time of The moisture content at that location; the resulting change in moisture content is ;set up This indicates a line whose endpoint is a production well. The streamlines; then, for any time... and production wells The following relationships hold true:

[0046] ;

[0047] in, It is the most recent streamline start time. To flow through the grid cells Location To the end of the streamline Flight time, For grid cells water saturation For grid cells Disturbance water saturation, For production wells In time The moisture content fitting error.

[0048] Furthermore, the specific process of step S4 is as follows:

[0049] Based on the linear relationship established in step S3, the equations for all times and all production wells are combined to construct a linear system. ,in:

[0050] This is the water saturation correction vector for the grid cells to be solved;

[0051] The coefficient matrix is ​​expressed as follows:

[0052] ;

[0053] in, For time point indexing, For the number of production wells, For production well index, row index Corresponding to production time and production wells Column index Corresponding grid cell , It is the most recent streamline start time. To flow through the grid cells Location To the end of the streamline Flight time, For grid cells water saturation For a line whose endpoint is a production well Streamlines;

[0054] Let the moisture content prediction error vector be expressed as follows:

[0055] ;

[0056] in, For the first production well In the Time points The error in moisture content prediction;

[0057] Establish and solve the following constrained least squares optimization problem to obtain the optimal correction amount. :

[0058] ;

[0059] in, This is the physical lower limit of water saturation. This represents the physical upper limit of water saturation. For grid cells The current water saturation value, The mesh element to be solved Water saturation correction vector This is a set of grid cells for which changes in water saturation are not permitted.

[0060] Furthermore, the specific process of step S5 is as follows:

[0061] The optimal correction vector Each element in With the corresponding grid cell The current water saturation values ​​are added together to update the water saturation field in the reservoir model. Then, based on the updated water saturation field and the current reservoir pressure field, streamline simulation is performed again to generate a new set of discrete streamlines that reflects the current seepage field state, providing a basis for the next iteration.

[0062] Furthermore, the specific process of step S6 is as follows:

[0063] Using the updated water saturation field and the new streamline distribution as input, steps S2, S3, S4, and S5 are repeated to form a closed iterative optimization loop. After each iteration, the historical fitting error of water cut for all production wells is calculated, and it is determined whether the preset accuracy convergence criterion has been met. If the accuracy requirement is met, the iteration terminates, and the final inverted water saturation field is output. If not, the next iteration continues until convergence.

[0064] Furthermore, the specific process of step S7 is as follows:

[0065] The Pearson correlation coefficient is used to represent the influence between two time series, and its expression is as follows:

[0066] ;

[0067] in, For two frequency domains and The correlation coefficient between them Representing the frequency domain and frequency domain covariance, For frequency domain The average value, For frequency domain The average value, For frequency domain standard deviation For frequency domain standard deviation Represents the expectation operator;

[0068] Calculation of water injection wells With production well Correlation coefficient between injection volume and product volume time series :

[0069] ;

[0070] in, For water injection wells With production wells The correlation coefficient, The lag coefficient, For injection well The injection volume time series, from time index 1 to Lagging data, For production wells The time series of liquid production, from the time index arrive ;

[0071] Combined with production wells Bottom hole pressure prediction error For connecting water injection wells With production wells streamline Grid cells passed through The penetration rate is corrected:

[0072] ;

[0073] in, This is the corrected permeability. The permeability before correction. and To adjust the coefficient, For production wells The error in bottom hole pressure prediction To inject from the well Departure and in time Arrival at the production well The streamline set, For grid cells Located on one of the streamlines.

[0074] The present invention has the following beneficial effects:

[0075] First, this invention exhibits strong robustness to the initial water saturation model. Even given a physically unreasonable initial field, it can automatically recover a reasonable distribution consistent with the production history through an iterative framework, significantly reducing the dependence on the accuracy of the initial geological model. Second, this invention achieves direct dynamic inversion of the water saturation field within a streamline simulation framework and performs co-optimization with the permeability field. By constructing a linear mapping relationship between saturation perturbation and production response, the history fitting problem is transformed into a constrained least squares problem for solution. This significantly improves computational efficiency while ensuring a clear physical mechanism, achieving high-precision and rapid fitting of dynamic data such as water cut. Furthermore, this invention can effectively handle abrupt changes in water cut data caused by water shut-off, profile control, and other measures, and has good adaptability and stability to complex noise in actual production data. In summary, this invention provides an efficient, stable, and practical automatic history fitting solution for large-scale water-drive reservoirs. Attached Figure Description

[0076] Figure 1 This is a flowchart illustrating the automatic history fitting method for water-drive reservoirs at the full grid scale based on perturbation pulse inversion according to the present invention.

[0077] Figure 2 This is a diagram of the initial water saturation model in an embodiment of the present invention;

[0078] Figure 3 This is a diagram of the water saturation model after historical fitting in an embodiment of the present invention;

[0079] Figure 4 This is a diagram of the initial permeability model in an embodiment of the present invention;

[0080] Figure 5 This is a permeability model diagram after historical fitting in an embodiment of the present invention;

[0081] Figure 6 This is an iterative fitting curve of reservoir water cut in an embodiment of the present invention;

[0082] Figure 7 This is a historical fitting effect diagram of reservoir fluid production and oil production in an embodiment of the present invention;

[0083] Figure 8 This is a historical fitting effect diagram of water cut, fluid production and oil production of production well 8666 under complex noise conditions in this embodiment of the invention;

[0084] Figure 9 This is a historical fitting effect diagram of water cut, fluid production and oil production of production well 8537 under complex noise conditions in this embodiment of the invention;

[0085] Figure 10This is a historical fitting effect diagram of water cut, fluid production and oil production of production well 8505 under complex noise conditions in an embodiment of the present invention;

[0086] Figure 11 This is a historical fitting effect diagram of water cut, fluid production and oil production of production well 8331-1 under complex noise conditions in an embodiment of the present invention;

[0087] Figure 12 This is a historical fitting effect diagram of water cut, fluid production and oil production of production well 8061 under complex noise conditions in an embodiment of the present invention;

[0088] Figure 13 This is a historical fitting effect diagram of water cut, fluid production and oil production of production well 8039 under complex noise conditions in an embodiment of the present invention. Detailed Implementation

[0089] The present invention will be further described in detail below with reference to the embodiments and accompanying drawings, but the embodiments of the present invention are not limited thereto.

[0090] like Figure 1 The method for automatic history fitting of water-drive reservoirs at the full grid scale based on perturbation pulse inversion, as shown, specifically includes the following steps:

[0091] Step S1: Based on the streamline simulation framework, the seepage field is discretized into multiple continuous time intervals and corresponding streamline sets;

[0092] Step S2: Solve the aqueous phase seepage equation using the method of characteristics to obtain the governing equation for the propagation of saturation along streamlines;

[0093] Step S3: Based on the governing equations, construct the perturbation transmission mechanism of water saturation micro-perturbation, establish its linear mapping relationship with the water cut error of the production well, and form a system of linear equations;

[0094] Step S4: Using the historical fitting error of water content as a constraint, construct and solve the least squares optimization problem to obtain the optimal correction amount of the water saturation field;

[0095] Step S5: Update the water saturation field using the optimal correction amount and generate a new streamline distribution;

[0096] Step S6: Repeat steps S2 to S5 for iteration until the moisture content fitting meets the accuracy requirements and the saturation field history fitting is achieved.

[0097] Step S7: Correct the permeability field based on the time series correlation coefficient between injection and production wells and the pressure prediction error of production wells.

[0098] Specifically, step S1 includes the following steps:

[0099] Define a continuous time discrete point sequence The entire simulation time is divided into a continuous time interval ;

[0100] In each time interval Internally, using streamline simulation, the transient flow field is represented as a set of discrete streamlines starting from the injection well and ending at the production well. ;

[0101] Along each streamline Using arc length coordinates The space is parameterized, where The total length of the streamline;

[0102] Using the streamlines as a framework, the reservoir seepage field is discretized into a series of dynamic flow tube units, and the fluid flow behavior within each flow tube unit is dominated by the streamlines it is attached to.

[0103] In subsequent iterations, the streamline distribution and the flow tube unit division are dynamically updated based on the updated reservoir parameter field.

[0104] Specifically, step S2 includes the following steps:

[0105] Taking oil-water two-phase flow as an example, traditional numerical simulation techniques simulate the evolution of the seepage field by solving the mass conservation equation that controls the underground flow. The conservation equation is expressed as follows:

[0106] ;

[0107] ;

[0108] in, Porosity For crude oil density, Oil saturation For absolute penetration rate, The relative permeability of the oil phase. For crude oil viscosity, For the flow rate of oil, The density of water, Water saturation The relative permeability of the aqueous phase. The viscosity of water, For water flow rate, For time, Let the coordinates be the arc length along the streamline. For pressure;

[0109] Assuming incompressible flow, and defining the Darcy velocity of water as... Simplified to:

[0110] ;

[0111] in, Darcy's velocity for water;

[0112] Solve using the method of characteristics, let The saturation propagation equation along the streamline is obtained as follows:

[0113] ;

[0114] The flight time of a streamline is defined as:

[0115] ;

[0116] in, Let these be the initial arc length coordinates along the streamline. Let the coordinates be the arc length along the streamline. For position arrive Flight time, For the total Darcy speed, This is the flow rate function of the water.

[0117] Integrating the saturation equation along the streamlines, we obtain the result at position and time The expression for water saturation:

[0118] ;

[0119] in, For position and time water saturation For position And the initial water saturation at time 0, From position Arrive at the location The time.

[0120] Specifically, step S3 includes the following steps:

[0121] Based on the aforementioned saturation propagation equation, when for mesh cells Apply perturbation to the initial water saturation hour:

[0122] ;

[0123] in, To disturb the water saturation;

[0124] By utilizing simplified discretization It can identify grid cells that affect the water cut of a production well within a specific time interval; thus, the impact of saturation disturbance in any grid cell on the dynamic response of the production well can be quantitatively tracked.

[0125] Consider a location that starts from the injection well and terminates at the production well. streamline On the grid cells ;set up , making ;Applied to grid cells Disturbance on initial water saturation It will propagate along the streamline and in time Then arrived at the production well :

[0126] ;

[0127] in, Indicates from position Arrive at the location Flight time;

[0128] Through calculation It can be determined that the saturation disturbance is in the production well. The arrival time; by using the error between the predicted value and the observed water content, the grid cell can be derived. Optimal correction for water saturation Flight time is:

[0129] ;

[0130] In time interval Inside, for those passing through grid cells And eventually reach the production well streamline , applied to disturbances on Will have an impact on time of The moisture content at that location; the resulting change in moisture content is ;set up This indicates a line whose endpoint is a production well. The streamlines; then, for any time... and production wells The following relationships hold true:

[0131] ;

[0132] in, It is the most recent streamline start time. To flow through the grid cells Location To the end of the streamline Flight time, For grid cells water saturation For grid cells Disturbance water saturation, For production wells In time The moisture content fitting error.

[0133] Specifically, step S4 includes the following steps:

[0134] Based on the linear relationship established in step S3, the equations for all times and all production wells are combined to construct a linear system. ,in:

[0135] This is the water saturation correction vector for the grid cells to be solved;

[0136] The coefficient matrix is ​​expressed as follows:

[0137] ;

[0138] in, For time point indexing, For the number of production wells, For production well index, row index Corresponding to production time and production wells Column index Corresponding grid cell , It is the most recent streamline start time. To flow through the grid cells Location To the end of the streamline Flight time, For grid cells water saturation For a line whose endpoint is a production well Streamlines;

[0139] Let the moisture content prediction error vector be expressed as follows:

[0140] ;

[0141] in, For the first production well In the Time points The error in moisture content prediction;

[0142] Establish and solve the following constrained least squares optimization problem to obtain the optimal correction amount. :

[0143] ;

[0144] in, This is the physical lower limit of water saturation. This represents the physical upper limit of water saturation. For grid cells The current water saturation value, The mesh element to be solved Water saturation correction vector This is a set of grid cells for which changes in water saturation are not permitted.

[0145] Specifically, step S5 includes the following steps:

[0146] The optimal correction vector Each element in With the corresponding grid cell The current water saturation values ​​are added together to update the water saturation field in the reservoir model. Then, based on the updated water saturation field and the current reservoir pressure field, streamline simulation is performed again to generate a new set of discrete streamlines that reflects the current seepage field state, providing a basis for the next iteration.

[0147] Specifically, step S6 includes the following steps:

[0148] Using the updated water saturation field and the new streamline distribution as input, steps S2, S3, S4, and S5 are repeated to form a closed iterative optimization loop. After each iteration, the historical fitting error of water cut for all production wells is calculated, and it is determined whether the preset accuracy convergence criterion has been met. If the accuracy requirement is met, the iteration terminates, and the final inverted water saturation field is output. If not, the next iteration continues until convergence.

[0149] Specifically, step S7 includes the following steps:

[0150] The Pearson correlation coefficient is used to represent the influence between two time series, and its expression is:

[0151] ;

[0152] in, For two frequency domains and The correlation coefficient between them Representing the frequency domain and frequency domain covariance, For frequency domain The average value, For frequency domain The average value, For frequency domain standard deviation For frequency domain standard deviation Represents the expectation operator;

[0153] Calculation of water injection wells With production wells Correlation coefficient between injection volume and product volume time series :

[0154] ;

[0155] in, For water injection wells With production wells The correlation coefficient, The lag coefficient, For injection well The injection volume time series, from time index 1 to Lagging data, For production wells The time series of liquid production, from the time index arrive ;

[0156] Combined with production wells Bottom hole pressure prediction error For connecting water injection wells With production wells streamline Grid cells passed through The penetration rate is corrected:

[0157] ;

[0158] in, This is the corrected permeability. The permeability before correction. and To adjust the coefficient, For production wells The error in bottom hole pressure prediction To inject from the well Departure and in time Arrival at the production well The streamline set, For grid cells Located on one of the streamlines.

[0159] The feasibility and advantages of the present invention are demonstrated by the following specific embodiments.

[0160] The embodiments of this invention employ a realistic three-dimensional geological model, such as... Figure 2 As shown, the model's grid system is 370×252×12, with a total of 1,118,880 grid cells, and the simulated study area is approximately 90.26 km². 2 The simulation includes 755 production wells and 251 water injection wells, totaling 1006 wells, resulting in a complex well network. The reservoir's physical properties exhibit strong heterogeneity, with porosity ranging from 2.3% to 18.9% and an average porosity of 14.4%. Permeability ranges from 0.015 mD to 179.2 mD, with an average permeability of 15.7 mD, classifying it as a typical low-porosity, low-permeability reservoir. This block has been under development since January 1995, with continuous water injection development. It possesses nearly thirty years of dynamic production history data from the start of production to February 2024, including monthly oil production, water production, water injection volume, and periodic pressure measurement data for all wells. The data is massive in scale and contains complex dynamic response information.

[0161] Step S1: Discretization of the seepage field based on the streamline simulation framework; The total simulation time is divided into monthly time intervals. At the beginning of each time interval, the pressure field is obtained through full reservoir simulation based on that time. The streamline simulation method is used to generate a streamline set from all 251 injection wells to all 755 production wells. Each streamline is parameterized by arc length, and the entire seepage field is dynamically discretized into associated flow tube units based on the generated streamlines, resulting in the streamline set and flow tube unit division results for subsequent analysis.

[0162] Step S2: Derivation of the control equation for saturation propagation along the streamline; considering the development characteristics of water-driven reservoirs, the influence of the gas phase is ignored, and the incompressible flow of the oil and water phases is taken into account; based on the streamlines and the total velocity distribution along the streamlines obtained in Step S1, the three-dimensional two-phase seepage partial differential equations are simplified by applying the method of characteristics; by calculating the flight time between any two points on the streamline through integration, a deterministic relationship is established between the water saturation at a certain point at the current moment and its upstream historical saturation, thus obtaining the control equation for saturation propagation along the streamline.

[0163] Step S3: Establish a linear mapping between saturation disturbance and production response; based on the saturation propagation mechanism along streamlines established in Step S2, perform disturbance analysis, with the saturation disturbance propagating downstream along all streamlines passing through the grid cell; by calculating the flight time of the disturbance along each streamline to the corresponding production well, determine the historical time point at which the disturbance affects the water cut of the production well; under the assumption that the disturbance is a small amount, there is a linear relationship between the saturation disturbance and the water cut response of the production well, that is, the change in the water cut of the production well is proportional to the amount of saturation disturbance in the grid cell; traverse the influence relationship of all 1,118,880 grid cells of the entire reservoir on 755 production wells at 349 time points, and construct a large sparse linear equation system reflecting the linear mapping relationship between saturation disturbance and water cut error of the production well.

[0164] Step S4: Solve for the optimal correction amount of the saturation field; Based on the linear equation system constructed in step S3, the nonlinear inversion problem with the objective function of minimizing the historical water content fitting error is transformed into a least squares optimization problem with constraints; By solving the large sparse matrix, the globally optimal saturation correction vector is obtained to characterize the adjustment direction and magnitude required for each grid cell to match the historical water content data.

[0165] Step S5: Update the saturation field and reconstruct streamlines; Apply the optimal correction vector obtained in step S4 to the water saturation field of the current reservoir model, update the saturation value of each grid cell, and obtain the updated water saturation field of the whole area; Recalculate the reservoir pressure distribution based on the updated water saturation field, and repeat step S1 based on the new pressure field to generate a new set of streamlines that matches the current seepage state.

[0166] Step S6: Iterate and optimize until the water content fitting converges; repeat steps S2 to S5 to form a closed iterative optimization loop; in this embodiment, even starting from a physically unreasonable initial field of full water saturation, the water saturation field obtained after several iterations quickly converges to a reasonable distribution; such as Figure 3 As shown, the saturation field obtained after 10 iterations clearly reveals the complex three-dimensional water-flooded structure formed after nearly 30 years of waterflooding. The distribution of high water-cut channels and remaining oil-rich areas is consistent with geological understanding; the overall water cut fitting curve for the entire area is as follows: Figure 6 As shown, the iterative fitting curve quickly approximates the actual observed data and reaches the preset convergence standard.

[0167] Step S7: Coordinated Correction of the Permeability Field; Based on the obtained water saturation field inversion results, the permeability field is optimized using pressure data and injection-production dynamic correlation; Pearson correlation coefficients of injection and production volume time series between injection and production wells are calculated to identify major seepage channels; combined with the long-term fitting error of bottom hole pressure in production wells, the permeability of grid cells traversed by highly correlated streamline paths is fine-tuned under physical constraints; The corrected permeability field, while maintaining geostatistical characteristics, shows higher consistency between its high-permeability bands and the dominant seepage channels revealed by dynamic inversion, such as... Figure 4 The initial permeability model shown is... Figure 5 The comparison of historical fitting penetration models is shown below.

[0168] Figure 7 The historical fitting effect of the simulated liquid production and oil production in the entire region was demonstrated, and the simulated curves matched the observed data well. This study also verified the ability of the method of this invention to handle complex noise and abrupt changes in actual data. Figures 8 to 13 The fitting results of water cut, fluid production and oil production of six production wells with different dynamic characteristics are presented respectively. Even when there are fluctuations and step-like abrupt changes in the data, the method of the present invention can still accurately track the changing trend and achieve high-precision historical matching.

[0169] This embodiment demonstrates that, for actual oil reservoirs with over a million grid cells, over a thousand wells, a production history of nearly thirty years, and complex geological conditions, the automatic history fitting method for water-driven oil reservoirs at the full grid scale based on perturbation pulse inversion provided by this invention can be implemented efficiently and stably. This method, through linearized iteration under a streamlined framework and multi-parameter collaborative inversion, successfully achieves high-precision automated fitting of the long-term production dynamics of large-scale oil reservoirs, breaking through the bottlenecks of traditional methods in terms of computational efficiency, initial model dependence, and inversion accuracy, and providing key technical support for the fine description and intelligent decision-making of oil fields.

[0170] The above description is merely a preferred embodiment of the present invention, but the scope of protection of the present invention is not limited thereto. Any variations or substitutions that can be easily conceived by those skilled in the art within the scope of the embodiments disclosed in the present invention should be included within the scope of protection of the present invention; therefore, the scope of protection of the present invention should be determined by the scope of the claims.

Claims

1. A method for automatic history fitting of water-drive reservoirs at the full grid scale based on perturbation pulse inversion, characterized in that, Includes the following steps: Step S1: Based on the streamline simulation framework, the seepage field is discretized into multiple continuous time intervals and corresponding streamline sets; Step S2: Solve the aqueous phase seepage equation using the method of characteristics to obtain the governing equation for the propagation of saturation along streamlines; The specific process of step S2 is as follows: The mass conservation equation is expressed as follows: ； ； in, Porosity For crude oil density, Oil saturation For absolute penetration rate, The relative permeability of the oil phase. For crude oil viscosity, For the flow rate of oil, The density of water, Water saturation The relative permeability of the aqueous phase. The viscosity of water, For water flow rate, For time, Let the coordinates be the arc length along the streamline. For pressure; Assuming incompressible flow, and defining the Darcy velocity of water as... Simplified to: ； in, Darcy's velocity for water; Solve using the method of characteristics, let The saturation propagation equation along the streamline is obtained as follows: ； The flight time along the streamline is defined as: ； in, Let these be the initial arc length coordinates along the streamline. Let the coordinates be the arc length along the streamline. For position arrive Flight time, For total Darcy speed, This is the flow rate function of the water. Integrating the saturation equation along the streamlines, we obtain the result at position and time The expression for water saturation: ； in, For position and time water saturation For position And the initial water saturation at time 0, From position Arrive at the location Time; Step S3: Based on the governing equations, construct the perturbation transmission mechanism of water saturation micro-perturbation, establish its linear mapping relationship with the water cut error of the production well, and form a system of linear equations; Step S4: Using the historical fitting error of water content as a constraint, construct and solve the least squares optimization problem to obtain the optimal correction amount of the water saturation field; Step S5: Update the water saturation field using the optimal correction amount and generate a new streamline distribution; Step S6: Repeat steps S2 to S5 for iteration until the moisture content fitting meets the accuracy requirements and the saturation field history fitting is achieved. Step S7: Correct the permeability field based on the time series correlation coefficient between injection and production wells and the pressure prediction error of production wells.

2. The automatic history fitting method for water-drive reservoirs at the full grid scale based on perturbation pulse inversion according to claim 1, characterized in that, The specific process of step S1 is as follows: Define a continuous time discrete point sequence The entire simulation time is divided into a continuous time interval ; In each time interval Internally, using streamline simulation, the transient flow field is represented as a set of discrete streamlines starting from the injection well and ending at the production well. ; Along each streamline Using arc length coordinates The space is parameterized, where This is the total length of the streamline; Using the streamlines as a framework, the reservoir seepage field is discretized into a series of dynamic flow tube units, and the fluid flow behavior within each flow tube unit is dominated by the streamlines it is attached to. In subsequent iterations, the streamline distribution and the flow tube unit division are dynamically updated based on the updated reservoir parameter field.

3. The automatic history fitting method for water-drive reservoirs at the full grid scale based on perturbation pulse inversion according to claim 2, characterized in that, The specific process of step S3 includes: Based on the aforementioned saturation propagation equation, when for mesh cells Apply perturbation to the initial water saturation hour: ； in, To disturb the water saturation; Consider a location that starts from the injection well and terminates at the production well. streamline On the grid cells ;set up , making ;Applied to grid cells Disturbance on initial water saturation It will propagate along the streamline and in time Then arrived at the production well : ； in, Indicates from position Arrive at the location Flight time; Through calculation Determine the saturation disturbance in the production well The arrival time is determined; the error between the predicted value and the observed water content is used to derive the grid cell. Optimal correction for water saturation Flight time is: ； In time interval Inside, for those passing through grid cells And eventually reach the production well streamline , applied to disturbances on Will have an impact on time of The moisture content at that location; the resulting change in moisture content is ;set up This indicates a line whose endpoint is a production well. The streamlines; then, for any time... and production wells The following relationships hold true: ； in, It is the most recent streamline start time. To flow through the grid cells Location To the end of the streamline Flight time, For grid cells water saturation For grid cells Disturbance water saturation, For production wells In time The moisture content fitting error.

4. The automatic history fitting method for water-drive reservoirs at the full grid scale based on perturbation pulse inversion according to claim 3, characterized in that, The specific process of step S4 is as follows: Based on the linear relationship established in step S3, the equations for all times and all production wells are combined to construct a linear system. ,in: This is the water saturation correction vector for the grid cells to be solved; The coefficient matrix is ​​expressed as follows: ； in, For time point indexing, For the number of production wells, For production well index, row index Corresponding to production time and production wells Column index Corresponding grid cell , It is the most recent streamline start time. To flow through the grid cells Location To the end of the streamline Flight time, For grid cells water saturation For a line whose endpoint is a production well Streamlines; Let the moisture content prediction error vector be expressed as follows: ； in, For the first production well In the Time points The error in moisture content prediction; Establish and solve the following constrained least squares optimization problem to obtain the optimal correction amount. : ； in, This is the physical lower limit of water saturation. This represents the physical upper limit of water saturation. For grid cells The current water saturation value, The mesh element to be solved Water saturation correction vector This is a set of grid cells for which changes in water saturation are not permitted.

5. The automatic history fitting method for water-drive reservoirs at the full grid scale based on perturbation pulse inversion according to claim 4, characterized in that, The specific process of step S5 is as follows: The optimal correction vector Each element in With the corresponding grid cell The current water saturation values ​​are added together to update the water saturation field in the reservoir model; Subsequently, based on the updated water saturation field and the current reservoir pressure field, streamline simulation was re-executed to generate a new set of discrete streamlines that reflects the current seepage field state, providing a basis for the next iteration.

6. The automatic history fitting method for water-drive reservoirs at the full grid scale based on perturbation pulse inversion according to claim 5, characterized in that, The specific process of step S6 is as follows: Using the updated water saturation field and the new streamline distribution as input, steps S2, S3, S4 and S5 are repeated to form a closed iterative optimization loop; After each iteration, the historical fitting error of water cut for all production wells is calculated, and it is determined whether the preset accuracy convergence criterion has been met. If the accuracy requirement is met, the iteration terminates and the final inverted water saturation field is output; if not, the next iteration continues until convergence.

7. The automatic history fitting method for water-drive reservoirs at the full grid scale based on perturbation pulse inversion according to claim 6, characterized in that, The specific process of step S7 is as follows: The Pearson correlation coefficient is used to represent the influence between two time series, and its expression is: ； in, For two frequency domains and The correlation coefficient between them Representing the frequency domain and frequency domain covariance, For frequency domain The average value, For frequency domain The average value, For frequency domain standard deviation For frequency domain standard deviation Represents the expectation operator; Calculation of water injection wells With production wells Correlation coefficient between injection volume and product volume time series : ； in, For water injection wells With production wells The correlation coefficient, The lag coefficient, For injection well The injection volume time series, from time index 1 to Lagging data, For production wells The time series of liquid production, from the time index arrive ; Combined with production wells Bottom hole pressure prediction error For connecting water injection wells With production wells streamline Grid cells passed through The penetration rate is corrected: ； in, This is the corrected permeability. The permeability before correction. and To adjust the coefficient, For production wells The error in bottom hole pressure prediction To inject from the well Departure and in time Arrival at the production well The streamline set, For grid cells Located on one of the streamlines.

Images