A carbonate rock pore structure identification method based on multi-modal deep learning

Abstract

The present application relates to the technical field of petroleum engineering drilling, and discloses a carbonate rock pore structure identification method based on multi-modal deep learning. Firstly, one-dimensional acoustic waveform data and porosity parameters obtained are normalized, noise is enhanced, and continuous wavelet transform is performed to obtain two types of inputs of one-dimensional acoustic curve and two-dimensional time-frequency diagram; then one-dimensional convolutional neural network is used to extract acoustic time domain local phase and amplitude features, global energy distribution features in the two-dimensional time-frequency diagram are extracted, and the porosity parameters are embedded as physical constraint features; three types of features are spliced to form a joint feature vector, three types of pore structure probabilities of cracks, intergranular pores and caves are output through a classification branch, and a continuous numerical value of a pore aspect ratio is output through a regression branch. Through multi-modal fusion of one-dimensional acoustic curve, two-dimensional time-frequency diagram and porosity physical parameters, the present application establishes the relationship between acoustic features and pore structure types, and realizes the identification of complex carbonate rock pore structure.

Description

Technical Field

[0001] This invention relates to the field of petroleum engineering drilling technology, and in particular to a method for identifying the pore structure of carbonate rocks based on multimodal deep learning. Background Technology

[0002] Carbonate reservoirs are an important type of oil and gas reservoir, and their internal pore spaces typically include various forms such as fractures, intergranular pores, and caverns. Compared with conventional clastic rocks, carbonate rocks exhibit stronger heterogeneity, and differences in pore morphology significantly affect the rock's equivalent elastic modulus, acoustic propagation velocity, and acoustic scattering attenuation characteristics. Therefore, accurately identifying the pore structure type of carbonate rocks and quantitatively inverting the pore aspect ratio is of great significance for reservoir evaluation and well logging interpretation.

[0003] In existing technologies, pore structure is typically characterized using experimental methods such as core mercury intrusion porosimetry, nuclear magnetic resonance (NMR), scanning electron microscopy (SEM), CT scans, and thin-section image analysis. While these methods can provide relatively intuitive information about pore structure, they are usually dependent on core samples and experimental conditions, making them difficult to directly apply to continuous well logging interpretation. Acoustic logging or core ultrasonic testing can reflect the influence of rock pore structure on elastic wave propagation, but traditional acoustic velocity analysis methods are prone to the problem of "same wave velocity, different pore structures," making it difficult to distinguish between different pore types such as fractures, intergranular pores, and caverns.

[0004] In recent years, intelligent algorithms such as convolutional neural networks, recurrent neural networks, and Transformers have been used for the interpretation of acoustic signals. A single one-dimensional acoustic wave curve model can capture local time-domain features such as the first arrival time, initial amplitude, and wake complexity, but it struggles to fully reflect the frequency components and global energy evolution of the acoustic wave. A single two-dimensional time-frequency diagram model can capture the distribution of acoustic energy across time and frequency scales, but the visualization process easily smooths out minute variations in the first arrival time. Therefore, existing single-modal intelligent models still suffer from incomplete feature information, high errors in identifying boundary pore structures such as intergranular pores, and insufficient accuracy in continuous inversion of pore aspect ratios. Summary of the Invention

[0005] The purpose of this invention is to address the aforementioned deficiencies in existing technologies by providing a method for identifying the pore structure of carbonate rocks based on multimodal deep learning. This invention establishes the relationship between acoustic features and pore structure types through multimodal fusion of one-dimensional acoustic curves, two-dimensional time-frequency diagrams, and porosity physical parameters, thereby enabling the identification of complex pore structures in carbonate rocks.

[0006] The present invention discloses a method for identifying the pore structure of carbonate rocks based on multimodal deep learning. The technical solution includes the following steps: 1. Acquiring acoustic wave data and porosity parameters: The well logging information collection module acquires the one-dimensional acoustic wave waveform data of the sample to be identified received by the downhole array acoustic wave receiving device, and obtains the corresponding porosity parameters. 2. Constructing multimodal input data: Normalize the one-dimensional acoustic waveform data to ensure that the amplitude scale of different samples is consistent, and obtain a standardized one-dimensional acoustic curve; add random noise to the training samples to perform data augmentation to improve the model's noise resistance; perform continuous wavelet transform on the one-dimensional acoustic waveform data to obtain a two-dimensional time-frequency diagram that can characterize the distribution law of acoustic energy on the time and frequency scales. III. Establishment of the Swin-1dCNN Multimodal Fusion Model: The multimodal fusion model includes a one-dimensional convolutional neural network branch, a Swin Transformer branch, and a porosity physical constraint embedding branch. The one-dimensional convolutional neural network branch processes the one-dimensional acoustic wave curve, extracting local temporal features of the acoustic wave, such as first-arrival time offset, first-wave amplitude variation, waveform period variation, and tail wave complexity. The Swin Transformer branch processes the two-dimensional time-frequency diagram, extracting global time-frequency energy features of the acoustic wave's main energy cluster location, low-frequency energy proportion, high-frequency energy retention, and cross-window time-frequency correlation. The porosity physical constraint embedding branch maps porosity parameters to porosity physical features. 4. Multimodal feature fusion: The acoustic wave temporal local features output from the one-dimensional convolutional neural network branch, the time-frequency domain global energy features output from the Swin Transformer branch, and the porosity physical features embedded in the porosity physical constraint branch are concatenated to form a multimodal joint feature vector. The multimodal joint feature vector is nonlinearly fused through a fully connected layer, a batch normalization layer, a ReLU activation function, and a Dropout layer to prevent the model from over-relying on a single modal feature. V. Identifying Pore Structure Types: The input classification branch of the fused multimodal joint feature vector is used to output a three-dimensional probability distribution through the Softmax function. The three-dimensional probabilities correspond to cracks, intergranular pores, and cavities, respectively. The pore type corresponding to the maximum probability is taken as the pore structure identification result of the sample to be identified. VI. Inversion of Pore Aspect Ratio: The input regression branch of the fused multimodal joint feature vector is modified by removing the Softmax function at the end of the classification branch and changing the output layer to a single scalar output, so that the model outputs the predicted value of the pore aspect ratio. During the training process, the mean square error loss function is used to optimize the model's ability to approximate the continuous pore aspect ratio. The inversion results are used to characterize the continuous geometric evolution relationship between cracks, intergranular pores and caverns.

[0007] Preferably, the aforementioned one-dimensional acoustic waveform data is a time-varying sequence of displacement, velocity, pressure, or amplitude data at the receiving point.

[0008] Preferably, in step two, the normalization process is as follows: To ensure consistency in sound wave amplitude scale, sampling length, and time-frequency representation at the model input, the original sound wave sequence is denoted as... x ={ x 1, x 2, ..., x L The standardized one-dimensional acoustic wave curve is obtained using the following normalization expression: , in, This represents the normalized sound wave sequence. Represents the original sound wave sequence. Represents a numerical index in a sound wave sequence. L Indicates the length of the sound wave sequence. μ x This represents the mean of the sound wave sequence. σ x This represents the standard deviation of the sound wave sequence. ε To prevent stable terms with zero denominators, this formula ensures that the amplitude distribution of different samples is on a comparable scale, thereby reducing the impact of differences in sensor gain or logging acquisition conditions on model training.

[0009] Preferably, in step two, the continuous wavelet transform process is as follows: The two-dimensional time-frequency graph is generated by continuous wavelet transform, and its mathematical expression is as follows: , in, For scale parameters, b For time shift parameters, ψ * It is the conjugate function of the mother wavelet function. Indicates the scale of the sound wave signal and time b The time-frequency response at the location; according to After constructing a two-dimensional time-frequency diagram from its energy spectrum, the model can simultaneously utilize the phase abrupt change near the arrival time, the location of the dominant frequency energy cluster, and the tailwave energy attenuation characteristics.

[0010] Preferably, the above-mentioned one-dimensional convolutional neural network branch uses multiple one-dimensional convolutional layers with different kernel sizes. Some of the convolutional kernels are used to reduce the dimension of the input waveform and capture waveform changes within a wide time window, while other convolutional kernels are used to extract first arrival time, phase change and tail wave morphology features. The input to a branch of a one-dimensional convolutional neural network is: , in, L0 represents the number of sampling points. C 0 represents the number of input channels; The following convolution operation slides along the time sampling direction to extract local temporal response features of acoustic wave first arrival, amplitude abrupt change, phase reversal, and wake complexity; Layer j The convolutional kernel performs sliding local weighting on the features of the previous layer to obtain the th convolutional kernel. j The output channel is in the first Feature values ​​at each time position : , in, j The output channel number also corresponds to the [number]. j One convolutional kernel; Output time position number; r This is the index of the internal position of the convolution kernel; q Enter the channel number; For the first Kernel length; This represents the number of input channels for the previous layer. For the first Layer j The convolutional kernel at the _th ... q The first input channel, the first r The weights at each core position; For the first j The bias term for each output channel; σ(.) is the nonlinear activation function; , , The first The stride, fill length, and hole ratio of the convolutional layer; No. The output length of the convolutional layer is: , when When =1, the expression degenerates into the formula for calculating the output length of a normal one-dimensional convolution. By reasonably setting the kernel length, stride, padding and dilation, the output dimension of each convolutional layer is controlled, so that the waveform branch and the subsequent fusion layer maintain dimensionality matching. Subsequently, the max pooling layer is applied to the first... j Each convolutional feature map The output of local downsampling is expressed as follows: , in, k Number the output position after pooling; This is the pooling window length; The pooling step size; For the pooling processj The first channel in the k The output value at each position; max pooling operation can preserve strong response features within the local window and reduce the impact of sampling point position perturbation on the recognition result.

[0011] Preferably, the above-mentioned Swin Transformer branch uses a pre-trained Swin Transformer model as a feature extractor, removes the original classification head, freezes the front computation stage, and activates the two Swin Transformer Blocks in the back; the Swin Transformer Block includes a window multi-head self-attention module, a moving window multi-head self-attention module, a layer normalization module, a multilayer perceptron module, and a residual connection module. The input to the Swin Transformer branch is the feature token matrix Z. (0) ∈R N×d ,in N For the number of patches, d For the token embedding dimension; for the feature token matrix within a local window Z w , No. r The query matrix, key matrix, and value matrix for each attention head are as follows: , Window self-attention is calculated as follows: , Window multi-head self-attention calculation is as follows: , Two consecutive Swing Transformer Blocks are updated sequentially using regular window multi-head self-attention (W-MSA) and moving window multi-head self-attention (SW-MSA): , Among them, Z w The feature token matrix within the window; M The window's side length is [length], and the area inside the window contains [area]. M 2 One patch; r Number the attention head; Q r , K r , V r The first r The query matrix, key matrix, and value matrix of each attention head; W Q,r , W K,r ,W V,r This is the corresponding trainable mapping matrix; T Represents the transpose matrix; d h = d / H The feature dimension of a single attention head; H For the number of heads; B r For the first r The relative position bias matrix within the window corresponding to each attention head; W O To output the mapping matrix; t Number the SwinTransformer Block; Z t-1 For the first t Input features of each block; Z represents the intermediate features after multi-head self-attention via a rule window; t For the process MLP Output characteristics after module; Z represents the intermediate features after multi-head self-attention via a moving window. t+1 The output characteristics of the Swing Transformer Block for moving windows; LN For layer normalization; MLP It is a multilayer perceptron.

[0012] Preferably, the porosity physical constraint embedding branch mentioned above takes a porosity scalar as input, and the porosity scalar is embedded and mapped by a multilayer perceptron (MLP). The MLP includes a linear mapping layer and a ReLU activation layer, which are used to output physical feature vectors.

[0013] Preferably, in step four, the multimodal feature fusion process is as follows: Let the waveform feature vector output by a branch of a one-dimensional convolutional neural network be... The time-frequency eigenvectors output by the SwinTransformer branch are: The physical feature vector output by the porosity physical constraint embedded branch is Then the multimodal joint feature vector is expressed as: , in, That is, the joint feature vector is obtained by splicing 64-dimensional waveform features, 768-dimensional time-frequency plot features and 32-dimensional physical constraint features along the feature dimensions; the joint feature vector can simultaneously characterize the local phase change, amplitude change and tail wave morphology features in a one-dimensional acoustic sequence, the global energy distribution and local time-frequency structure features in a two-dimensional time-frequency plot, and the physical constraint information of porosity for pore structure identification in the same feature space; Subsequently, the 864-dimensional joint feature vector is input into the fusion classifier; the fusion classifier first performs dimensionality reduction and cross-modal feature interaction on the joint features through a linear fully connected layer: , in, W 1∈R 64×864 , b 1∈R 64 , ∈R 64 ; To improve training stability and enhance nonlinear expressive power, the dimensionality-reduced features are sequentially subjected to one-dimensional batch normalization, ReLU activation, and random deactivation: , in, BN (.) indicates a batch normalization operation, used to reduce the impact of differences in feature scales between different modalities on the training process; ReLU (.) represents a nonlinear activation function, used to enhance the nonlinear expressive power of fused features; Dropout (.) represents a random deactivation operation, used to suppress overfitting of local feature nodes, thereby reducing the risk of overfitting and improving the model's generalization ability.

[0014] Preferably, the types of pore structures identified in step five include the following: The fused discriminative features output by the fusion classifier The input classification output layer yields unnormalized classification scores for the three types of pore structures: , in, o =[ o 1, o 2, o [3] represents the classification score vectors for the three types of porous structures. W c and b c These are the weight matrix and bias term of the classification output layer, respectively; the three types of pore structures correspond to cracks, intergranular pores, and caverns, respectively. Subsequently, the unnormalized classification scores are converted into a class probability distribution using the Softmax function: , in, Indicates that the sample belongs to the first Predicted probability of porous structures; p 1. p 2. p 3 corresponds to the classification probabilities of cracks, intergranular pores, and caverns, respectively; The pore structure category of the sample to be identified is determined by the category corresponding to the highest predicted probability: , in, The type of pore structure predicted by the model; The training loss for the classification branch uses the cross-entropy loss function. For a single sample, the classification loss is expressed as: , in, L cls For classification loss; For the uniquely thermally encoded label of the true category of pore structure, when the sample belongs to the first When class, =1, otherwise =0.

[0015] Preferably, the aspect ratio of the inverted pores mentioned in step six includes the following: Let the fusion discriminative features output by the fusion classifier be... Then the predicted pore aspect ratio output by the regression branch is expressed as: , in, To predict the aspect ratio of pores, W r and b r These are the weight matrix and bias term of the regression output layer, respectively; since the porosity is a continuous variable, the Softmax function is not used at the end of the regression branch, but a single linear output node is used to directly output continuous predicted values; During training, the regression branch uses the mean squared error loss function as the optimization objective: , in, L reg For regression loss; B This refers to the batch sample size. AR q For the first q The true aspect ratio of the pores in each sample; For the first q The true aspect ratio of the pores in each sample.

[0016] Compared with the prior art, the beneficial effects of the present invention are as follows: The present invention emits sound waves downhole using an ultrasonic transmitting device. After the sound waves propagate in the formation, they are received by an array of sound wave receiving devices. The one-dimensional sound wave waveform data and porosity parameter data are transmitted to the well logging information acquisition device through a logging cable. Subsequently, the sound wave multimodal data processing module first normalizes the one-dimensional sound wave waveform data and generates two-dimensional time-frequency diagram data through continuous wavelet transform, which is then transmitted to the intelligent recognition system module. The intelligent recognition system module then performs recognition and prediction, and visualizes the pore structure type or pore aspect ratio through the recognition result output module. This invention can simultaneously utilize the local time-domain characteristics of the acoustic curve, the global frequency-domain characteristics of the acoustic time-frequency diagram, and the physical constraint characteristics of porosity to overcome the multiple solutions problem in single acoustic velocity interpretation and single-mode intelligent models. This method can automatically classify and identify fractures, intergranular pores, and caverns, and further output continuous predicted values ​​of pore aspect ratio, providing more refined pore structure parameters for carbonate reservoir evaluation and acoustic logging interpretation. Attached Figure Description

[0017] Figure 1 This is a flowchart of the present invention; Figure 2 This is a schematic diagram of a carbonate rock pore structure identification system based on multimodal deep learning; Figure 3 This is a schematic diagram of the Swin-1dCNN multimodal fusion model structure; Figure 4 This is a schematic diagram of the pore structure recognition results of the Swin-1dCNN multimodal fusion model in the test set; Figure 5 This is a schematic diagram of the pore structure identification results of a one-dimensional convolutional neural network model in a blind test set; Figure 6 This is a schematic diagram of the pore structure identification results of the Swin Transformer model in the blind test set; Figure 7 This is a schematic diagram of the pore structure recognition results of the Swin-1dCNN multimodal fusion model in a blind test set; Figure 8 This is a schematic diagram of the aperture aspect ratio prediction results of the Swin-1dCNN multimodal fusion model in the test set; Figure 9 This is a schematic diagram of the aperture aspect ratio prediction results of a one-dimensional convolutional neural network model in a blind test set; Figure 10 This is a schematic diagram of the porosity prediction results of the Swin Transformer model in the blind test set; Figure 11This is a schematic diagram of the aperture aspect ratio prediction results of the Swin-1dCNN multimodal fusion model in the blind test set; In the diagram: 1-Drilling platform; 2-Drill pipe; 3-Casing; 4-Open hole section; 5-Array acoustic wave receiver; 6-Ultrasonic wave transmitter; 7-Logging cable; 8-Logging information acquisition device; 9-Acoustic multimodal data processing module; 10-Intelligent recognition system module; 11-Recognition result output module. Detailed Implementation

[0018] The preferred embodiments of the present invention will be described below with reference to the accompanying drawings. It should be understood that the preferred embodiments described herein are for illustration and explanation only and are not intended to limit the present invention.

[0019] Example 1: This invention provides a carbonate rock pore structure identification system based on multimodal deep learning. The system includes a well logging information collection module, an acoustic multimodal data processing module 9, an intelligent identification system module 10, and an identification result output module 11. The well logging information collection module includes an array acoustic receiver 5, an ultrasonic transmitter 6, a well logging cable 7, and a well logging information acquisition device 8. On the drilling platform 1, the drill bit is lowered into the open hole section 4 below the casing 3 via the drill pipe 2. The ultrasonic transmitter 6 is located in the lower middle part of the drill pipe 2, and one or more array acoustic receivers 5 are located in the upper middle part of the drill pipe 2. The array acoustic receivers 5 are connected to the well logging information acquisition device 8 installed on the ground via the well logging cable 7. The input end of the acoustic multimodal data processing module 9 is connected to the well logging information acquisition device 8 via a wire, and the output end is connected to the intelligent identification system module 10. The output end of the intelligent identification system module 10 is connected to the identification result output module 11. The well logging information collection module is used to acquire one-dimensional acoustic waveform data and porosity parameters of the carbonate rock formation to be identified; the acoustic multimodal data processing module 9 is used to normalize and perform continuous wavelet transform on the one-dimensional acoustic waveform data to form one-dimensional acoustic curve data and two-dimensional time-frequency diagram data; the intelligent identification system module 10 is used to extract the local waveform features of the one-dimensional acoustic curve, the global energy features and porosity physical features of the two-dimensional time-frequency diagram, and fuse the three types of features for identification; the identification result output module 11 is used to output the three types of pore structure types, namely fractures, intergranular pores and caverns, or to output the predicted value of the pore aspect ratio.

[0020] This invention uses an ultrasonic transmitter 6 to emit sound waves downhole. After the sound waves propagate in the formation, they are received by an array of acoustic wave receivers 5. The one-dimensional acoustic waveform data and porosity parameter data are transmitted to the logging information acquisition device 8 via a logging cable 7. Subsequently, the acoustic multimodal data processing module 9 first normalizes the one-dimensional acoustic waveform data and generates two-dimensional time-frequency graph data through continuous wavelet transform, which is then transmitted to the intelligent identification system module 10. The intelligent identification system module 10 performs identification and prediction, and the identification result output module 11 visualizes the pore structure type or pore aspect ratio.

[0021] The present invention discloses a method for identifying the pore structure of carbonate rocks based on multimodal deep learning, comprising the following steps: 1. Acquiring acoustic wave data and porosity parameters: The well logging information collection module acquires the one-dimensional acoustic wave waveform data of the sample to be identified received by the downhole array acoustic wave receiver 5, and obtains the corresponding porosity parameters. The aforementioned one-dimensional acoustic waveform data is a time-varying sequence of displacement, velocity, pressure, or amplitude data at the receiving point.

[0022] 2. Constructing multimodal input data: Normalize the one-dimensional acoustic waveform data to ensure that the amplitude scale of different samples is consistent, and obtain a standardized one-dimensional acoustic curve; add random noise to the training samples to perform data augmentation to improve the model's noise resistance; perform continuous wavelet transform on the one-dimensional acoustic waveform data to obtain a two-dimensional time-frequency diagram that can characterize the distribution law of acoustic energy on the time and frequency scales. The normalization process is as follows: To ensure consistency in sound wave amplitude scale, sampling length, and time-frequency representation at the model input, the original sound wave sequence is denoted as... x ={ x 1, x 2, ..., x L The standardized one-dimensional acoustic wave curve is obtained using the following normalization expression: , in, This represents the normalized sound wave sequence. Represents the original sound wave sequence. Represents a numerical index in a sound wave sequence. L Indicates the length of the sound wave sequence. μ x This represents the mean of the sound wave sequence. σ x This represents the standard deviation of the sound wave sequence. εTo prevent stable terms with zero denominators, this formula ensures that the amplitude distribution of different samples is on a comparable scale, thereby reducing the impact of differences in sensor gain or logging acquisition conditions on model training.

[0023] Preferably, in step two, the continuous wavelet transform process is as follows: The two-dimensional time-frequency graph is generated by continuous wavelet transform, and its mathematical expression is as follows: , in, For scale parameters, b For time shift parameters, ψ * It is the conjugate function of the mother wavelet function. Indicates the scale of the sound wave signal and time b The time-frequency response at the location; according to After constructing a two-dimensional time-frequency diagram from its energy spectrum, the model can simultaneously utilize the phase abrupt change near the arrival time, the location of the dominant frequency energy cluster, and the tailwave energy attenuation characteristics.

[0024] III. Establishment of the Swin-1dCNN Multimodal Fusion Model: The multimodal fusion model includes a one-dimensional convolutional neural network branch, a Swin Transformer branch, and a porosity physical constraint embedding branch. The one-dimensional convolutional neural network branch processes the one-dimensional acoustic wave curve, extracting local temporal features of the acoustic wave, such as first-arrival time offset, first-wave amplitude variation, waveform period variation, and tail wave complexity. The Swin Transformer branch processes the two-dimensional time-frequency diagram, extracting global time-frequency energy features of the acoustic wave's main energy cluster location, low-frequency energy proportion, high-frequency energy retention, and cross-window time-frequency correlation. The porosity physical constraint embedding branch maps porosity parameters to porosity physical features. Among them, the above-mentioned one-dimensional convolutional neural network branch uses multiple one-dimensional convolutional layers with different kernel sizes. A larger portion of the convolutional kernels is used to reduce the dimension of the input waveform and capture waveform changes within a wide time window, while a smaller portion of the convolutional kernels is used to extract first arrival time, phase transition, and tail wave morphology features. The input to a branch of a one-dimensional convolutional neural network is: , in, L 0 represents the number of sampling points. C 0 represents the number of input channels; The following convolution operation slides along the time sampling direction to extract local temporal response features of acoustic wave first arrival, amplitude abrupt change, phase reversal, and wake complexity; Layer j The convolutional kernel performs sliding local weighting on the features of the previous layer to obtain the th convolutional kernel. jThe output channel is in the first Feature values ​​at each time position : , in, j The output channel number also corresponds to the [number]. j One convolutional kernel; Output time position number; r This is the index of the internal position of the convolution kernel; q Enter the channel number; For the first Kernel length; This represents the number of input channels for the previous layer. For the first Layer j The convolutional kernel at the _th ... q The first input channel, the first r The weights at each core position; For the first j The bias term for each output channel; σ(.) is the nonlinear activation function; , , The first The stride, fill length, and hole ratio of the convolutional layer; No. The output length of the convolutional layer is: , when When =1, the expression degenerates into the formula for calculating the output length of a normal one-dimensional convolution. By reasonably setting the kernel length, stride, padding and dilation, the output dimension of each convolutional layer is controlled, so that the waveform branch and the subsequent fusion layer maintain dimensionality matching. Subsequently, the max pooling layer is applied to the first... j Each convolutional feature map The output of local downsampling is expressed as follows: , in, k Number the output position after pooling; This is the pooling window length; The pooling step size; For the pooling process j The first channel in the k The output value at each position; max pooling operation can preserve strong response features within the local window and reduce the impact of sampling point position perturbation on the recognition result.

[0025] Preferably, the above-mentioned Swin Transformer branch uses a pre-trained Swin Transformer model as a feature extractor, removes the original classification head, freezes the front computation stage, and activates the two Swin Transformer Blocks in the back; the Swin Transformer Block includes a window multi-head self-attention module, a moving window multi-head self-attention module, a layer normalization module, a multilayer perceptron module, and a residual connection module. The input to the Swin Transformer branch is the feature token matrix Z. (0) ∈R N×d ,in N For the number of patches, d For the token embedding dimension; for the feature token matrix within a local window Z w , No. r The query matrix, key matrix, and value matrix for each attention head are as follows: , Window self-attention is calculated as follows: , Window multi-head self-attention calculation is as follows: , Two consecutive Swing Transformer Blocks are updated sequentially using regular window multi-head self-attention (W-MSA) and moving window multi-head self-attention (SW-MSA): , Among them, Z w The feature token matrix within the window; M The window's side length is [length], and the area inside the window contains [area]. M 2 One patch; r Number the attention head; Q r , K r , V r The first r The query matrix, key matrix, and value matrix of each attention head; W Q,r , W K,r , W V,r This is the corresponding trainable mapping matrix; T Represents the transpose matrix; d h = d / H The feature dimension of a single attention head; H For the number of heads; B r For the first r The relative position bias matrix within the window corresponding to each attention head; W O To output the mapping matrix; t Number the SwinTransformer Block; Z t-1 For the first t Input features of each block; Z represents the intermediate features after multi-head self-attention via a rule window; t For the process MLP Output characteristics after module; Z represents the intermediate features after multi-head self-attention via a moving window. t+1 The output characteristics of the Swing Transformer Block for moving windows; LN For layer normalization; MLP It is a multilayer perceptron.

[0026] Preferably, the porosity physical constraint embedding branch mentioned above takes a porosity scalar as input, and the porosity scalar is embedded and mapped by a multilayer perceptron (MLP). The MLP includes a linear mapping layer and a ReLU activation layer, which are used to output physical feature vectors.

[0027] 4. Multimodal feature fusion: The acoustic wave temporal local features output from the one-dimensional convolutional neural network branch, the time-frequency domain global energy features output from the Swin Transformer branch, and the porosity physical features embedded in the porosity physical constraint branch are concatenated to form a multimodal joint feature vector. The multimodal joint feature vector is nonlinearly fused through a fully connected layer, a batch normalization layer, a ReLU activation function, and a Dropout layer to prevent the model from over-relying on a single modal feature. Preferably, in step four, the multimodal feature fusion process is as follows: Let the waveform feature vector output by a branch of a one-dimensional convolutional neural network be... The time-frequency eigenvectors output by the SwinTransformer branch are: The physical feature vector output by the porosity physical constraint embedded branch is Then the multimodal joint feature vector is expressed as: , in, That is, the joint feature vector is obtained by splicing 64-dimensional waveform features, 768-dimensional time-frequency plot features and 32-dimensional physical constraint features along the feature dimensions; the joint feature vector can simultaneously characterize the local phase change, amplitude change and wake morphology features in a one-dimensional acoustic wave sequence, the global energy distribution and local time-frequency structure features in a two-dimensional time-frequency plot, and the physical constraint information of porosity for pore structure identification in the same feature space.

[0028] Subsequently, the 864-dimensional joint feature vector is input into the fusion classifier; the fusion classifier first performs dimensionality reduction and cross-modal feature interaction on the joint features through a linear fully connected layer: , in, W 1∈R 64×864 , b 1∈R 64 , ∈R 64 ; To improve training stability and enhance nonlinear expressive power, the dimensionality-reduced features are sequentially subjected to one-dimensional batch normalization, ReLU activation, and random deactivation: , in, BN (.) indicates a batch normalization operation, used to reduce the impact of differences in feature scales between different modalities on the training process; ReLU (.) represents a nonlinear activation function, used to enhance the nonlinear expressive power of fused features; Dropout (.) represents a random deactivation operation, used to suppress overfitting of local feature nodes, thereby reducing the risk of overfitting and improving the model's generalization ability.

[0029] V. Identifying Pore Structure Types: The input classification branch of the fused multimodal joint feature vector is used to output a three-dimensional probability distribution through the Softmax function. The three-dimensional probabilities correspond to cracks, intergranular pores, and cavities, respectively. The pore type corresponding to the maximum probability is taken as the pore structure identification result of the sample to be identified. Preferably, the types of pore structures identified in step five include the following: The fused discriminative features output by the fusion classifier The input classification output layer yields unnormalized classification scores for the three types of pore structures: , in, o =[ o 1, o 2, o [3] represents the classification score vectors for the three types of porous structures. W c and bc These are the weight matrix and bias term of the classification output layer, respectively; the three types of pore structures correspond to cracks, intergranular pores, and caverns, respectively. Subsequently, the unnormalized classification scores are converted into a class probability distribution using the Softmax function: , in, Indicates that the sample belongs to the first Predicted probability of porous structures; p 1. p 2. p 3 corresponds to the classification probabilities of cracks, intergranular pores, and caverns, respectively; The pore structure category of the sample to be identified is determined by the category corresponding to the highest predicted probability: , in, The type of pore structure predicted by the model; The training loss for the classification branch uses the cross-entropy loss function. For a single sample, the classification loss is expressed as: , in, L cls For classification loss; For the uniquely thermally encoded label of the true category of pore structure, when the sample belongs to the first When class, =1, otherwise =0.

[0030] VI. Inversion of Pore Aspect Ratio: The input regression branch of the fused multimodal joint feature vector is modified by removing the Softmax function at the end of the classification branch and changing the output layer to a single scalar output, so that the model outputs the predicted value of the pore aspect ratio. During the training process, the mean square error loss function is used to optimize the model's ability to approximate the continuous pore aspect ratio. The inversion results are used to characterize the continuous geometric evolution relationship between cracks, intergranular pores and caverns.

[0031] Specifically, for the porosity inversion task, a regression branch is set after the multimodal fusion module. This regression branch shares the front-end feature extraction module and the multimodal feature fusion structure with the porosity structure classification branch. That is, the one-dimensional convolutional neural network branch, the Swing Transformer branch, the porosity physical constraint embedding branch, the branch-level random deactivation, the 864-dimensional feature splicing, and the fusion dimensionality reduction operation remain unchanged. Only after the fusion of discriminative features is an independent single scalar regression output layer set to output the predicted porosity. The aspect ratio of the inverted pores mentioned in step six includes the following: Let the fusion discriminative features output by the fusion classifier be... Then the predicted pore aspect ratio output by the regression branch is expressed as: , in, To predict the aspect ratio of pores, W r and b r These are the weight matrix and bias term of the regression output layer, respectively; since the porosity is a continuous variable, the Softmax function is not used at the end of the regression branch, but a single linear output node is used to directly output continuous predicted values; During training, the regression branch uses the mean squared error loss function as the optimization objective: , in, L reg For regression loss; B This refers to the batch sample size. AR q For the first q The true aspect ratio of the pores in each sample; For the first q The true aspect ratio of the pores in each sample.

[0032] The system and method provided by this invention can simultaneously utilize the local time-domain characteristics of acoustic curves, the global frequency-domain characteristics of acoustic time-frequency maps, and the physical constraints of porosity to overcome the multiple-solution problem in single acoustic velocity interpretation and single-mode intelligent models. This method can automatically classify and identify fractures, intergranular pores, and caverns, and further output continuous predicted values ​​of porosity aspect ratio, providing more refined pore structure parameters for carbonate reservoir evaluation and acoustic logging interpretation.

[0033] Figure 4 This is a schematic diagram of the pore structure recognition results of the Swin-1dCNN multimodal fusion model on the test set. After training, the Swin-1dCNN model achieved a classification accuracy of 98.2% on the test set for the identification of three types of pore structures. For easily confused intergranular pore samples, the Swin-1dCNN model achieved an intergranular pore identification precision of 97.0% and an intergranular pore identification recall of 97.3%. For crack samples, both the identification precision and recall reached 100.0%. These results indicate that the one-dimensional acoustic curve branch can compensate for the smoothing effect of the two-dimensional time-frequency map on the first arrival details of the first wave, and the two-dimensional time-frequency map branch can compensate for the lack of frequency dimension in the one-dimensional acoustic curve, thereby reducing misclassification among cracks, intergranular pores, and cavities.

[0034] Figure 5 This is a schematic diagram of the pore structure identification results of a one-dimensional convolutional neural network model in a blind test set; Figure 6 This is a schematic diagram of the pore structure identification results of the Swin Transformer model in the blind test set; Figure 7 This is a schematic diagram of the pore structure recognition results of the Swin-1dCNN multimodal fusion model in a blind test set; In an independent blind test experiment, blind test samples were constructed using porosity and pore aspect ratio conditions not found in the training and test sets. For example, the porosity was set to 12%, and pore aspect ratios of 0.22 and 0.25, which were not found in the training set, were used to verify the model's generalization ability when faced with unseen pore geometric distributions. The Swin-1dCNN model achieved a recognition accuracy of 94.7% on this blind test set. Compared with single-modal models, in this blind test implementation, the recognition accuracy using only a one-dimensional convolutional neural network (1dCNN) model was 88.9%, and the recognition accuracy using only a Swin Transformer model was 92.1%. The Swin-1dCNN fusion model improved the accuracy by 5.8% and 2.6%, respectively, indicating that the multimodal fusion structure can reduce the model's dependence on a single data distribution.

[0035] Figure 8 This is a schematic diagram showing the aperture aspect ratio prediction results of the Swin-1dCNN multimodal fusion model in the test set. In one implementation of pore aspect ratio inversion, the terminal classification output layer of the aforementioned Swin-1dCNN model is replaced with a single scalar regression output layer, and the model is trained using the mean squared error loss function to output the predicted pore aspect ratio. Test results show that the coefficient of determination R0 of the Swin-1dCNN model's predictions is [value missing]. 2 It is 0.941.

[0036] Figure 9 This is a schematic diagram of the aperture aspect ratio prediction results of a one-dimensional convolutional neural network model in a blind test set; Figure 10 This is a schematic diagram of the porosity prediction results of the Swin Transformer model in the blind test set; Figure 11 This is a schematic diagram of the aperture aspect ratio prediction results of the Swin-1dCNN multimodal fusion model in the blind test set; In one blind test inversion implementation, the Swin-1dCNN fusion model's pore aspect ratio prediction determination coefficient R in the independent blind test set... 2The value remained around 0.88, still superior to the single-mode baseline model. This result indicates that the model can not only output three discrete pore structure types—fractures, intergranular pores, and caverns—but also quantitatively invert the aspect ratio parameters of continuous pores, thus providing support for the implementation of pore geometry characterization and rock physical parameter prediction in carbonate reservoir evaluation.

[0037] The above-described implementation results demonstrate that the present invention can improve classification accuracy and blind test generalization ability in pore structure identification tasks, and improve the continuous parameter fitting accuracy in pore aspect ratio inversion tasks. It can effectively reduce the ambiguity of single acoustic velocity interpretation and single-mode intelligent models in the identification of complex carbonate rock pore structures.

[0038] In practical applications, after the acoustic logging data to be interpreted is input into the system, it undergoes preprocessing to obtain a one-dimensional acoustic curve and a two-dimensional time-frequency diagram, along with porosity parameters. The system automatically outputs the predicted pore structure type and pore aspect ratio. These results can be used to further determine the development characteristics of fractures, intergranular pores, and caverns in the reservoir, providing a basis for carbonate reservoir evaluation, effective reservoir space identification, and rock physical modeling.

[0039] This invention reduces the ambiguity in single-sonic velocity interpretation and improves the accuracy of intelligent identification of carbonate rock pore structure and quantitative inversion of pore aspect ratio. This invention enables automatic classification and identification of complex pore structures in carbonate rocks and inversion of pore aspect ratio based on acoustic signals, providing technical support for carbonate reservoir evaluation, acoustic logging interpretation, and rock physical parameter prediction.

[0040] The above description is merely a partial preferred embodiment of the present invention. Any person skilled in the art can modify the above-described technical solutions or modify them into equivalent technical solutions. Therefore, any simple modifications or equivalent transformations made based on the technical solutions of the present invention fall within the scope of protection claimed by the present invention.

Claims

1. A method for identifying the pore structure of carbonate rocks based on multimodal deep learning, characterized by: Includes the following steps:

1. Acquiring acoustic wave data and porosity parameters: The one-dimensional acoustic wave waveform data of the sample to be identified received by the array acoustic wave receiving device (5) downhole is obtained through the well logging information collection module, and the corresponding porosity parameters are obtained.

2. Constructing multimodal input data: Normalize the one-dimensional acoustic waveform data to ensure that the amplitude scale of different samples is consistent, and obtain a standardized one-dimensional acoustic curve; Add random noise to the training samples to perform data augmentation to improve the model's noise resistance. A continuous wavelet transform is performed on one-dimensional acoustic waveform data to obtain a two-dimensional time-frequency diagram that can characterize the distribution of acoustic energy on the time and frequency scales. III. Establishment of the Swin-1dCNN Multimodal Fusion Model: The multimodal fusion model includes a one-dimensional convolutional neural network branch, a Swin Transformer branch, and a porosity physical constraint embedding branch. The one-dimensional convolutional neural network branch processes the one-dimensional acoustic wave curve, extracting local temporal features of the acoustic wave, such as first-arrival time offset, first-wave amplitude variation, waveform period variation, and tail wave complexity. The Swin Transformer branch processes the two-dimensional time-frequency diagram, extracting global time-frequency energy features of the acoustic wave's main energy cluster location, low-frequency energy proportion, high-frequency energy retention, and cross-window time-frequency correlation. The porosity physical constraint embedding branch maps porosity parameters to porosity physical features. The Swing Transformer branch uses a pre-trained Swing Transformer model as a feature extractor, removes the original classification head, freezes the front computation stage, and activates the two Swing Transformer Blocks in the back. The Swing Transformer Block includes a window multi-head self-attention module, a moving window multi-head self-attention module, a layer normalization module, a multilayer perceptron module, and a residual connection module. The input to the Swin Transformer branch is the feature token matrix Z. (0) ∈R N×d ,in N For the number of patches, d For the token embedding dimension; for the feature token matrix within the local window Z w , No. r The query matrix, key matrix, and value matrix for each attention head are as follows: ， Window self-attention is calculated as follows: ， Window multi-head self-attention calculation is as follows: ， Two consecutive Swing Transformer Blocks are updated sequentially using regular window multi-head self-attention (W-MSA) and moving window multi-head self-attention (SW-MSA): ， Among them, Z w The feature token matrix within the window; M The window's side length is [length], and the area inside the window contains [areas]. M 2 One patch; r Number the attention head; Q r , K r , V r The first r The query matrix, key matrix, and value matrix of each attention head; W Q,r , W K,r , W V,r This is the corresponding trainable mapping matrix; T Represents the transpose matrix; d h = d / H The feature dimension of a single attention head; H For the number of heads; B r For the first r The relative position bias matrix within the window corresponding to each attention head; W O To output the mapping matrix; t Number the SwinTransformer Block; Z t-1 For the first t Input features of each block; Z represents the intermediate features after multi-head self-attention via a rule window; t For the process MLP Output characteristics after module; Z represents the intermediate features after multi-head self-attention via a moving window. t+1 The output characteristics of the Swing Transformer Block for moving windows; LN For layer normalization; MLP It is a multilayer perceptron; 4. Multimodal feature fusion: The acoustic wave temporal local features output from the one-dimensional convolutional neural network branch, the time-frequency domain global energy features output from the SwinTransformer branch, and the porosity physical features embedded in the porosity physical constraint branch are concatenated to form a multimodal joint feature vector. The multimodal joint feature vector is nonlinearly fused through a fully connected layer, a batch normalization layer, a ReLU activation function, and a Dropout layer to prevent the model from over-relying on a single modal feature. V. Identifying Pore Structure Types: The input classification branch of the fused multimodal joint feature vector is used to output a three-dimensional probability distribution through the Softmax function. The three-dimensional probabilities correspond to cracks, intergranular pores, and cavities, respectively. The pore type corresponding to the maximum probability is taken as the pore structure identification result of the sample to be identified. VI. Inversion of Pore Aspect Ratio: The input regression branch of the fused multimodal joint feature vector is modified by removing the Softmax function at the end of the classification branch and changing the output layer to a single scalar output, so that the model outputs the predicted value of the pore aspect ratio. During the training process, the mean square error loss function is used to optimize the model's ability to approximate the continuous pore aspect ratio. The inversion results are used to characterize the continuous geometric evolution relationship between cracks, intergranular pores and caverns.

2. The method for identifying the pore structure of carbonate rocks based on multimodal deep learning according to claim 1, characterized in that: The one-dimensional acoustic waveform data is a time-varying sequence of displacement, velocity, pressure, or amplitude data at the receiving point.

3. The method for identifying the pore structure of carbonate rocks based on multimodal deep learning according to claim 2, characterized in that: In step two, the normalization process is as follows: To ensure consistency in the acoustic wave amplitude scale, sampling length, and time-frequency representation at the model input, the original acoustic wave sequence is denoted as... x ={ x 1, x 2, ..., x L The standardized one-dimensional acoustic wave curve is obtained using the following normalization expression: ， in, This represents the normalized sound wave sequence. Represents the original sound wave sequence. Represents a numerical index in a sound wave sequence. L Indicates the length of the sound wave sequence. μ x This represents the mean of the sound wave sequence. σ x This represents the standard deviation of the sound wave sequence. ε To prevent stable terms with a denominator of zero, this formula ensures that the amplitude distribution of different samples is on a comparable scale, thereby reducing the impact of differences in sensor gain or logging acquisition conditions on model training.

4. The method for identifying the pore structure of carbonate rocks based on multimodal deep learning according to claim 3, characterized in that: in step two, the continuous wavelet transform process is as follows: The two-dimensional time-frequency graph is generated by continuous wavelet transform, and its mathematical expression is as follows: ， in, For scale parameters, b For time shift parameters, ψ * It is the conjugate function of the mother wavelet function. Indicates the scale of the sound wave signal and time b The time-frequency response at the location; according to After constructing a two-dimensional time-frequency diagram from its energy spectrum, the model can simultaneously utilize the phase abrupt change near the arrival time, the location of the dominant frequency energy cluster, and the tailwave energy attenuation characteristics.

5. The method for identifying the pore structure of carbonate rocks based on multimodal deep learning according to claim 4, characterized in that: The one-dimensional convolutional neural network branch uses multiple one-dimensional convolutional layers with different kernel sizes. Some kernels are used to reduce the dimension of the input waveform and capture waveform changes within a wide time window, while others are used to extract first arrival time, phase transition, and tail wave morphology features. The input to a branch of a one-dimensional convolutional neural network is: ， in, L 0 represents the number of sampling points. C 0 represents the number of input channels; The following convolution operation slides along the time sampling direction to extract local temporal response features of acoustic wave first arrival, amplitude abrupt change, phase reversal, and wake complexity; Layer j The convolutional kernel performs sliding local weighting on the features of the previous layer to obtain the th convolutional kernel. j The output channel is in the first Feature values ​​at each time position : ， in, j The output channel number also corresponds to the [number]. j One convolutional kernel; Output the time position number; r This is the index of the internal position of the convolution kernel; q Enter the channel number; For the first Layer convolution kernel length; This represents the number of input channels for the previous layer. For the first Layer j The convolutional kernel at the _th ... q The first input channel, the first r The weights at each core position; For the first j The bias term for each output channel; σ(.) is the nonlinear activation function; , , The first The stride, fill length, and hole ratio of the convolutional layer; No. The output length of the convolutional layer is: ， when When =1, the expression degenerates into the formula for calculating the output length of a normal one-dimensional convolution. By reasonably setting the kernel length, stride, padding and dilation, the output dimension of each convolutional layer is controlled, so that the waveform branch and the subsequent fusion layer maintain dimensionality matching. Subsequently, the max pooling layer is applied to the first... j Each convolutional feature map The output of local downsampling is expressed as follows: ， in, k Number the output position after pooling; This is the pooling window length; The pooling step size; For the pooling process j The first channel in the k The output value at each position; max pooling operation can preserve strong response features within the local window and reduce the impact of sampling point position perturbation on the recognition result.

6. The method for identifying the pore structure of carbonate rocks based on multimodal deep learning according to claim 5, characterized in that: The porosity physical constraint embedding branch takes a porosity scalar as input, and the porosity scalar is embedded and mapped by a multilayer perceptron (MLP). The MLP includes a linear mapping layer and a ReLU activation layer, which are used to output physical feature vectors.

7. The method for identifying the pore structure of carbonate rocks based on multimodal deep learning according to claim 6, characterized in that: In step four, the multimodal feature fusion process is as follows: Let the waveform feature vector output by a branch of a one-dimensional convolutional neural network be... The time-frequency eigenvectors output by the Swing Transformer branch are: The physical feature vector output by the porosity physical constraint embedded branch is Then the multimodal joint feature vector is expressed as: ， in, That is, the joint feature vector is obtained by splicing 64-dimensional waveform features, 768-dimensional time-frequency plot features and 32-dimensional physical constraint features along the feature dimensions; the joint feature vector can simultaneously characterize the local phase change, amplitude change and tail wave morphology features in a one-dimensional acoustic sequence, the global energy distribution and local time-frequency structure features in a two-dimensional time-frequency plot, and the physical constraint information of porosity for pore structure identification in the same feature space; Subsequently, the 864-dimensional joint feature vector is input into the fusion classifier; the fusion classifier first performs dimensionality reduction and cross-modal feature interaction on the joint features through a linear fully connected layer: ， in, W 1∈R 64×864 , b 1∈R 64 , ∈R 64 ; To improve training stability and enhance nonlinear expressive power, the dimensionality-reduced features are sequentially subjected to one-dimensional batch normalization, ReLU activation, and random deactivation: ， in, BN (.) indicates a batch normalization operation, used to reduce the impact of differences in feature scales between different modalities on the training process; ReLU (.) represents a nonlinear activation function, used to enhance the nonlinear expressive power of fused features; Dropout (.) represents a random deactivation operation, used to suppress overfitting of local feature nodes, thereby reducing the risk of overfitting and improving the model's generalization ability.

8. The method for identifying the pore structure of carbonate rocks based on multimodal deep learning according to claim 7, characterized in that: The types of pore structures mentioned in step five include the following: The fused discriminative features output by the fusion classifier The input classification output layer yields unnormalized classification scores for the three types of pore structures: ， in, o =[ o 1, o 2, o [3] represents the classification score vectors for the three types of porous structures. W c and b c These are the weight matrix and bias term of the classification output layer, respectively; the three types of pore structures correspond to cracks, intergranular pores, and caverns, respectively. Subsequently, the unnormalized classification scores are converted into a class probability distribution using the Softmax function: ， in, Indicates that the sample belongs to the first Predicted probability of porous structures; p 1. p 2. p 3 corresponds to the classification probabilities of cracks, intergranular pores, and caverns, respectively; The pore structure category of the sample to be identified is determined by the category corresponding to the highest predicted probability: ， in, The type of pore structure predicted by the model; The training loss for the classification branch uses the cross-entropy loss function. For a single sample, the classification loss is expressed as: ， in, L cls For classification loss; For the uniquely thermally encoded label of the true category of pore structure, when the sample belongs to the first When class, =1, otherwise =0.

9. The method for identifying the pore structure of carbonate rocks based on multimodal deep learning according to claim 8, characterized in that: The aspect ratio of the inverted pores mentioned in step six includes the following: Let the fusion discriminative features output by the fusion classifier be... Then the predicted pore aspect ratio output by the regression branch is expressed as: ， in, To predict the aspect ratio of pores, W r and b r These are the weight matrix and bias term of the regression output layer, respectively; since the porosity is a continuous variable, the Softmax function is not used at the end of the regression branch, but a single linear output node is used to directly output continuous predicted values; During training, the regression branch uses the mean squared error loss function as the optimization objective: ， in, L reg For regression loss; B This refers to the batch sample size. AR q For the first q The true aspect ratio of the pores in each sample; For the first q The true aspect ratio of the pores in each sample.

Images