Frequency Domain Divide and Conquer and Physical Constraints: Multi-Field Intelligent Prediction Method and System for CO2 Oil Flooding

Abstract

This invention belongs to the field of CO2 enhanced oil recovery technology and discloses a multi-physics intelligent prediction method and system for CO2 enhanced oil recovery based on frequency domain divide-and-conquer and physical constraints. The method includes: extracting three-dimensional time-series data samples of multiple physics fields from a constructed numerical simulation model of a CO2 enhanced oil recovery reservoir; obtaining the original time-series compressed feature vector of each physics field from the data samples through a time-series processing branch; mapping the time-series compressed feature vector to the frequency domain space; performing high-low frequency division and dynamic fusion of the features of different physics fields to obtain global frequency-domain fusion features, which are then inverted into time-domain dynamic incremental features. Adaptive weighting is then performed to obtain the fusion features of each physics field. A prediction model is constructed and trained, ensuring the physical consistency of the predicted values ​​during the training process; parameters are input into the trained prediction model, and the CO2 enhanced oil recovery multi-physics field prediction results are output. This invention can simultaneously meet the engineering requirements of high precision, high efficiency, strong robustness, physical compliance, and multi-physics field collaboration.

Description

Technical Field

[0001] This invention belongs to the field of CO2 flooding technology, and specifically relates to a multi-field intelligent prediction method and system for CO2 flooding based on frequency domain divide-and-conquer and physical constraints. Background Technology

[0002] The global climate problem caused by excessive industrial CO2 emissions is becoming increasingly severe, and carbon capture, utilization and storage (CCUS) technology is one of the key technologies for achieving carbon neutrality. CO2 enhanced oil recovery, as an important application of CCUS technology, injects CO2 into oil and gas reservoirs, utilizes the displacement effect of CO2 to improve oil and gas recovery, and simultaneously achieves geological CO2 sequestration, thus possessing the dual value of emission reduction and resource development.

[0003] After CO2 is injected into an oil reservoir, it undergoes complex multiphase flow in a multi-scale pore-fracture system. The spatiotemporal dynamic changes of its oil displacement effect exhibit significant characteristics of strong nonlinearity, multi-field coupling, and multi-scale correlation. Accurately predicting the spatiotemporal dynamic evolution of CO2 oil displacement effect is a core prerequisite for evaluating oil displacement efficiency, optimizing injection parameters, and supporting real-time decision-making and scheme optimization in CO2 oil displacement projects.

[0004] Current methods for dynamically predicting the effects of CO2 flooding mainly rely on traditional reservoir numerical simulation. This method faces severe grid dependency and computational efficiency bottlenecks, with a single complete simulation taking several hours to several days, which cannot meet the rapid computational needs for optimizing flooding schemes and real-time dynamic control.

[0005] In recent years, deep learning-based surrogate model technology has begun to be applied in oil and gas reservoir engineering and CO2 enhanced oil recovery, which can significantly improve computational efficiency. However, existing technical solutions have obvious drawbacks: 1. The use of a single feature extraction module without constructing independent links for the distribution characteristics and evolution laws of different physical fields results in the submergence of the differentiated features of different physical fields and insufficient feature capture capability for key areas such as CO2 transport front and permeability mutation zone.

[0006] 2. Feature fusion is mostly completed in the time and spatial domains. No frequency domain feature representation and processing method is constructed for the nonlinear evolution process of CO2 oil displacement. It is impossible to transform the strong nonlinear problem in the time domain into a linear problem in the frequency domain for solution. Long-term prediction is prone to error accumulation and numerical oscillation.

[0007] 3. The lack of separate processing and adaptive fusion of high and low frequency features makes it impossible to simultaneously address the differentiated extraction and fusion requirements of displacement front mutation features and reservoir background stability features, resulting in insufficient prediction robustness.

[0008] 4. Physical constraints are only introduced in the loss function without constructing a physical residual fusion link at the feature level, and no collaborative prediction mechanism is built for the coupling mechanism between multiple physical fields. The prediction results are prone to non-physical interpretations that violate the laws of seepage mechanics, and cannot guarantee the physical consistency between multiple fields such as pressure, saturation, and concentration, thus limiting field applications.

[0009] The following are representative solutions of existing related technologies: Existing technology 1: The invention patent with publication number CN118917355A discloses a method for constructing a recurrent neural network surrogate model for predicting CO2 sequestration in deep saline aquifers. Based on the multiphase seepage numerical simulation of CO2 sequestration in deep saline aquifers, a recurrent neural network surrogate model is established using RU-Net network to simulate CO2 saturation and reservoir temperature during the CO2 geological sequestration process, thereby realizing the regional sequestration potential assessment and analysis of CO2 sequestration in deep saline aquifers.

[0010] The existing technology has the following problems: it does not consider the multi-scale heterogeneity and multiphase flow displacement mechanism of oil and gas reservoirs, it can only capture time-series characteristics, its spatial characterization ability is insufficient, and it only supports single-field prediction, which cannot meet the requirements of simultaneous prediction of the entire CO2 oil recovery process and multiple physical fields.

[0011] For example, in the second prior art, "Deep learning-assisted optimization for enhanced oil recovery and CO2 sequestration considering gas channeling constraints" (XYZhuang, WD Wang, YL Su, ZX Dai, Petroleum Science, 2025), a method for CO2 flooding and sequestration assisted by a proxy model based on graph attention networks is proposed. This method constructs an inter-well correlation graph structure based on numerical simulation data of CO2 flooding multiphase flow, integrates reservoir geological attributes, well network spatial distribution and development strategies, and analyzes the temporal dependency of the CO2 flooding process through a multi-layer stacked structure of graph attention layers and Transformers, thereby achieving dynamic prediction of development indicators such as cumulative oil production and CO2 sequestration.

[0012] The existing technology 2 has the following problems: it only focuses on macroscopic indicators at the well network scale, and does not perform fine modeling of the displacement front and saturation spatial distribution at the three-dimensional grid scale, which makes it impossible to achieve full-domain long-term dynamic prediction and difficult to support fine real-time control.

[0013] In summary, existing technical solutions cannot simultaneously meet the engineering requirements of high precision, high efficiency, strong robustness, physical compliance, and multi-physics field collaboration. Summary of the Invention

[0014] To address the above problems, this invention provides a multi-field intelligent prediction method for CO2 flooding based on frequency domain divide-and-conquer and physical constraints, comprising the following steps: A numerical simulation model of CO2-driven oil reservoir was constructed, and multiphysics three-dimensional time-series data samples were extracted and preprocessed from the CO2-driven oil reservoir numerical simulation model. An independent and parallel temporal processing branch is constructed for each physical field. The temporal features of the preprocessed data samples are extracted and compressed through the temporal processing branch to obtain the original temporal compressed feature vector of each physical field. The original time-series compressed feature vectors are mapped to the frequency domain space, and frequency domain noise filtering is performed. The features of different physical fields in the filtered frequency domain space are divided into high and low frequencies, and the high-frequency abrupt feature and the low-frequency background feature are dynamically fused to obtain the global frequency domain fusion feature. The global frequency domain fusion features are inverted into time-domain dynamic incremental features; The original time-series compressed feature vector and the time-domain dynamic incremental feature are adaptively weighted to obtain the fusion feature of each physical field. A multi-physics collaborative prediction model is constructed, which takes the fusion feature of each physical field as input and outputs the multi-field feature of the next time step. The training process of the multi-physics collaborative prediction model ensures the physical consistency of the predicted value through the fusion loss function embedded with mechanistic knowledge. Input the CO2 flooding engineering parameters of the target reservoir into the trained multiphysics collaborative prediction model, and output the CO2 flooding multiphysics prediction results.

[0015] Furthermore, a numerical simulation model of a CO2-enhanced oil reservoir is constructed, including the following steps: Based on the CO2 flooding engineering parameters of the target reservoir, a numerical simulation model of the CO2 flooding reservoir with multiple grid nodes in the target dimension is established.

[0016] Furthermore, multiphysics three-dimensional time-series data samples are extracted from the numerical simulation model of CO2-enhanced reservoirs and preprocessed, including the following steps: Obtain the output file of the numerical simulation model of CO2 oil displacement reservoir, and analyze and extract the multiphysics three-dimensional matrix data of each time step in the entire cycle of CO2 oil displacement. A time-series sliding window is used to construct a supervised learning training sample set. The multi-physics three-dimensional matrix data of consecutive preset time steps are used as the model input X. Each input sample is simultaneously labeled with the displacement stage label of the corresponding time window. The full multi-physics three-dimensional matrix data of the next time step is used as the model prediction target Y, thus completing the basic construction of the supervised sample. The supervised samples are divided into training, validation and test sets according to a set ratio. Each physics field is independently standardized to convert the three-dimensional matrix data of each physics field into a standard normal distribution. The completed training data is stored to obtain preprocessed multi-physics three-dimensional time series data samples.

[0017] Furthermore, the multi-physics three-dimensional matrix data includes pressure field, oil phase saturation field, water phase saturation field, CO2 gas phase saturation field, CO2 concentration field, permeability field, porosity field, and crude oil viscosity field.

[0018] Furthermore, the preprocessed multi-physics 3D time-series data samples are subjected to time-series feature extraction and compression through a time-series processing branch to obtain the original time-series compressed feature vector for each physics field, including the following steps: The three-dimensional tensor of the preset time step input to each time series processing branch is embedded with position information using a learnable three-dimensional position encoding matrix. The dimension of the encoding matrix is ​​completely consistent with the input tensor. The three-dimensional spatial coordinate information of the reservoir grid I, J, K and the time series information of the time step are synchronously embedded into the input features to obtain the position-encoded three-dimensional time series tensor. The three-dimensional temporal tensor after position encoding is divided into fixed-size blocks along the spatial dimension, resulting in multiple spatial blocks. Each spatial block corresponds to the temporal features of a preset time step. Each block is flattened and mapped to a token vector of a preset dimension through a linear mapping layer. The single physics data of each time step is converted into token vectors of multiple preset dimensions to obtain the token sequence of the preset time step. The token sequence at a preset time step is input into the VIT encoder, which is a multi-layer TransformerBlock stack, to complete the adaptive extraction of temporal-spatial coupling features. The VIT encoder outputs a token sequence of the same length as the input. The token sequence output by the VIT encoder is compressed into a single vector of a set dimension through global average pooling and a linear compression layer, which is used as the temporal compression feature of the physical field. Each temporal processing branch outputs an independent temporal compression feature vector of a set dimension, thus obtaining the original temporal compression feature vector of each physical field.

[0019] Furthermore, the original time-series compressed feature vector is mapped to the frequency domain space, and frequency domain noise filtering is performed, including the following steps: An independent frequency domain mapping branch is constructed for each physical field's temporal compressed feature vector; the temporal compressed feature vector of each physical field with a set dimension is increased in dimension through a linear layer, and then converted into a two-dimensional feature tensor through a reshape operation to obtain the preprocessed two-dimensional feature tensor. The preprocessed two-dimensional feature tensor is input into a multi-layered stacked FNO layer. Each FNO layer includes a Fast Fourier Transform (FFT) module, a frequency domain filtering and fusion module, and an Inverse Fast Fourier Transform (IFFT) module. The FFT module transforms the feature tensor in the time domain to the Fourier frequency domain. In the frequency domain, the frequency domain filtering and fusion module decomposes the original frequency domain features into high-frequency features, low-frequency features, and environmental noise features. A learnable weight matrix is ​​used to perform a linear transformation, converting the nonlinear evolution relationship in the time domain into a linear representation in the frequency domain. Finally, the IFFT module transforms the feature tensor back to the spatial domain, completing one frequency domain feature mapping. In the frequency domain transformation stage of each FNO layer, a learnable Gaussian filter kernel is embedded to smooth the frequency domain features and automatically remove high-frequency noise components in the frequency domain. After multiple FNO layers and noise filtering, the frequency domain feature vector corresponding to the set dimension of each physical field is output. Each physical field outputs an independent frequency domain feature vector, completing the full mapping from time domain features to frequency domain space.

[0020] Furthermore, the features of different physical fields in the filtered frequency domain are divided into high and low frequencies, and the high-frequency abrupt change features and low-frequency background features are dynamically fused to obtain global frequency domain fused features, including the following steps: For each physical field, the frequency domain feature vector of a given dimension is automatically decomposed into high-frequency feature components and low-frequency feature components through a learnable adaptive frequency domain partitioning threshold. The high-frequency feature components are converted into a token sequence through convolution transformation. Each high-frequency component of a physical field is an independent token. The token sequence is input into the VIT encoder. The VIT encoder automatically captures the correlation between high-frequency features of different physical fields through a self-attention mechanism, completes the cross-physical field deep fusion of key high-frequency features, and outputs the fused global high-frequency feature vector. The low-frequency feature components are converted into a token sequence, with each low-frequency component of a physical field being an independent token. The token sequence is then input into the ViT encoder to complete the cross-physical field deep fusion of the background low-frequency features and output the fused global low-frequency feature vector. The fused global high-frequency feature vector and the global low-frequency feature vector are concatenated along the channel dimension. The concatenated features are then reduced to a global frequency domain fused feature vector of a set dimension by a multilayer perceptron, thus obtaining the global frequency domain fused features.

[0021] Furthermore, the global frequency domain fusion features are inverted into time domain dynamic incremental features, including the following steps: The frequency domain features are inverted to the time domain features through the fast inverse Fourier transform of the Fourier neural operator, and then the dynamic incremental features of each physical field are generated through nonlinear transformation.

[0022] Furthermore, the original temporal compressed feature vector is adaptively weighted with the temporal dynamic incremental feature to obtain the fused feature of each physical field, including the following steps: Construct an independent gated weighting unit for each physical field; The gated weighting unit adaptively weights and fuses the original temporal compressed feature vector and the temporal dynamic incremental feature through the output weights to obtain the fused feature of each physical field.

[0023] Furthermore, the physical consistency constraint loss function is determined based on the data fitting loss and the inter-field coupling physical constraint loss.

[0024] Furthermore, the inter-field coupling physical constraint loss is constructed based on the seepage mechanics mass conservation equation and the saturation constraint equation, including saturation constraint, pressure-saturation coupling constraint, and concentration conservation constraint.

[0025] Furthermore, the saturation constraint includes ensuring that the sum of the saturations of the oil phase, water phase, and gas phase is always equal to 1 at each grid node, and penalizing saturation predictions that exceed the range of 0-1. The pressure and saturation coupling constraint includes constructing a coupling relationship between pressure gradient and fluid saturation change based on Darcy's law for multiphase seepage, and penalizing non-physical prediction results that violate seepage laws; Concentration conservation constraints include ensuring the total mass of CO2 is conserved within a closed reservoir system and penalizing predicted anomalous changes in CO2 concentration.

[0026] Furthermore, it also includes the following steps: CO2 enhanced oil recovery engineering parameters include geological parameters, fluid parameters, and operating condition parameters.

[0027] This invention also provides a frequency-domain divide-and-conquer and physically constrained CO2 enhanced oil recovery multi-field intelligent prediction system, comprising: The dataset construction module is used to build a numerical simulation model of CO2 flooded reservoirs, extract multiphysics three-dimensional time-series data samples from the CO2 flooded reservoir numerical simulation model and perform preprocessing. The parallel temporal feature extraction module is used to construct an independent parallel temporal processing branch for each physical field. The temporal processing branch extracts and compresses the temporal features of the preprocessed data samples to obtain the original temporal compressed feature vector of each physical field. The feature mapping module is used to map the original time-series compressed feature vectors to the frequency domain space and perform frequency domain noise filtering. The frequency division and fusion module is used to divide the features of different physical fields in the filtered frequency domain space into high and low frequencies, and dynamically fuse the high-frequency abrupt features and low-frequency background features to obtain global frequency domain fusion features. The frequency domain inversion and incremental generation module is used to invert global frequency domain fusion features into time domain dynamic incremental features; The model training module is used to adaptively weight the original time-series compressed feature vector and the time-domain dynamic incremental feature to obtain the fusion feature of each physical field, and construct a multi-physics collaborative prediction model. The fusion feature of each physical field is used as input to output the multi-field feature of the next time step. The training process of the multi-physics collaborative prediction model ensures the physical consistency of the predicted value through the fusion loss function embedded with mechanistic knowledge. The multiphysics prediction module is also used to input the CO2 flooding engineering parameters of the target reservoir into the trained multiphysics collaborative prediction model and output the CO2 flooding multiphysics prediction results.

[0028] The beneficial effects of this invention are: 1. This invention effectively addresses the core shortcomings of existing CO2 enhanced oil recovery (EOR) dynamic prediction technologies, achieving significant improvements in prediction accuracy, computational efficiency, robustness, and physical compliance. It also possesses strong scenario adaptability and application value. By constructing independent and parallel temporal feature extraction branches for each physical field, this invention can accurately capture the differentiated evolutionary characteristics of different physical fields, effectively avoiding feature overload. Combined with a high- and low-frequency divide-and-conquer frequency domain fusion architecture, it can accurately capture subtle changes in key regions such as the CO2 migration front and permeability abrupt change zones. Compared to traditional methods, the prediction error in key regions is reduced by more than 48%, and the average prediction R² (coefficient of determination) of core EOR dynamic features (pressure field, CO2 vapor saturation field, etc.) reaches over 0.96. The prediction results highly match the actual values ​​in reservoir numerical simulations, meeting the needs of refined engineering predictions.

[0029] 2. This invention uses frequency domain mapping technology to transform the strongly nonlinear evolution process of CO2 flooding in the time domain into a linear representation in the frequency domain. This fundamentally solves the problems of error accumulation and numerical oscillation that are prone to occur in long-term predictions of traditional time domain prediction models. Even for long-term predictions of CO2 injection over many years, it can still maintain stable prediction accuracy, and its robustness and generalization ability are greatly improved.

[0030] 3. This invention achieves adaptive coupling between the original time-domain features and the dynamically incremental features retrieved from the frequency domain through a gated physical residual fusion link. At the same time, it embeds the coupling mechanism between seepage mechanics fields into the loss function and constructs a physical consistency constraint loss function to achieve dual physical constraints. This ensures that the prediction results strictly follow the basic laws such as mass conservation, saturation constraint, and pressure-saturation coupling. The physical consistency compliance rate of the prediction results reaches 100%, which completely solves the pain point of traditional deep learning models being prone to non-physical interpretation and greatly improves the adaptability to field applications.

[0031] 4. This invention adopts a highly efficient feature extraction and frequency domain processing architecture. The time required for a single three-dimensional full reservoir multiphysics prediction is less than 4 seconds. Compared with traditional reservoir numerical simulation methods, the computational efficiency is improved by more than 4 orders of magnitude. It completely solves the bottleneck of long calculation cycle and high cost of traditional methods. It can realize batch optimization of CO2 oil displacement schemes and real-time dynamic control of injection parameters, providing strong support for rapid decision-making on the engineering site.

[0032] 5. The parallel branch architecture of this invention enables independent feature extraction of each physical field. Combined with high- and low-frequency divide-and-conquer processing in the frequency domain, it adapts to the fusion requirements of key abrupt change features and background stationary features, and can fully adapt to the complex features of multi-physical field coupling and strong spatiotemporal nonlinearity in CO2 oil recovery. Compared with the traditional shared weight feature extraction method, the targeting and effectiveness of feature representation are greatly improved. In addition, this method does not require adjustment of the core architecture. It only needs to adapt the input and output physical field types and corresponding physical constraint equations, and can be extended to various CCUS scenarios such as shale oil / gas reservoirs, conventional oil and gas reservoirs, deep saline aquifers, and coalbed methane reservoirs. At the same time, it can adapt to CO2 oil recovery and storage models with different grid dimensions and different injection methods, and has strong scenario adaptability and technical scalability, with broad application prospects.

[0033] Other features and advantages of the invention will be set forth in the following description, and will be apparent in part from the description, or may be learned by practicing the invention. The objects and other advantages of the invention may be realized and obtained by means of the structures pointed out in the description and the drawings. Attached Figure Description

[0034] To more clearly illustrate the technical solutions in the embodiments of the present invention or the prior art, the drawings used in the description of the embodiments or the prior art will be briefly introduced below. Obviously, the drawings described below are some embodiments of the present invention. For those skilled in the art, other drawings can be obtained based on these drawings without creative effort.

[0035] Figure 1 A schematic diagram of the main process of a multi-field intelligent prediction method for CO2 oil displacement based on frequency domain divide-and-conquer and physical constraints according to an embodiment of the present invention is shown. Figure 2 A detailed flowchart of a multi-field intelligent prediction method for CO2 flooding based on frequency domain divide-and-conquer and physical constraints according to an embodiment of the present invention is shown. Figure 3a A comparison chart of CO2 gas phase saturation prediction effects according to embodiments of the present invention is shown; Figure 3b A comparison chart of the oil phase saturation prediction effects according to an embodiment of the present invention is shown; Figure 3cA comparison chart of the water phase saturation prediction effects according to an embodiment of the present invention is shown; Figure 3d A comparison chart showing the predictive effect of pressure changes during CO2 oil displacement according to an embodiment of the present invention is shown; Figure 4 A comparison diagram of the predicted CO2 saturation distribution at different depths according to an embodiment of the present invention is shown; Figure 5 A comparison chart showing the prediction effect of concentration distribution of oil, gas and water at the same depth location according to an embodiment of the present invention is shown. Figure 6a This diagram illustrates the comparison between actual production data and predicted values ​​of the production oil-gas ratio time series of the validation set samples according to an embodiment of the present invention. Figure 6b This diagram illustrates the comparison between actual and predicted gas production time series data of the validation set samples according to an embodiment of the present invention. Figure 6c This diagram illustrates the comparison between actual production data and predicted values ​​of the production volume time series of the validation set sample according to an embodiment of the present invention; Figure 7 A schematic diagram of the structure and flow of a CO2 flooding multi-field intelligent prediction system based on frequency domain divide-and-conquer and physical constraints according to an embodiment of the present invention is shown. Detailed Implementation

[0036] To make the objectives, technical solutions, and advantages of the embodiments of the present invention clearer, the technical solutions of the embodiments of the present invention will be clearly and completely described below with reference to the accompanying drawings. Obviously, the described embodiments are only some embodiments of the present invention, not all embodiments. Based on the embodiments of the present invention, all other embodiments obtained by those skilled in the art without creative effort are within the scope of protection of the present invention.

[0037] It should be noted that the terms "first," "second," etc., used in this application are used to distinguish similar objects and are not necessarily used to describe a specific order or sequence. It should be understood that such data can be interchanged where appropriate for the embodiments of this application described herein.

[0038] This invention provides a multi-field intelligent prediction method and system for CO2 flooding based on frequency domain divide-and-conquer and physical constraints. It constructs a novel technical architecture encompassing time-domain parallel feature extraction, frequency domain divide-and-conquer fusion, time-domain inversion, and physical constraint prediction. By extracting independent features from multiple physics fields, performing frequency-domain linearization mapping, and fusing high and low frequencies, it overcomes the technical limitations of existing methods in feature capture and nonlinear processing, improving the precision of dynamic prediction for CO2 flooding. Simultaneously, it addresses the problems of strong grid dependence and low computational efficiency in traditional reservoir numerical simulation methods, achieving accurate, efficient, and robust prediction of multi-physics field dynamics during CO2 flooding. This provides technical support for scheme optimization, safety assessment, and real-time control in CCUS engineering.

[0039] like Figure 1 and Figure 2 As shown, a multi-field intelligent prediction method for CO2 flooding based on frequency domain divide-and-conquer and physical constraints includes the following steps: S1. Construct a numerical simulation model of CO2-driven oil reservoir, extract multiphysics three-dimensional time-series data samples from the CO2-driven oil reservoir numerical simulation model and perform preprocessing.

[0040] The construction of a CO2 flooding reservoir numerical simulation model includes the following steps: Based on the CO2 flooding engineering parameters of the target reservoir, establish a CO2 flooding reservoir numerical simulation model with a target dimension of 10×128×128 (K×I×J) and multiple grid nodes (e.g., 163,840) to ensure the computational stability and convergence of the numerical model and stably output the full-cycle time series data required for neural network training.

[0041] The process of extracting multiphysics three-dimensional time-series data samples from the numerical simulation model of CO2-enhanced oil reservoirs and preprocessing them includes the following steps: S11. Multiphysics 3D Data Extraction: Obtain the output file of the CO2 flooding reservoir numerical simulation model, and parse and extract the multiphysics 3D matrix data of each time step in the entire CO2 flooding cycle, including 8 core physical fields such as pressure field, oil phase saturation field, water phase saturation field, CO2 gas phase saturation field, CO2 concentration field, permeability field, porosity field, and crude oil viscosity field. Each physical field is a 10×128×128 three-dimensional tensor that matches the mesh of the CO2 flooding reservoir numerical simulation model.

[0042] S12. Construction of supervised samples: A time-series sliding window is used to construct a supervised learning training sample set. The multi-physics three-dimensional matrix data of consecutive preset time steps (e.g., 6 time steps) are used as the model input X. Each input sample is simultaneously labeled with the displacement stage label of the corresponding time window. The full multi-physics three-dimensional matrix data of the next time step is used as the model prediction target Y, thus completing the basic construction of supervised samples.

[0043] It should be noted that the length of the timing sliding window can be adjusted to other lengths such as 3 steps or 10 steps, depending on the forecasting requirements.

[0044] S13. Dataset Partitioning and Preprocessing: The supervised samples are divided into training, validation, and test sets according to a set ratio (e.g., a 7:2:1 ratio). Z-score normalization is used to perform independent normalization on each physics field, converting the 3D matrix data of each physics field into a standard normal distribution (mean 0, standard deviation 1), eliminating the numerical magnitude differences between different physics fields while preserving the independent distribution characteristics of each physics field. The trained data is stored in HDF5 format with a memory optimization mechanism to avoid memory overflow during reading due to the large size of the high-dimensional tensor dataset, thus obtaining preprocessed multi-physics 3D time-series data samples.

[0045] It should be noted that the standardization method can be replaced by other data normalization methods such as Min-Max normalization and logarithmic normalization; the types of high-dimensional features in the training dataset can be added or deleted according to the characteristics of the target reservoir scenario; and the displacement stage annotation can be refined into more stages according to the development process.

[0046] S2. Construct independent and parallel time-series processing branches for each physics field. Through these branches, extract and compress the time-series features of the preprocessed multi-physics 3D time-series data samples to obtain the original time-series compressed feature vector for each physics field.

[0047] In this embodiment of the invention, an independent and parallel temporal processing branch is constructed for each physical field. A VisionTransformer with embedded reservoir 3D mesh location coding is used to adaptively extract and compress temporal features, accurately capturing the correlation between reservoir temporal evolution and spatial heterogeneity without introducing any physical constraints or prior knowledge. The process includes the following steps: S21. Parallel temporal processing branch design: For the multiple input physical fields (e.g., the 8 physical fields mentioned in step S1), construct completely independent parallel temporal processing branches with non-shared weights. Each temporal processing branch has a completely identical structure and performs feature extraction only on the temporal 3D data of the corresponding physical field. There is no information interaction between branches.

[0048] S22. 3D Mesh Position Encoding Embedding: The 3D tensor (dimension 6×10×128×128) of the preset time step (e.g., 6 time steps) input to each time series processing branch is embedded with position information using a learnable 3D position encoding matrix. The dimension of the encoding matrix is ​​completely consistent with the input tensor. The I, J, K 3D spatial coordinate information of the reservoir mesh and the time series information of the time step are synchronously embedded into the input features, solving the problem of Transformer's insensitivity to position information. After encoding, the output dimension remains 6×10×128×128, and the position-encoded 3D time series tensor is obtained.

[0049] S23. Tensor Blocking and Linear Mapping: The position-encoded 3D temporal tensor is divided into fixed-size blocks (patches) along the spatial dimension, resulting in multiple spatial blocks. Each spatial block corresponds to the temporal features of a preset time step (e.g., 6 time steps). After flattening each block, it is mapped to a token vector of a preset dimension through a linear mapping layer, converting the single physics data of each time step into token vectors of multiple preset dimensions, thus obtaining the token sequence of the preset time steps.

[0050] For example, the position-encoded 3D temporal tensor is divided into fixed-size patches along the spatial dimension. The patch size is set to 2×16×16 (K×I×J), corresponding to the original 10×128×128 spatial dimension. After patching, 5×8×8=320 spatial patch patches are obtained, and each patch patch corresponds to the temporal features of 6 time steps. After flattening each patch, it is mapped to a 256-dimensional token vector through a linear mapping layer. Finally, the single physics data of each time step is converted into 320 256-dimensional tokens, and the total token sequence length of 6 time steps is 1920.

[0051] S24. Adaptive extraction of temporal features using VIT (Vision Transformer): The token sequence at a preset time step is input into a VIT encoder consisting of multiple stacked TransformerBlocks. Each time step corresponds to one TransformerBlock. Each TransformerBlock includes a multi-head self-attention mechanism (MSA) layer, a layer normalization (LayerNorm) layer, and a feedforward neural network with a hidden layer dimension of 1024. For example, the number of heads in the multi-head self-attention mechanism layer is the same as the number of physical fields (e.g., 8 physical fields, and the number of heads is also set to 8), or they can be different. The self-attention mechanism automatically captures the long-range dependencies in the temporal dimension and the heterogeneous associations in the spatial dimension of the token sequence, completing the adaptive extraction of temporal-spatial coupling features. The VIT encoder outputs a token sequence of the same length as the input.

[0052] It should be noted that the VIT architecture of the temporal feature extraction module can be replaced with other attention variant architectures such as Swin, Transformer, and TimeSformer.

[0053] S25. Feature Compression and Output: The token sequence output by the VIT encoder is compressed into a single vector of a set dimension (e.g., dimension 512) through global average pooling (GAP) and a linear compression layer, which is used as the temporal compression feature of the physical field. Each temporal processing branch outputs an independent temporal compression feature vector of a set dimension, thus obtaining the original temporal compression feature vector of each physical field.

[0054] S3. Map the original time-series compressed feature vectors to the frequency domain and perform frequency domain noise filtering, including: increasing the dimensionality of the original time-series compressed feature vectors of each physical field, mapping them to the frequency domain through the Fast Fourier Transform module in the Fourier Neural Operator (FNO), transforming the time-domain nonlinear evolution into a frequency-domain linear representation, and embedding a noise filtering module to remove environmental noise in the frequency domain of different physical field features. Specifically, this includes the following steps: S31. Feature Upscaling and Preprocessing: For each physical field's temporal compressed feature vector, an independent frequency domain mapping branch is constructed. The weights of the frequency domain mapping branches are not shared and there is no information exchange. The temporal compressed feature vector of each physical field with a set dimension is upscaled through a linear layer (e.g., from 512 dimensions to 1024 dimensions), and then converted into a 32×32 two-dimensional feature tensor through a reshape operation to adapt to the input dimension requirements of FNO. The upscaling process only uses linear mapping and does not introduce nonlinear activation, thus preserving the linear correlation of the original features.

[0055] S32. FNO (Fourier Neural Operator) Frequency Domain Transformation: The preprocessed two-dimensional feature tensor is input into a multi-layer (e.g., 4-layer) stacked FNO layer. Each FNO layer includes a Fast Fourier Transform (FFT) module, a frequency domain filtering and fusion module, and an Inverse Fast Fourier Transform (IFFT) module. The Fourier Transform module completely transforms the feature tensor in the time domain to the Fourier frequency domain. In the frequency domain, the frequency domain filtering and fusion module decomposes the original frequency domain features into high-frequency features, low-frequency features, and environmental noise features. A learnable weight matrix is ​​used to complete the linear transformation, converting the nonlinear evolution relationship in the time domain into a linear representation in the frequency domain. Finally, the Inverse Fourier Transform module transforms it back to the spatial domain, completing one frequency domain feature mapping. The stacking of multiple FNO layers enables the layer-by-layer extraction of multi-scale frequency domain features.

[0056] S33, Frequency Domain Noise Filtering Module: In the frequency domain transformation stage of each FNO layer, a learnable Gaussian filter kernel is embedded to smooth the frequency domain features and automatically remove high-frequency noise components in the frequency domain. The weights of the filter kernel are adaptively updated as the model is trained. No fixed filtering threshold is set to avoid human intervention in the feature distribution. The filtered frequency domain features then enter the subsequent fusion transformation stage.

[0057] S34. Frequency Domain Feature Output: After multiple FNO layers and noise filtering, the frequency domain feature vector of the set dimension (e.g., 32×32) corresponding to each physical field is output. Each physical field outputs an independent frequency domain feature vector, completing the full mapping from time domain features to frequency domain space.

[0058] It should be noted that the Fourier transform in the frequency domain mapping and inversion module can be replaced by other frequency domain transformation methods such as wavelet transform and discrete cosine transform.

[0059] S4. Divide the features of different physical fields in the filtered frequency domain space into high and low frequencies, and dynamically fuse the high-frequency abrupt change features and the low-frequency background features to obtain the global frequency domain fusion features.

[0060] This step designs a multi-frequency domain hierarchical adaptive filtering module, which separates the features of different physical fields in the frequency domain into high-frequency and low-frequency components for separate processing. Key high-frequency features such as the CO2 migration front and permeability abrupt change regions are dynamically fused using VIT, and background low-frequency features such as oil, gas, and water distribution characteristics and pressure distribution are also dynamically fused using VIT. This achieves adaptive deep coupling of multi-physical field frequency domain features, specifically including the following steps: S41. Adaptive High- and Low-Frequency Decomposition and Filtering: For the frequency domain feature vector of each physical field with a set dimension (32×32=1024), it is automatically decomposed into high-frequency feature components and low-frequency feature components through a learnable adaptive frequency domain partitioning threshold. The partitioning threshold is adaptively updated. After decomposition, the high-frequency components and low-frequency components are subjected to independent adaptive filtering. The high-frequency components are bandpass filtered to retain edge abrupt features, and the low-frequency components are lowpass filtered to smooth global background features. The filter kernel weights are all learnable parameters.

[0061] S42. Cross-field fusion of high-frequency features based on VIT: The high-frequency feature components are converted into token sequences through convolution transformation. Each high-frequency component of a physical field is an independent token with a token dimension of 512 and a sequence length of 8. The token sequence is input into a VIT encoder composed of 3 stacked TransformerBlocks, with an 8-head attention mechanism and an FFN hidden layer dimension of 2048. The correlation between high-frequency features of different physical fields is automatically captured through the self-attention mechanism to complete the deep fusion of key high-frequency features across physical fields and output a fused 1024-dimensional global high-frequency feature vector.

[0062] S43. Cross-field fusion of low-frequency features based on VIT: The low-frequency feature components are converted into token sequences. Each low-frequency component of a physical field is an independent token with a token dimension of 512 and a sequence length of 8. The token sequence is input into a ViT encoder consisting of 3 layers of TransformerBlocks with the same high-frequency branch structure and completely independent weights to complete the deep fusion of background low-frequency features across physical fields and output a fused 1024-dimensional global low-frequency feature vector.

[0063] S44. Feature Fusion: The fused global high-frequency feature vector and the global low-frequency feature vector are concatenated in the channel dimension. The concatenated features are reduced to a global frequency domain fusion feature vector of a set dimension (e.g., 1024 dimensions) through a multilayer perceptron (MLP), thus completing the adaptive deep coupling of high and low frequency division and control of multi-physics fields and obtaining global frequency domain fusion features.

[0064] S5. Inverting the global frequency domain fusion features into time domain dynamic incremental features includes: completing the inversion of frequency domain features to time domain features through the inverse IFFT transform of the FNO Fourier neural operator, and then generating the dynamic incremental features of each physical field through nonlinear transformation. Specifically, this includes the following steps: S51. Frequency domain feature preprocessing: The global frequency domain fusion features are converted into feature vectors of dimension 1024 through a linear mapping layer, and then converted into two-dimensional feature tensors of 32×32 through a reshape operation to adapt to the input dimension requirements of the FNO inverse transform layer.

[0065] S52, FNO frequency domain inverse transform: The preprocessed two-dimensional feature tensor is input into the FNO inverse transform layer to complete the inverse transform from the frequency domain to the time domain, and outputs a high-dimensional feature tensor in the time domain space.

[0066] S53. Nonlinear Transformation and Dynamic Incremental Feature Generation: The high-dimensional feature tensor in the time domain output by the inverse FNO transform is processed through three consecutive convolutional nonlinear transformation layers. The convolutional kernel size is set to 3×3, the stride is 1, the padding is 1, and the activation function is GELU. Finally, through a linear mapping layer, the feature tensor is transformed into an 8×10×128×128 three-dimensional tensor that matches the dimension of the input physical field, which serves as the dynamic incremental feature of each physical field, completing the complete inversion from frequency domain features to time domain dynamic incremental features.

[0067] S6. Adaptively weight the original temporal compressed feature vector with the temporal dynamic incremental feature to obtain the fusion feature of each physical field. Use convolutional layers to construct a multi-physics collaborative prediction model. Take the fusion feature of each physical field as input and output the multi-field feature of the next time step. During the training process of the multi-physics collaborative prediction model, the fusion loss function embedded by the mechanism knowledge ensures the physical consistency of the predicted values ​​of pressure, saturation, concentration and other parameters.

[0068] This step constructs a gated dynamic weighted residual link for each physics field, adaptively weighting and fusing the original temporal compressed feature vector extracted in step S2 with the temporal dynamic incremental feature generated in step S5; finally, through a multi-physics field collaborative prediction model, the inter-field coupling mechanism is embedded in the loss, and the next time step result of each physics field is output synchronously to ensure the physical consistency of pressure, saturation, and concentration.

[0069] The process of adaptively weighting the original temporal compressed feature vector with the temporal dynamic incremental feature to obtain the fused feature of each physical field includes the following steps: S61. Construct an independent gated weighting unit for each physical field. The gated unit is implemented using a 2-layer fully connected network with the activation function being Sigmoid. The output is a dynamic weighting weight between [0,1].

[0070] S62. The gated weighting unit adaptively weights and fuses the original temporal compressed feature vector extracted in step S2 and the temporal dynamic incremental feature generated in step S5 using the output weights to obtain the fused feature of each physical field, as follows:

[0071] in, For the first i The fusion characteristics of individual physical fields The dynamic weights output by the gating unit. For the first extracted in step S2 i Original temporal compressed feature vectors of physical fields The first generated in step S5 i The time-domain dynamic incremental characteristics of a physical field.

[0072] The instruction manual is required; the gated weighting unit can be replaced with other gated structures such as GRU or LSTM.

[0073] The construction of the multiphysics collaborative prediction model includes the following steps: constructing a collaborative prediction head consisting of two 3D convolutional layers, with the kernel size set to 2×2×2, stride of 1, padding of 1, and GELU activation function. The output of the gated weighted unit is used as the input of the collaborative prediction head to obtain the multiphysics collaborative prediction model. The final output of the multiphysics collaborative prediction model is a three-dimensional tensor with dimensions of 8×10×128×128, corresponding to the prediction results of the eight physics fields in the next time step.

[0074] The physical consistency constraint loss function is determined based on the data fitting loss and the inter-field coupling physical constraint loss, as follows:

[0075] in, This is the data fitting loss, which is the mean square error loss between the predicted and the true values, used to ensure the accuracy of data fitting. The weighting coefficient for the physical constraint loss is set to 0.2; This represents the physical constraint loss due to inter-field coupling.

[0076] Among them, the inter-field coupling physical constraint loss is constructed based on the seepage mechanics mass conservation equation and saturation constraint equation, specifically including saturation constraint, pressure-saturation coupling constraint, and concentration conservation constraint.

[0077] The saturation constraint includes ensuring that the sum of the saturation of the oil phase, water phase, and gas phase is always equal to 1 at each grid node, and penalizing saturation prediction results that exceed the range of 0-1.

[0078] The pressure and saturation coupling constraint includes constructing a coupling relationship between pressure gradient and fluid saturation change based on Darcy's law for multiphase seepage, and penalizing non-physical prediction results that violate seepage laws; Concentration conservation constraints include ensuring the total mass of CO2 is conserved within a closed reservoir system and penalizing predicted anomalous changes in CO2 concentration.

[0079] By embedding the physical constraint loss into the inter-field coupling mechanism, the physical consistency of pressure, saturation, and concentration in the prediction results is ensured, thus avoiding the generation of non-physical interpretations.

[0080] The process of inputting the fused features of each physics field into the multiphysics collaborative prediction model for training includes the following steps: Training employs the RMSprop optimizer (base learning rate 2e-4, decay factor 0.85, eps=1.5e-8), and enables FP16 mixed precision and data parallel training mode, combined with gradient accumulation and gradient clipping to ensure training stability.

[0081] The learning rate adopts a cosine annealing plus 5 rounds of linear warm-up strategy, with a total training of 100 epochs and an early stopping mechanism enabled (training is terminated when the total loss of the validation set does not improve for 10 consecutive rounds); the training batch is 32 and the validation batch is 64. Data loading is enhanced by 8 threads and multi-processing, and overfitting is suppressed by combining Dropout (inactivation rate of 0.1) and L2 weight decay (coefficient 1e-5).

[0082] The system maintains real-time logs during training, and after each training epoch, it saves the model weights with the best accuracy.

[0083] It should be noted that the model optimizer can be replaced with other deep learning optimizers such as SGD, Adam, and AdamW; the learning rate scheduling strategy can be replaced with other scheduling strategies such as StepLR; the training hyperparameters can be adapted to the hardware environment and dataset size; and the weight coefficients of the physical constraint loss can be adapted to the scenario requirements.

[0084] S7. Input the CO2 flooding engineering parameters of the target reservoir into the trained multiphysics collaborative prediction model and output the CO2 flooding multiphysics prediction results. The CO2 flooding engineering parameters include geological parameters, fluid parameters and operating condition parameters.

[0085] By inputting the geological parameters, fluid parameters, and operating condition parameters of the target reservoir into the trained model, the reservoir saturation field and other characteristics under different injection times can be quickly predicted, thus fully depicting the dynamic evolution process of CO2 oil displacement.

[0086] S8. Evaluation of the accuracy and progress of CO2 flooding multiphysics prediction results, as follows: The accuracy of dynamic prediction is evaluated using indicators such as MAE, MSE, and R²; the CO2 flooding prediction results are plotted using visualization methods to provide a predictive basis for subsequent flooding potential assessment, storage safety analysis, and injection parameter optimization.

[0087] The deep learning architecture design of this invention enables rapid prediction of the dynamics of multiple physics fields in CO2 flooding, meeting the high-efficiency computing requirements for real-time engineering decision-making and scheme optimization. It solves the problems of insufficient capture of the differentiated features of multiple physics fields and low prediction accuracy of key areas in existing deep learning proxy models. Through independent and parallel VIT time series feature extraction branches, it accurately captures the evolution law of each physics field and the key regional features such as the CO2 migration front and permeability mutation zone, solving the problems of error accumulation and numerical oscillation in long-term prediction of existing methods. By using FNO frequency domain mapping, the strong nonlinear evolution in the time domain is transformed into a linear representation in the frequency domain, which improves the stability and robustness of the model's long-term predictions and solves the problems of existing models' prediction results being prone to non-physical interpretations and difficulty in ensuring physical consistency between multiple physics fields.

[0088] By using gated physical residual fusion and inter-field coupling mechanism constraints, the prediction results are ensured to conform to the basic laws of seepage mechanics, improving the engineering adaptability of the model. Finally, through the innovative design of parallel VIT time-series feature extraction architecture, FNO frequency domain feature mapping and divide-and-conquer fusion link, combined with gated physical residual fusion and multi-physics collaborative prediction head, dual optimization at the feature level and result level is achieved. This enables accurate, efficient and robust prediction of multi-physics dynamics during CO2 oil displacement, providing technical support for CCUS engineering scheme optimization, safety assessment and real-time control.

[0089] Figure 3a , 3b Figures 3c and 3d are verification diagrams of the CO2 oil displacement prediction effect, which verify the prediction accuracy of the method of the present invention for the dynamic evolution of reservoir physical characteristics during CO2 oil displacement. Figure 3a The graph shows a comparison of the prediction results of CO2 gas phase saturation. The horizontal axis represents the actual values ​​on different grids obtained from reservoir numerical simulation, and the vertical axis represents the predicted values ​​output by the method in this embodiment of the invention. Figure 3b The graph shows a comparison of the prediction results for oil phase saturation. The horizontal axis represents the actual value obtained from the reservoir numerical simulation, and the vertical axis represents the predicted value output by the method of this invention. Figure 3c The graph shows a comparison of the prediction results of water phase saturation. The horizontal axis represents the actual value obtained from the reservoir numerical simulation, and the vertical axis represents the predicted value output by the method of this invention. Figure 3d The graph shows a comparison of the predicted pressure changes during CO2 flooding. The horizontal axis represents the actual pressure value obtained from reservoir numerical simulation, and the vertical axis represents the predicted value output by the method of this invention. The black dashed line in the graph is the ideal fitting line of y=x, and the colored scatter points in each subgraph are the corresponding data points of the predicted results and actual values ​​of the validation set.

[0090] from Figure 3a , 3b As can be seen from 3c and 3d, the data points in all subgraphs closely fit the ideal fitting line of y=x, with no obvious system deviation. This indicates that the dual-branch fusion deep learning model constructed in this embodiment of the invention has high prediction accuracy for the CO2 oil displacement process. The prediction results are consistent with the actual values ​​of numerical simulation, and can accurately depict the dynamic process of CO2 displacing crude oil in the oil reservoir.

[0091] like Figure 4 As shown, Figure 4The comparison chart shows the prediction results of CO2 saturation distribution at different depths, verifying the accuracy of the method in this embodiment of the invention in predicting gas phase distribution during CO2 flooding. It displays the comparison results of the actual values, predicted values, and relative errors of the planar distribution of gas phase saturation in the target reservoir at K=2 and K=4 layers. Each group contains three sub-figures: the left sub-figure of each group shows the actual distribution value obtained from the reservoir numerical simulation; the middle sub-figure of each group shows the prediction results of the method of this invention. It can be seen that the output prediction is highly consistent with the actual value of the numerical simulation, with no obvious visual difference; the maximum relative error at K=2 layer is less than 5.341%, and the maximum relative error at K=4 layer is less than 3.918%. The overall error is only concentrated in the local grid area at the CO2 displacement front, and the prediction error of the main reservoir area is close to 0.

[0092] pass Figure 4 The results of the comparison show that the method of the present invention can accurately characterize the flow pattern and spatial distribution characteristics of CO2 in a plane, and has excellent prediction accuracy for CO2 distribution at different depths.

[0093] like Figure 5 As shown, Figure 5 This chart compares the predicted concentration distribution of the three phases (oil, gas, and water) at the same depth, verifying the accuracy of the method in this invention for predicting the physical properties of multiple reservoirs during CO2 flooding. It displays the comparison of the actual values, model predictions, and relative errors of the target reservoir model at depth K=3. The chart is divided into three groups: the first row contains three sub-figures showing the actual gas, water, and oil phase saturation distributions generated by numerical simulation; the second row contains three sub-figures showing the gas, water, and oil phase saturation distributions predicted by this method; and the third row contains three sub-figures showing the relative error distributions of the predicted values ​​in terms of the gas, water, and oil phase saturation distributions. The maximum relative error for predicting gas and water phase saturation is less than 3.5%, and the maximum relative error for predicting oil phase saturation does not exceed 4.85%.

[0094] Figure 5 The results show that the saturation distribution field and local enrichment regions of crude oil predicted by the method of the present invention are in high agreement with the actual values ​​of numerical simulation, with no significant visual difference. The prediction error is concentrated only in the local grid regions of the migration front and the heterogeneous boundary of the reservoir where CO2 concentration changes abruptly, and the prediction relative error of the main reservoir region is close to 0. The results indicate that the method of the present invention can accurately fit the dynamic characteristics of CO2-driven oil displacement.

[0095] Figure 6a , 6bThe figure 6c shows a comparison of the prediction results of CO2 flooding production response time series data, which verifies the prediction accuracy of the method of the present invention on the real reservoir production dynamic time series data. It also shows the comparison results of production oil-gas ratio, production fluid rate and production gas rate in one of the validation set samples. Figure 6a In the diagram, the solid line represents the actual value of the production oil-gas ratio time series obtained from reservoir numerical simulation calculations, and the scatter points represent the predicted values ​​output by the method of this invention. Figure 6b In the diagram, the solid line represents the actual gas production value obtained from the numerical simulation calculation of the reservoir, and the scatter points represent the predicted values ​​output by the method in this embodiment of the invention. Figure 6c In the diagram, the solid line represents the actual production volume obtained from reservoir numerical simulation calculations, and the scatter points represent the predicted values ​​output by the method in this embodiment of the invention.

[0096] Figure 6a , 6b The results in section 6c show that the dynamic prediction of the production capacity response data for production wells is in high agreement with the actual values ​​from the numerical simulation, with no significant difference. These results demonstrate that the method described in this embodiment of the invention is effective for accurately predicting dynamic data of CO2-enhanced oil production in the field.

[0097] Based on the above frequency domain divide-and-conquer and physical constraint-based multi-field intelligent prediction method for CO2 enhanced oil recovery, such as... Figure 7 As shown, this embodiment of the invention also provides a CO2 enhanced oil recovery multi-field intelligent prediction system based on frequency domain divide-and-conquer and physical constraints, including a dataset construction module, a parallel time-series feature extraction module, a feature mapping module, a frequency division and fusion module, a frequency domain inversion and incremental generation module, a model training module, and a multi-physics prediction module.

[0098] The dataset construction module is used to build a numerical simulation model of CO2 flooding reservoirs, extract multiphysics three-dimensional time-series data samples from the CO2 flooding reservoir numerical simulation model and perform preprocessing.

[0099] The parallel temporal feature extraction module is used to construct an independent parallel temporal processing branch for each physical field. The temporal processing branch extracts and compresses the temporal features of the preprocessed data samples to obtain the original temporal compressed feature vector of each physical field.

[0100] The feature mapping module is used to map the original time-series compressed feature vectors to the frequency domain space and perform frequency domain noise filtering.

[0101] The frequency division and fusion module is used to divide the features of different physical fields in the filtered frequency domain space into high and low frequencies, and dynamically fuse the high-frequency abrupt change features and the low-frequency background features to obtain the global frequency domain fusion features.

[0102] The frequency domain inversion and incremental generation module is used to invert global frequency domain fusion features into time domain dynamic incremental features.

[0103] The model training module is used to adaptively weight the original temporal compressed feature vector and the temporal dynamic incremental feature to obtain the fusion feature of each physical field, and to construct a multi-physics collaborative prediction model. The fusion feature of each physical field is input into the multi-physics collaborative prediction model for training. The loss function of the prediction model adopts the physical consistency constraint loss function.

[0104] The multiphysics prediction module is also used to input the CO2 flooding engineering parameters of the target reservoir into the trained multiphysics collaborative prediction model and output the CO2 flooding multiphysics prediction results.

[0105] For example, the frequency domain divide-and-conquer and physically constrained CO2 flooding multi-field intelligent prediction method and system of this invention are applied to CO2 flooding and storage in a low-permeability oil reservoir in a certain area. The specific implementation process is as follows: 1. Data Acquisition and Preprocessing: Based on the CO2 flooding engineering parameters (including geological parameters, fluid parameters, and operating condition parameters) of the target shale oil reservoir, a 10×128×128 grid numerical simulation model of CO2 flooding in the shale oil reservoir was constructed. Numerical simulation calculations were completed using a reservoir numerical simulator, and the output result file and input feature parameters were summarized. High-dimensional feature data of eight core physical fields were parsed and extracted from the 3D grid. Data missing issues were checked and verified, outliers were removed, and inverse distance weighted interpolation was used to supplement missing data, thus completing data preprocessing.

[0106] 2. Construction and Standardization of Time-Series Datasets: A 6-step time-series sliding window was used to construct multi-physics 3D time-series data samples, with simultaneous labeling during the displacement phase. A total of 20,420 valid samples were generated, of which 4,084 were randomly selected as the test set, 2,042 as the validation set, and the remaining 14,294 as the training set. HDF5 format was used for training data storage. During data parsing, memory was cleaned up and released in real time to avoid memory overflow due to the large size of the model data. The Z-score standardization method was used to perform independent standardization processing for each physics field.

[0107] 3. Model Construction and Training: A multi-physics collaborative prediction model based on frequency domain divide-and-conquer and physical consistency constraints was constructed. Hyperparameter configuration was performed according to step S6 to complete the training setup, and an NVIDIA A100 GPU was used for model training. During training, a test set validation was performed at each epoch, and the average prediction accuracy for the dynamic fields of CO2 concentration and oil saturation was monitored. Early termination was triggered when the validation set metrics did not improve for 10 consecutive epochs. Training finally terminated at the 58th epoch, and the model weights with the optimal dynamic field average MSE and R² were saved.

[0108] 4. Model Accuracy Validation: Based on the test set data, the model achieves an average R² of 0.9821 for core dynamic features of oil displacement, including CO2 concentration field, oil phase saturation field, water phase saturation field, and pressure field. Specifically, the predicted R² for the pressure field is 0.9837, for the CO2 gas phase saturation field it is 0.9895, for the oil phase saturation field it is 0.9742, and for the water phase saturation field it is 0.9680. The prediction accuracy meets the requirements of engineering applications. For fracture regions with abrupt permeability changes and the CO2 displacement front, the prediction error is reduced by 48.3% compared to the traditional spatiotemporal fusion UNet model, with no significant numerical oscillations.

[0109] Qualitative verification results: Professional visualization processing was performed on the test set samples, and cloud maps comparing the actual and predicted values ​​of CO2 gas phase saturation, oil phase saturation field, and pressure field were plotted with fixed K layers. The results show that the model can accurately capture the spatial distribution of the CO2 migration front, and there is no visually discernible difference between the predicted results and the actual values ​​of numerical simulation. The percentage error is generally less than 3.5%, with only a small error appearing in the crack area where permeability changes abruptly, demonstrating excellent spatial distribution prediction capabilities. In the long-term 30-year injection cycle prediction, there was no error accumulation or non-physical explanation, and the physical consistency compliance rate reached 100%.

[0110] 5. Dynamic Prediction Application of CO2 Displacement Effect: Based on a trained multi-physics collaborative prediction model, and inputting CO2 flooding engineering parameters of the target reservoir, dynamic prediction of CO2 injection process over 30 years is completed. A single full-field prediction takes less than 4 seconds, representing an efficiency improvement of over 12,000 times compared to traditional numerical simulation. Based on the prediction results, core indicators such as the spatial distribution of CO2, sweep efficiency, oil recovery rate, and proportion of stored crude oil are clarified at different CO2 injection times. Oil displacement potential assessment and storage safety analysis are completed, and injection pressure and injection volume parameters are optimized, providing decision support for on-site engineering implementation.

[0111] The technical solutions of this invention can be extended to other underground porous media multiphase seepage scenarios such as dynamic prediction of oil and gas reservoir development, prediction of groundwater pollutant transport, and multi-field dynamic prediction of geothermal development. Only the physical field types and corresponding physical constraint equations of the input and output need to be adjusted, while the core method architecture remains unchanged. All of these are alternative solutions to this invention.

[0112] Although the present invention has been described in detail with reference to the foregoing embodiments, those skilled in the art should understand that modifications can still be made to the technical solutions described in the foregoing embodiments, or equivalent substitutions can be made to some of the technical features; and these modifications or substitutions do not cause the essence of the corresponding technical solutions to deviate from the spirit and scope of the technical solutions of the embodiments of the present invention.

Claims

1. A multi-field intelligent prediction method for CO2 enhanced oil recovery based on frequency domain divide-and-conquer and physical constraints, characterized in that, Includes the following steps: A numerical simulation model of CO2-driven oil reservoir was constructed, and multiphysics three-dimensional time-series data samples were extracted and preprocessed from the CO2-driven oil reservoir numerical simulation model. An independent and parallel temporal processing branch is constructed for each physical field. The temporal features of the preprocessed data samples are extracted and compressed through the temporal processing branch to obtain the original temporal compressed feature vector of each physical field. The original time-series compressed feature vectors are mapped to the frequency domain space, and frequency domain noise filtering is performed. The features of different physical fields in the filtered frequency domain space are divided into high and low frequencies, and the high-frequency abrupt feature and the low-frequency background feature are dynamically fused to obtain the global frequency domain fusion feature. The global frequency domain fusion features are inverted into time-domain dynamic incremental features; The original time-series compressed feature vector and the time-domain dynamic incremental feature are adaptively weighted to obtain the fusion feature of each physical field. A multi-physics collaborative prediction model is constructed, which takes the fusion feature of each physical field as input and outputs the multi-field feature of the next time step. The training process of the multi-physics collaborative prediction model ensures the physical consistency of the predicted value through the fusion loss function embedded with mechanistic knowledge. Input the CO2 flooding engineering parameters of the target reservoir into the trained multiphysics collaborative prediction model, and output the CO2 flooding multiphysics prediction results.

2. The frequency domain divide-and-conquer and physically constrained multi-field intelligent prediction method for CO2 oil displacement according to claim 1, characterized in that, Constructing a numerical simulation model for CO2-enhanced oil reservoirs includes the following steps: Based on the CO2 flooding engineering parameters of the target reservoir, a numerical simulation model of the CO2 flooding reservoir with multiple grid nodes in the target dimension is established.

3. The frequency domain divide-and-conquer and physically constrained multi-field intelligent prediction method for CO2 oil displacement according to claim 1, characterized in that, Multiphysics three-dimensional time-series data samples were extracted from the numerical simulation model of CO2-enhanced oil reservoirs and preprocessed, including the following steps: Obtain the output file of the numerical simulation model of CO2 oil displacement reservoir, and analyze and extract the multiphysics three-dimensional matrix data of each time step in the entire cycle of CO2 oil displacement. A time-series sliding window is used to construct a supervised learning training sample set. The multi-physics three-dimensional matrix data of consecutive preset time steps are used as the model input X. Each input sample is simultaneously labeled with the displacement stage label of the corresponding time window. The full multi-physics three-dimensional matrix data of the next time step is used as the model prediction target Y, thus completing the basic construction of the supervised sample. The supervised samples are divided into training, validation and test sets according to a set ratio. Each physics field is independently standardized to convert the three-dimensional matrix data of each physics field into a standard normal distribution. The completed training data is stored to obtain preprocessed multi-physics three-dimensional time series data samples.

4. The frequency domain divide-and-conquer and physically constrained multi-field intelligent prediction method for CO2 oil displacement according to claim 3, characterized in that, Multiphysics three-dimensional matrix data, including pressure field, oil phase saturation field, water phase saturation field, CO2 gas phase saturation field, CO2 concentration field, permeability field, porosity field and crude oil viscosity field.

5. The frequency domain divide-and-conquer and physically constrained multi-field intelligent prediction method for CO2 oil displacement according to claim 1, characterized in that, The temporal feature extraction and compression of the preprocessed multiphysics 3D temporal data samples are performed through a temporal processing branch to obtain the original temporal compressed feature vector for each physics field. This includes the following steps: The three-dimensional tensor of the preset time step input to each time series processing branch is embedded with position information using a learnable three-dimensional position encoding matrix. The dimension of the encoding matrix is ​​completely consistent with the input tensor. The three-dimensional spatial coordinate information of the reservoir grid I, J, K and the time series information of the time step are synchronously embedded into the input features to obtain the position-encoded three-dimensional time series tensor. The three-dimensional temporal tensor after position encoding is divided into fixed-size blocks along the spatial dimension, resulting in multiple spatial blocks. Each spatial block corresponds to the temporal features of a preset time step. Each block is flattened and mapped to a token vector of a preset dimension through a linear mapping layer. The single physics data of each time step is converted into token vectors of multiple preset dimensions to obtain the token sequence of the preset time step. The token sequence at a preset time step is input into the VIT encoder, which is a multi-layer TransformerBlock stack, to complete the adaptive extraction of temporal-spatial coupling features. The VIT encoder outputs a token sequence of the same length as the input. The token sequence output by the VIT encoder is compressed into a single vector of a set dimension through global average pooling and a linear compression layer, which is used as the temporal compression feature of the physical field. Each temporal processing branch outputs an independent temporal compression feature vector of a set dimension, thus obtaining the original temporal compression feature vector of each physical field.

6. The frequency domain divide-and-conquer and physically constrained multi-field intelligent prediction method for CO2 oil displacement according to claim 1, characterized in that, The original time-series compressed feature vectors are mapped to the frequency domain space, and frequency domain noise filtering is performed, including the following steps: An independent frequency domain mapping branch is constructed for each physical field's temporal compressed feature vector; the temporal compressed feature vector of each physical field with a set dimension is increased in dimension through a linear layer, and then converted into a two-dimensional feature tensor through a reshape operation to obtain the preprocessed two-dimensional feature tensor. The preprocessed two-dimensional feature tensor is input into a multi-layered stacked FNO layer. Each FNO layer includes a Fast Fourier Transform (FFT) module, a frequency domain filtering and fusion module, and an Inverse Fast Fourier Transform (IFFT) module. The FFT module transforms the feature tensor in the time domain to the Fourier frequency domain. In the frequency domain, the frequency domain filtering and fusion module decomposes the original frequency domain features into high-frequency features, low-frequency features, and environmental noise features. A learnable weight matrix is ​​used to perform a linear transformation, converting the nonlinear evolution relationship in the time domain into a linear representation in the frequency domain. Finally, the IFFT module transforms the feature tensor back to the spatial domain, completing one frequency domain feature mapping. In the frequency domain transformation stage of each FNO layer, a learnable Gaussian filter kernel is embedded to smooth the frequency domain features and automatically remove high-frequency noise components in the frequency domain. After multiple FNO layers and noise filtering, the frequency domain feature vector corresponding to the set dimension of each physical field is output. Each physical field outputs an independent frequency domain feature vector, completing the full mapping from time domain features to frequency domain space.

7. The frequency domain divide-and-conquer and physically constrained multi-field intelligent prediction method for CO2 oil displacement according to claim 1, characterized in that, The features of different physical fields in the filtered frequency domain space are divided into high and low frequencies, and the high-frequency abrupt change features and low-frequency background features are dynamically fused to obtain the global frequency domain fused features, including the following steps: For each physical field, the frequency domain feature vector of a given dimension is automatically decomposed into high-frequency feature components and low-frequency feature components through a learnable adaptive frequency domain partitioning threshold. The high-frequency feature components are converted into a token sequence through convolution transformation. Each high-frequency component of a physical field is an independent token. The token sequence is input into the VIT encoder. The VIT encoder automatically captures the correlation between high-frequency features of different physical fields through a self-attention mechanism, completes the cross-physical field deep fusion of key high-frequency features, and outputs the fused global high-frequency feature vector. The low-frequency feature components are converted into a token sequence, with each low-frequency component of a physical field being an independent token. The token sequence is then input into the ViT encoder to complete the cross-physical field deep fusion of the background low-frequency features and output the fused global low-frequency feature vector. The fused global high-frequency feature vector and the global low-frequency feature vector are concatenated along the channel dimension. The concatenated features are then reduced to a global frequency domain fused feature vector of a set dimension by a multilayer perceptron, thus obtaining the global frequency domain fused features.

8. The frequency domain divide-and-conquer and physical constraint-based multi-field intelligent prediction method for CO2 oil displacement according to claim 1, characterized in that, The process of inverting global frequency domain fused features into time domain dynamic incremental features includes the following steps: The frequency domain features are inverted to the time domain features through the fast inverse Fourier transform of the Fourier neural operator, and then the dynamic incremental features of each physical field are generated through nonlinear transformation.

9. The frequency domain divide-and-conquer and physical constraint-based multi-field intelligent prediction method for CO2 oil displacement according to claim 1, characterized in that, The original time-series compressed feature vector is adaptively weighted with the time-domain dynamic incremental feature to obtain the fused feature of each physical field, including the following steps: Construct an independent gated weighting unit for each physical field; The gated weighting unit adaptively weights and fuses the original temporal compressed feature vector and the temporal dynamic incremental feature through the output weights to obtain the fused feature of each physical field.

10. The frequency domain divide-and-conquer and physically constrained multi-field intelligent prediction method for CO2 oil displacement according to claim 1, characterized in that, The physical consistency constraint loss function is determined based on the data fitting loss and the inter-field coupling physical constraint loss.

11. The frequency domain divide-and-conquer and physically constrained multi-field intelligent prediction method for CO2 oil displacement according to claim 10, characterized in that, The interfield coupling physical constraint loss is constructed based on the seepage mechanics mass conservation equation and saturation constraint equation, including saturation constraint, pressure-saturation coupling constraint, and concentration conservation constraint.

12. The frequency domain divide-and-conquer and physically constrained multi-field intelligent prediction method for CO2 oil displacement according to claim 11, characterized in that, The saturation constraint includes ensuring that the sum of the saturations of the oil phase, water phase, and gas phase is always equal to 1 at each grid node, and penalizing saturation predictions that exceed the range of 0-1. The pressure and saturation coupling constraint includes constructing a coupling relationship between pressure gradient and fluid saturation change based on Darcy's law for multiphase seepage, and penalizing non-physical prediction results that violate seepage laws; Concentration conservation constraints include ensuring the total mass of CO2 is conserved within a closed reservoir system and penalizing predicted anomalous changes in CO2 concentration.

13. The frequency domain divide-and-conquer and physically constrained multi-field intelligent prediction method for CO2 oil displacement according to any one of claims 1-12, characterized in that, CO2 enhanced oil recovery engineering parameters include geological parameters, fluid parameters, and operating condition parameters.

14. A frequency-domain divide-and-conquer and physically constrained CO2 enhanced oil recovery multi-field intelligent prediction system, characterized in that, include: The dataset construction module is used to build a numerical simulation model of CO2 flooded reservoirs, extract multiphysics three-dimensional time-series data samples from the CO2 flooded reservoir numerical simulation model and perform preprocessing. The parallel temporal feature extraction module is used to construct an independent parallel temporal processing branch for each physical field. The temporal processing branch extracts and compresses the temporal features of the preprocessed data samples to obtain the original temporal compressed feature vector of each physical field. The feature mapping module is used to map the original time-series compressed feature vectors to the frequency domain space and perform frequency domain noise filtering. The frequency division and fusion module is used to divide the features of different physical fields in the filtered frequency domain space into high and low frequencies, and dynamically fuse the high-frequency abrupt features and low-frequency background features to obtain global frequency domain fusion features. The frequency domain inversion and incremental generation module is used to invert global frequency domain fusion features into time domain dynamic incremental features; The model training module is used to adaptively weight the original time-series compressed feature vector and the time-domain dynamic incremental feature to obtain the fusion feature of each physical field, and construct a multi-physics collaborative prediction model. The fusion feature of each physical field is used as input to output the multi-field feature of the next time step. The training process of the multi-physics collaborative prediction model ensures the physical consistency of the predicted value through the fusion loss function embedded with mechanistic knowledge. The multiphysics prediction module is also used to input the CO2 flooding engineering parameters of the target reservoir into the trained multiphysics collaborative prediction model and output the CO2 flooding multiphysics prediction results.

Images