A high-resolution dipole flexural wave velocity inversion method

Abstract

The present invention discloses a high-resolution dipole flexural wave velocity inversion method for oil and gas exploration in unconventional reservoirs such as shale and tight sandstone. The method comprises collecting conventional well logging data and array acoustic logging data from wells and performing curve preprocessing on them; processing the preprocessed monopole array and dipole array waveform data using a correlation method to obtain low-resolution P-wave velocity, P-wave arrival time, and low-resolution flexural wave velocity; establishing a numerical model based on the acquired data and performing forward simulation to obtain dipole flexural waveforms, which are then analyzed in the frequency-wavenumber domain to determine the high-order flexural wave frequency range; filtering the measured dipole array waveform data within this range and performing correlation analysis constrained by the flexural wave arrival time to obtain high-resolution dipole flexural wave velocity. The method can improve the vertical resolution of flexural wave velocity to 0.3 m, making the logging interpretation of thin interbeds more accurate and better understanding the characteristics of oil reservoirs.

Description

Technical Field

[0001] The present invention relates to the field of oil and gas exploration in unconventional reservoirs such as shale and tight sandstone, and more specifically to a high-resolution dipole bending wave velocity inversion method. Background Art

[0002] With the continuous development and deepening of oil and gas field exploration and development, the requirements for logging technology and interpretation are becoming higher and higher. High-resolution acoustic logging plays an important role in the division of thin layers and subdivision of thick layers. It also plays a key role in the evaluation of future logging data and the re-evaluation of logging data from old wells. When evaluating formations, it is often necessary to estimate the acoustic parameters and rock physical properties of horizontal thin layers in order to better understand the characteristics and reserves of oil reservoirs. The conventional processing of array acoustic logging data is limited by the structure of the instrument and usually gives an average slowness curve within the span of the receiving transducer array. Its highest resolution is the span between the first receiving transducer and the last receiving transducer. This often ignores the characteristics of formations with a span smaller than the array span, which is not conducive to the logging interpretation of thin interbeds.

[0003] In 1987, Hsu and Chang proposed using a multi-source waveform similarity superposition method to process repeated subarray data to improve the resolution of slowness curves. In 2000, Tang XM et al. proposed a signal processing method for evaluating thin formations using high-resolution acoustic slowness logging, improving the resolution from 3.5 ft to 0.5 ft. They later proposed a dual-waveform matching method that uses subarrays to combine different formation spans, achieving compressional wave velocity extraction with a resolution of 0.5 ft. In 2005, Li Pengju et al. proposed a correlation spectrum combined with multi-shot processing method. This method uses the correlation spectrum method to extract the monopole compressional wave velocity of the shortest co-transmitted and co-received subarrays in multipole array acoustic logging. The compressional wave velocities of all the shortest co-transmitted or co-received subarrays across the same depth are then averaged to obtain the co-transmitted or co-received compressional wave velocity, respectively. In 2019, T. Lei et al. proposed a new interpretation algorithm. The algorithm, based on a powerful downscaling technique, extracts all convolution relationships between acoustic well logs of different resolutions to generate an overdetermined matrix. Velocity is estimated using the Moore-Penrose pseudo-inverse method to improve the vertical resolution of the logs. In 2020, Zhou Haoyi proposed combining bisection and simulated annealing with the slowness-time method to extract component wave velocities. This method improves upon the traditional slowness-time correlation method in terms of processing efficiency and accuracy.

[0004] The aforementioned research primarily focused on extracting acoustic wave velocity from array acoustic wave data, but failed to consider how to leverage the depth-offset multidimensional information and high-frequency components in waveform data for high-resolution acoustic wave velocity extraction. In particular, it failed to consider extracting high-resolution dipole flexural wave velocities, a modal wave close to formation shear wave velocity that has practical applications in analyzing soft elastic reservoirs such as shale. No related techniques have been identified in the prior art. In response, the present invention proposes a high-resolution dipole flexural wave velocity inversion method that first determines the frequency range of high-order flexural waves using forward modeling and then uses this high-order flexural wave frequency range to perform high-resolution flexural wave velocity inversion. Summary of the Invention

[0005] The purpose of the present invention is to provide a high-resolution dipole flexural wave velocity inversion method. This method makes the logging interpretation of thin interbeds more accurate, so as to better understand the characteristics and reserves of oil reservoirs. It provides more effective data support for measuring reservoirs with soft elastic properties and has practical application value.

[0006] The present invention provides a high-resolution dipole bending wave velocity inversion method, comprising the following steps:

[0007] S1. Collect conventional logging data and array acoustic logging data of the well;

[0008] S2. Performing curve preprocessing on conventional logging data and array acoustic logging data to obtain monopole array waveform data and dipole array waveform data;

[0009] S3, using a velocity extraction method to process the monopole array waveform data to obtain low-resolution longitudinal wave velocity and longitudinal wave arrival time, and using a velocity extraction method to process the dipole array waveform data to obtain low-resolution flexural wave velocity;

[0010] S4. Based on the wellbore diameter, density, compressional wave velocity, and flexural wave velocity, a numerical model of the well is established and a forward simulation is performed to obtain a simulated dipole flexural waveform;

[0011] S5. Perform frequency-wavenumber domain analysis on the simulated dipole bending waveform to obtain the frequency distribution range of the high-order bending wave;

[0012] S6. Selecting a high-order bending wave frequency distribution range to perform filtering processing on the measured dipole array waveform data to obtain high-order bending wave array waveform data;

[0013] S7. Perform correlation analysis of the high-order flexural wave array waveform data with flexural wave arrival time constraints to obtain high-resolution dipole flexural wave velocities.

[0014] In S1, the conventional logging data of the well mainly include natural gamma, well diameter, and density curve; the array acoustic logging data of the well mainly include the original waveform, gain curve and delay curve of the monopole array acoustic logging and the original waveform, gain curve and delay curve of the dipole array acoustic logging.

[0015] In S2, the curve preprocessing of array acoustic logging data mainly includes gain recovery and delay recovery. Gain recovery is performed first, and then delay recovery. During the waveform acquisition process, in order to reduce the storage volume of recorded data, waveform data for a period of time before the arrival of the first wave is often not collected. The purpose of delay recovery is to fill zero before the zero moment of the monopole array acoustic waveform, so as to obtain waveform data with accurate time.

[0016] In S3, the velocity extraction method is a correlation method. The correlation method is selected for extracting a P-wave velocity curve from the preprocessed monopole array acoustic waveform data. The correlation function of the array waveform is calculated in two dimensions, time and velocity. The velocity at which the maximum function value is obtained is the low-resolution P-wave velocity. Based on the low-resolution P-wave velocity, the P-wave arrival time curve can be directly calculated taking into account the well diameter and the length of the logging tool.

[0017] The correlation method is used to extract the shear wave velocity curve from the preprocessed dipole array acoustic waveform data. The correlation function of the array waveform is calculated in the two dimensions of time and velocity. The velocity with the maximum function value is the low-resolution dipole bending wave velocity.

[0018] In said S4, a comprehensive analysis is performed on the wellbore diameter, density, P-wave velocity and dipole bending wave velocity curves of the well to find a layer section with no wellbore diameter collapse and no sudden changes in density, P-wave velocity and dipole bending wave velocity;

[0019] The density, wellbore diameter, compressional wave velocity and dipole flexural wave velocity of the interval are taken as input parameters, and a numerical simulation method is used to perform forward simulation and calculate the dipole flexural wave array waveform.

[0020] The numerical simulation method adopts real axis integration method or finite difference method.

[0021] In said S5, said frequency-wavenumber domain analysis mainly involves performing a two-dimensional Fourier transform on the simulated dipole bending waveform data to obtain a frequency-wavenumber two-dimensional spectrum corresponding to the simulated dipole bending wave array waveform;

[0022] In the frequency-wavenumber two-dimensional spectrum, the horizontal axis is set to wavenumber and the vertical axis is set to frequency; from it, a first-order bending wave mode with a lower frequency and a second-order bending wave mode with a higher frequency can be observed; the frequency band distribution range corresponding to the second-order bending wave mode is the frequency distribution range of the higher-order bending wave.

[0023] In S6, the filtering process refers to setting a digital bandpass filter and processing the dipole bending wave array waveform data measured in S4 with the digital bandpass filter to filter out low-order bending wave information, and the remaining waveform data is high-order bending wave array waveform data.

[0024] The frequency band passing range of the digital bandpass filter is set to the frequency band distribution range of the high-order dipole bending wave obtained in S5.

[0025] In S7, the correlation analysis of the bending wave arrival time constraint is an improvement on the existing correlation method. Specifically, the windowing start time is no longer a fixed value, but is instead based on the dipole bending wave arrival time curve. This is done to reduce the impact of noise.

[0026] Compared with the prior art, the advantages of the present invention are:

[0027] 1. This paper proposes a new method for determining the frequency range of high-order flexural waves using forward modeling, and then using this high-order flexural wave frequency range to perform high-resolution flexural wave velocity inversion. Compared with existing acoustic velocity extraction methods, this method can improve the vertical resolution of flexural wave velocity to 0.3 m, enabling more accurate logging interpretation of thin interbeds, and better understanding of oil reservoir characteristics and reserves. It also provides more effective data support for measuring the elastic properties of soft reservoirs, and has practical application value. BRIEF DESCRIPTION OF THE DRAWINGS

[0028] Figure 1 The original well logging data of Well Y and the pre-processed array acoustic logging data provided by the embodiment of the present invention are shown in FIG.

[0029] Figure 2 The simulated dipole bending wave waveform corresponding to the forward model constructed based on the Y well logging curve provided in the embodiment of the present invention;

[0030] Figure 3 A two-dimensional spectrum in the frequency-wavenumber domain corresponding to the simulated dipole bending wave waveform provided by an embodiment of the present invention;

[0031] Figure 4 A comparison diagram of the high-resolution dipole bending wave velocity and the low-resolution dipole bending wave velocity provided by an embodiment of the present invention. DETAILED DESCRIPTION

[0032] In order to make the technical solution of the present invention more clear, the present invention is further described in detail below with reference to the accompanying drawings and through examples.

[0033] A high-resolution dipole bending wave velocity inversion method comprises the following steps:

[0034] S1. Collect conventional logging data and array acoustic logging data of the well, mainly including natural gamma, wellbore diameter, and density curves; array acoustic logging data of the well, mainly including the original waveform, gain curve and delay curve of monopole array acoustic logging and the original waveform, gain curve and delay curve of dipole array acoustic logging.

[0035] Figure 1 The first track is the natural gamma of Well Y; the second track is the well diameter; the third track is the density; the fourth track is the depth track; the fifth track is the original monopole waveform density map WVPO, from which the longitudinal wave information can be observed; the eighth track is the original dipole bending waveform density map WVFLO, from which the shear wave information can be observed.

[0036] S2. Curve preprocessing is performed on conventional logging data and array acoustic logging data to obtain caliper curves, density curves, monopole array waveform data, and dipole array waveform data. Monopole array acoustic data preprocessing primarily involves gain recovery and delay recovery, with gain recovery performed first, followed by delay recovery. During waveform acquisition, waveform data for a period before the arrival of the first wave is often not collected to reduce the amount of recorded data stored. The purpose of delay recovery is to fill in zeros before the zero time of the monopole array acoustic waveform, thereby obtaining waveform data with accurate time.

[0037] Figure 1 The sixth track is the monopolar waveform density map WVP after preprocessing, from which it can be observed that the longitudinal wave signal becomes clearer after waveform preprocessing; the ninth track is the dipole bending wave density map WVFL after preprocessing, from which it can be observed that the dipole bending wave signal becomes clearer after waveform preprocessing.

[0038] S3. Use the velocity extraction method to process the monopole array waveform data to obtain low-resolution longitudinal wave velocity and longitudinal wave arrival time, and use the velocity extraction method to process the dipole array waveform data to obtain low-resolution bending wave velocity.

[0039] The correlation method is used to extract the P-wave velocity curve from the preprocessed monopole array acoustic waveform data. This method calculates the correlation function of the array waveform in both time and velocity dimensions, with the velocity at the maximum function value representing the low-resolution P-wave velocity. Based on the low-resolution P-wave velocity, the P-wave arrival time curve can be directly calculated, taking into account the wellbore diameter and logging tool length. Similarly, the correlation method is used to extract the S-wave velocity curve from the preprocessed dipole array acoustic waveform data. This method calculates the correlation function of the array waveform in both time and velocity dimensions, with the velocity at the maximum function value representing the low-resolution dipole flexural wave velocity.

[0040] Figure 1Track 7 shows a low-resolution longitudinal wave velocity curve, which is the result of processing the monopole array waveform data using the velocity extraction method. Track 10 shows a low-resolution dipole flexural wave velocity curve, which is the result of processing the dipole flexural wave array waveform data using the velocity extraction method. It can be seen from the curve that the longitudinal resolution of the dipole flexural wave velocity is low.

[0041] S4. Based on the wellbore diameter, density, and compressional and flexural wave velocities, a numerical model of the well is established and forward simulation is performed to obtain a simulated dipole flexural waveform.

[0042] Taking into account various conditions such as no expansion of the wellbore diameter, no sudden changes in density, P-wave velocity and dipole flexural wave velocity, the wellbore diameter, density, P-wave velocity and dipole flexural wave velocity at a depth of 2322 m were selected to establish a numerical model of the wellbore formation. This embodiment uses the real axis integration method to implement forward simulation. Usually, a real array acoustic logging instrument will measure eight dipole flexural wave waveforms at each depth position, which have different source distances. In the numerical simulation process, these source distances should also be used as input parameters. The simulation results are shown in Figure 2. Figure 2 As shown in the figure, the horizontal axis is time, and the vertical axis is the corresponding source distance of the eight dipole bending wave array waveforms, with the minimum source distance being 10.25 ft and the maximum source distance being 13.75 ft. The figure shows pure dipole bending waveforms, whose arrival time is delayed as the source distance increases.

[0043] S5. Perform frequency-wavenumber domain analysis on the simulated dipole bending waveform to obtain the frequency distribution range of the high-order bending wave.

[0044] against Figure 2 The simulated dipole bending wave array waveform shown in the figure is subjected to double Fourier transform to obtain a two-dimensional spectrum in the frequency-wavenumber domain, as shown in Figure 3 As shown. The horizontal axis is the wave number and the vertical axis is the frequency. The oblique line VFL in the figure represents the projection of the shear wave velocity in the two-dimensional spectrum. The box in the figure is the first-order dipole flexural wave mode, from which it can be observed that its intersection with the shear wave velocity curve is at around 5kHz. Therefore, if the second-order dipole flexural wave (such as Figure 3 The position indicated by the arrow in the middle) or higher-order dipole bending waves need to be filtered out.

[0045] S6. Select a high-order bending wave frequency distribution range to perform filtering processing on the measured dipole array waveform data to obtain high-order bending wave array waveform data.

[0046] Figure 4 The first channel is the depth channel; the second channel is the pre-processed dipole bending wave density map; the third channel is the frequency corresponding to the second waveform, from which we can see that the sound wave frequency band is relatively wide, with response characteristics in the range of 3-8kHz. Figure 3The analysis results show that to extract high-order dipole flexural waves, waveform signals below 5 kHz must be filtered out. Taking into account the bandwidth broadening of the measured original dipole flexural wave waveform, a digital bandpass filter with a passband frequency range of 6-8 kHz is designed to obtain pure high-order dipole flexural waves. Figure 4 The fourth track is the dipole bending wave after digital bandpass filtering based on the second track; the fifth track is the corresponding spectrum of the fourth track. It can be observed that the spectrum components below 5kHz are indeed filtered out, and only the spectrum energy above 5kHz is retained, which corresponds to the high-order dipole bending wave.

[0047] S7. Perform correlation analysis of the high-order flexural wave array waveform data with flexural wave arrival time constraints to obtain high-resolution dipole flexural wave velocities.

[0048] A correlation analysis with first-order flexural wave arrival time constraints is performed on the waveform data of the high-order flexural wave array. In particular, the arrival time constraint correlation analysis method used is an improvement on the existing correlation method. Specifically, the window start time is no longer a fixed value, but is based on the dipole flexural wave arrival time curve. The purpose of this implementation is to reduce the influence of noise. After analysis, a high-resolution dipole flexural wave velocity curve is obtained, such as Figure 4 As shown in the seventh track, compared with the low-resolution velocity curve corresponding to the first-order flexural wave shown in the sixth track, the high-resolution flexural wave velocity curve shows more detailed information about the formation changes, especially in Figure 4 Arrows indicate depths of 2318.5m, 2332.5m and 2336m.

[0049] At this point, the entire high-resolution dipole bending wave velocity inversion method process is completed.

Claims

1. A high-resolution dipole bending wave velocity inversion method, characterized in that: The following steps are involved: S1. Collect conventional logging data and array acoustic logging data of the well; S2. Performing curve preprocessing on conventional logging data and array acoustic logging data to obtain monopole array waveform data and dipole array waveform data; S3, using a velocity extraction method to process the monopole array waveform data to obtain low-resolution longitudinal wave velocity and longitudinal wave arrival time, and using a velocity extraction method to process the dipole array waveform data to obtain low-resolution flexural wave velocity; S4. Based on the wellbore diameter, density, compressional wave velocity, and flexural wave velocity, a numerical model of the well is established and a forward simulation is performed to obtain a simulated dipole flexural waveform; S5. Perform frequency-wavenumber domain analysis on the simulated dipole bending waveform to obtain the frequency distribution range of the high-order bending wave; S6. Selecting a high-order bending wave frequency distribution range to perform filtering processing on the measured dipole array waveform data to obtain high-order bending wave array waveform data; S7. Perform correlation analysis of the high-order flexural wave array waveform data with flexural wave arrival time constraints to obtain high-resolution dipole flexural wave velocities.

2. A high-resolution dipole bending wave velocity inversion method according to claim 1, characterized in that: In S1, the conventional logging data of the well mainly include natural gamma, well diameter, and density curve; the array acoustic logging data of the well mainly include the original waveform, gain curve and delay curve of the monopole array acoustic logging and the original waveform, gain curve and delay curve of the dipole array acoustic logging.

3. The high-resolution dipole bending wave velocity inversion method according to claim 1, characterized in that: In S2, the curve preprocessing of array acoustic logging data mainly includes gain recovery and delay recovery. Gain recovery is performed first, and then delay recovery. During the waveform acquisition process, in order to reduce the storage volume of recorded data, waveform data for a period of time before the arrival of the first wave is often not collected. The purpose of delay recovery is to fill zero before the zero moment of the monopole array acoustic waveform, so as to obtain waveform data with accurate time.

4. The high-resolution dipole bending wave velocity inversion method according to claim 1, characterized in that: In S3, the velocity extraction method is a correlation method. The correlation method is selected for extracting a P-wave velocity curve from the preprocessed monopole array acoustic waveform data. The correlation function of the array waveform is calculated in two dimensions, time and velocity. The velocity at which the maximum function value is obtained is the low-resolution P-wave velocity. Based on the low-resolution P-wave velocity, the P-wave arrival time curve can be directly calculated taking into account the well diameter and the length of the logging tool. The correlation method is used to extract the shear wave velocity curve from the preprocessed dipole array acoustic waveform data. The correlation function of the array waveform is calculated in the two dimensions of time and velocity. The velocity with the maximum function value is the low-resolution dipole bending wave velocity.

5. The high-resolution dipole bending wave velocity inversion method according to claim 1, characterized in that: In said S4, a comprehensive analysis is performed on the wellbore diameter, density, P-wave velocity and dipole bending wave velocity curves of the well to find a layer section with no wellbore diameter collapse and no sudden changes in density, P-wave velocity and dipole bending wave velocity; The density, wellbore diameter, compressional wave velocity and dipole flexural wave velocity of the interval are taken as input parameters, and a numerical simulation method is used to perform forward simulation and calculate the dipole flexural wave array waveform.

6. A high-resolution dipole bending wave velocity inversion method according to claim 5, characterized in that: The numerical simulation method adopts real axis integration method or finite difference method.

7. The high-resolution dipole bending wave velocity inversion method according to claim 1, characterized in that: In said S5, said frequency-wavenumber domain analysis mainly involves performing a two-dimensional Fourier transform on the simulated dipole bending waveform data to obtain a frequency-wavenumber two-dimensional spectrum corresponding to the simulated dipole bending wave array waveform; In the frequency-wavenumber two-dimensional spectrum, the horizontal axis is set to wavenumber and the vertical axis is set to frequency; from it, a first-order bending wave mode with a lower frequency and a second-order bending wave mode with a higher frequency can be observed; the frequency band distribution range corresponding to the second-order bending wave mode is the frequency distribution range of the higher-order bending wave.

8. The high-resolution dipole bending wave velocity inversion method according to claim 1, characterized in that: In S6, the filtering process refers to setting a digital bandpass filter, whose frequency band pass range is set to the frequency band distribution range of the high-order dipole bending wave obtained in S5, and using this digital bandpass filter to process the dipole bending wave array waveform data measured in S4, thereby filtering out low-order bending wave information, and the waveform data retained is the high-order bending wave array waveform data.

9. The high-resolution dipole bending wave velocity inversion method according to claim 1, characterized in that: In S7, the correlation analysis of the bending wave arrival time constraint is performed by setting the window start time to the time of the dipole bending wave arrival time curve based on the existing correlation method to reduce the impact of noise on the results.

Images