A data physically dual-driven crack propagation prediction method

Abstract

The application discloses a kind of data physical double drive's crack propagation prediction method, belong to petroleum engineering field, including the following steps: step 1, construct real observation dataset and sampling dataset;Step 2, construct the hybrid architecture of fusion Transformer encoder and graph attention network;Step 3, construct three independent and parallel decoders to map shared features into the geometric size and mechanical parameters of crack;Step 4, based on linear elastic fracture mechanics and material balance principle to establish physical loss function, and combined with data loss function to construct hybrid loss function for model training;Step 5, use the trained model to predict the geometric size and mechanical parameters of the fracture.The application embeds physical priori knowledge into the multi-task deep learning framework, and embeds the physical loss into the loss function, which significantly improves the accuracy and interpretability of the model prediction.

Description

TECHNICAL FIELD

[0001] The application belongs to the field of petroleum engineering, and particularly relates to a data-physical dual-driven fracture propagation prediction method based on a Transformer and a graph attention network. BACKGROUND

[0002] In the field of oil and gas field exploration and development, large-scale fracturing of unconventional reservoirs has become the main way for the development of unconventional oil and gas reservoirs, and the fracture geometry after fracturing directly determines the effect of fracturing operation. Accurate prediction of these parameters is crucial for optimizing fracturing operation design, improving recovery efficiency and reducing development cost. In order to accurately predict the geometric parameters of the fracturing fracture, engineers use mechanism-based numerical simulation methods to predict these fracture parameters. The mechanism-based numerical simulation method is to establish a partial differential equation set according to rock mechanics, fluid mechanics and fracture mechanics, and then solve it by using the finite element method, boundary element method, etc. Although this method has clear physical meaning and strong interpretability, it is complex to calculate, time-consuming and cannot meet the needs of rapid prediction in oilfield site.

[0003] To solve the above problems, a pure data-driven method based on machine learning begins to replace the numerical simulation process for prediction. The pure data-driven method based on machine learning is to use a machine learning model to learn the mapping relationship between input parameters and output parameters from a large amount of historical data. In the training process, a loss function is first defined to calculate the error between the model prediction value and the true observation value, and then an optimizer is used to iteratively optimize the loss function. In each iteration, all learnable parameters in the model are updated using the backpropagation algorithm. After multiple rounds of training, a trained model is obtained which can better represent the nonlinear relationship between input and output. Although this method is fast in calculation, it is very dependent on the quality and quantity of historical data, and the model as a "black box" leads to a lack of interpretability of the prediction results and low accuracy. SUMMARY

[0004] To solve the above problems, the application proposes a data-physical dual-driven fracture propagation prediction method, which introduces a hybrid network architecture combining a Transformer encoder and a graph attention network to capture the complex interaction between features and establish a feature map based on physical priori, and then outputs the geometric parameters of the fracture through three parallel and independent decoders. This model structure can accurately quantify the complex coupling effect between input parameters, improve prediction accuracy, and introduce a hybrid loss function combining data loss and physical loss to ensure that the prediction results can fit the true observation data while following the physical law, solving the problem of poor interpretability of pure data-driven models.

[0005] The technical scheme of the application is as follows:

[0006] A data-physical dual-driven crack propagation prediction method, comprising the following steps:

[0007] Step 1, construct a real observation data set with labels and a physical sampling data set without labels, and perform normalization processing;

[0008] Step 2, construct a hybrid architecture that fuses a Transformer encoder with a graph attention network, perform deep interaction on the input features, and fuse with physical priors to generate physically information-enhanced shared features;

[0009] Step 3, construct three independent and parallel decoders to map the shared features to the geometric size and mechanical parameters of the crack;

[0010] Step 4, establish a physical loss function based on linear elastic fracture mechanics and material balance principle, and combine it with a data loss function to construct a hybrid loss function with physical meaning for model training;

[0011] Step 5, use the trained model to predict the geometric size and mechanical parameters of the fracturing crack.

[0012] Further, the specific process of step 1 is:

[0013] Step 1.1, obtain the reservoir rock geomechanical parameters, fluid parameters, fracturing construction parameters, and fracturing crack geometric parameters interpreted by microseismic, and after processing the missing values and abnormal values, construct a real observation data set with labels , then divide the training set and test set in proportion;

[0014] Step 1.2, determine the minimum and maximum values of each feature based on the numerical range of each input feature in the real observation data set, thereby defining the boundary of the physical parameter space, then use Latin hypercube sampling to generate a large number of sampling points in the multi-dimensional parameter space, and establish a physical sampling data set without labels , use the entire physical sampling data set as the training set;

[0015] Step 1.3, use the geomechanical parameters, fluid parameters, fracturing construction parameters, and fracturing crack geometric parameters in the real observation data set and the physical sampling data set as samples to construct a sample set, and perform normalization processing; specifically:

[0016] Arrange the geomechanical parameters, fluid parameters, and fracturing construction parameters in the real observation data set and the physical sampling data set as samples respectively, and correspondingly obtain two input feature sample sets and ; the fracture geometry parameters in the real observation dataset are arranged as output sample sets ; wherein, is the th sample set in the real observation dataset ; is the th sample set in the real observation dataset ; is the th sample set in the real observation dataset ;

[0017] Each sample set in each sample set is normalized, and the two input feature sample sets after normalization are , ; and the output sample set after normalization is .

[0018] Further, the specific process of step 2 is as follows:

[0019] Step 2.1, the two input feature sample sets after normalization , are classified according to their physical properties into geomechanical parameters, fluid parameters and fracturing construction parameters; a fully connected neural network is used to linearly project the classified parameters, and they are mapped to a high-dimensional vector space to obtain three high-dimensional feature vectors; then the high-dimensional feature vectors are spliced, and a learnable position code is introduced to form a high-dimensional embedding feature with a unified dimension;

[0020] Step 2.2, a Transformer encoder is constructed to capture the complex relationship between features;

[0021] Step 2.3, a graph attention network is constructed to generate an enhanced feature matrix combining physical priors and data-driven features, and the graph attention network includes a first graph attention layer, a Dropout layer, a second graph attention layer and an activation layer;

[0022] Step 2.4, the output of the Transformer encoder and the output of the graph attention network are spliced to obtain shared features enhanced by physical information.

[0023] Further, the specific process of step 2.1 is as follows:

[0024] Step 2.1.1, the two input feature sample sets after normalization , are classified according to their physical properties into three subsets of geomechanical parameters, fluid parameters and fracturing construction parameters; then an independent linear layer is used to project each type of parameter:

[0025] ;

[0026] ;

[0027] ;

[0028] wherein, , , are three subsets of geomechanical parameters, fluid parameters and fracturing operation parameters respectively; , , are the learnable weight matrices corresponding to the three subsets; , , are the high-dimensional feature vectors obtained after linear projection of the three subsets, , , are the high-dimensional embedding dimensions after projection of the three subsets; is the batch size;

[0029] Step 2.1.2, splice the high-dimensional feature vectors obtained by projection to form a preliminary fusion feature, the formula is as follows:

[0030] ;

[0031] wherein, is the preliminary fusion feature, is the total embedding dimension;

[0032] Step 2.1.3, introduce a learnable position encoding to form a high-dimensional embedding feature of a unified dimension; the formula is as follows:

[0033] ;

[0034] wherein, is the high-dimensional embedding feature; is the learnable position encoding of the embedding feature.

[0035] Further, the specific process of step 2.2 is as follows:

[0036] Step 2.2.1, convert the high-dimensional embedding feature into a query matrix , a key matrix and a value matrix through linear transformation; then generate the context vector of the attention head through , , , the formula is as follows:

[0037] ;

[0038] wherein, is the context vector of the i-th attention head; , , , are the query matrix, key matrix and value matrix of the i-th attention head, respectively; is the transpose operation; is the scaling factor; is the softmax function;

[0039] Step 2.2.2, the context vectors of the plurality of attention heads are spliced, and the output of the multi-head attention module is obtained through a linear projection matrix , the formula is:

[0040] ;

[0041] wherein, is a learnable linear projection matrix; is the splicing operation; is the context vector of the i-th attention head; Step 2.2.3, the output of the multi-head attention module is connected in residual with the high-dimensional embedding feature

[0042] and layer normalization operation is performed to obtain the embedding feature :

[0043] ;

[0044] wherein, is the layer normalization operation;

[0045] Step 2.2.4, the nonlinear transformation feature is obtained by performing nonlinear transformation on based on the feedforward neural network:

[0046] ;

[0047] wherein, is the ReLU activation function; , are learnable weights; , are biases;

[0048] Step 2.2.5, the second residual connection and layer normalization operation are performed on and to obtain the output of the first encoder​​ :

[0049] ;

[0050] Step 2.2.6, superimpose two encoder layers, i.e. as Input Step 2.2.1 again until the final output of the encoder part is obtained .

[0051] Further, the specific process of Step 2.3 is as follows:

[0052] Step 2.3.1, input the high-dimensional embedding features obtained in Step 2.1 into the first graph attention layer for linear transformation, attention calculation and multi-head attention output; the specific process is as follows: First, linear transformation is performed:

[0053]

[0054] ;

[0055] wherein, is a learnable linear transformation matrix; is the transformed node feature;

[0056] Then, attention calculation is performed; according to the prior knowledge in the field of hydraulic fracturing, the connection relationship between nodes is defined, an adjacency matrix based on physical prior is constructed, and the attention weight is calculated based on the adjacency matrix. When there is a physical connection between node and node , the attention coefficient is calculated:

[0057] ;

[0058] wherein, is an adjacency matrix;

[0059] ;

[0060] wherein, is an attention weight vector; , are the transformed node features of node and node , is feature splicing of node and node ; is the attention coefficient between node and node ; is a linear rectifier function with a leakage coefficient; ​

[0061] Normalization operation is performed on all neighbor nodes:

[0062] ;

[0063] wherein, is the attention weight between node and node ; is the exponential operation; is the element in the adjacency matrix in the i-th row and the j-th column, the i-th row corresponds to node , and the j-th column corresponds to node ; is the element in the adjacency matrix in the i-th row and the j-th column; is the attention coefficient between node and node ; is the number of neighbor nodes; Then, multi-head attention output is performed; two attention heads are used to calculate and integrate node neighborhood features respectively: ;

[0064] wherein, are the output features of node

[0065] of the 1st and 2nd attention heads respectively; is a nonlinear activation function; are the attention weights between node

[0066] and node in the 1st and 2nd attention heads respectively;

[0067] The obtained is spliced to obtain the final enhanced feature: ; is the enhanced feature of node

[0068] ; Step 2.3.2, input

[0069] into the Dropout layer to perform the Dropout operation:

[0070] ; is the enhanced feature of node

[0071] ; Step 2.3.3, input

[0072] into the Dropout layer to perform the Dropout operation: ​​​

[0073] wherein, is a random mask matrix consistent with the dimension of ; is a deactivation probability; is a corresponding element multiplication; is a node enhancement matrix after the Dropout operation; is a Dropout layer;

[0074] Step 2.3.3, inputting the shared feature into the first graph attention layer to obtain a first graph attention layer output as the input of the second graph attention layer, repeating the operation of the first graph attention layer to obtain a high-order node feature:

[0075] ;

[0076] wherein, is the output of the second graph attention layer; is the second graph attention layer;

[0077] Step 2.3.4, inputting the output of the second graph attention layer into the activation layer to obtain a final enhanced feature matrix:

[0078] ;

[0079] wherein, is the enhanced feature matrix.

[0080] Further, in the step 2.4, the calculation formula of the shared feature is:

[0081] ;

[0082] wherein, is the shared feature.

[0083] Further, the specific process of the step 3 is:

[0084] Step 3.1, constructing a physical information attention layer, inputting the shared feature obtained in the step 2.4 into the physical information attention layer to obtain a shared physical feature ; the specific process is: Step 3.1.1, generating a query matrix

[0085] , a key matrix , and a value matrix :

[0086] ;

[0087] ;

[0088] ​ ;

[0089] wherein, , , are three learnable weight coefficients;

[0090] Step 3.1.2, calculating the similarity of query matrix and key matrix , and then obtaining a set of attention weights after passing through the softmax function:

[0091] ;

[0092] Step 3.1.3, calculating the weighted sum of to obtain the shared physical features :

[0093] ;

[0094] Step 3.2, inputting into three parallel and independent multi-layer perceptron decoders at the same time, respectively obtaining six predicted values: crack area, crack length, crack width, crack height, net pressure, and stress intensity factor; wherein the crack area, crack length, crack width, and crack height are geometric dimensions, and the net pressure and stress intensity factor are mechanical parameters.

[0095] Further, the specific process of the step 4 is:

[0096] Step 4.1, constructing a physical loss based on the residual error of the equation of elastic mechanics , and the specific formula is:

[0097] ;

[0098] wherein, is the mean square error; is the predicted value of the stress intensity factor; is the predicted value of the crack width; is the predicted value of the net pressure; is a proportional coefficient; is the equivalent elastic modulus;

[0099] Step 4.2, constructing a physical loss based on the residual error of the crack propagation criterion :

[0100] ;

[0101] wherein, is the fracture toughness;

[0102] Step 4.3, Constructing physical loss based on residual of flow equation :

[0103] ;

[0104] wherein, is the predicted value of crack length; is the predicted value of crack height; is the displacement; is the time; is the shape factor;

[0105] Step 4.4, The total physical loss is :

[0106] ;

[0107] wherein, , , is the weight coefficient of different loss terms;

[0108] Step 4.5, Calculate data loss :

[0109] ;

[0110] wherein, is the number of samples; is the th predicted output of the th sample; is the th real value of the th sample;

[0111] Step 4.6, Add the data loss and the total physical loss in proportion to obtain the total hybrid loss function :

[0112] ;

[0113] wherein, , is the weight coefficient of different loss;

[0114] Step 4.7, Use the predicted value output by step 3 to calculate the hybrid loss function , and use the Adam optimizer to minimize the hybrid loss, then update all learnable parameters in the model based on the training set by back propagation, so that the prediction result of the model gradually approaches and meets the physical constraints, thus completing the model training.

[0115] Further, the specific process of step 5 is: inputting the samples in the step 1 test set and the adjacency matrix in step 2 into the trained model to predict the crack propagation, including crack area, crack length and crack width.

[0116] The beneficial technical effects brought by the present application: the present application divides the input features into three categories: geomechanical parameters, fluid parameters and fracturing operation parameters, so that the input features interact in their respective subspaces first, capturing the strong correlation between the features, while reducing the noise generated between weakly correlated features, enhancing the generalization of the model, so that the model can accurately predict the geometry of the crack even when it faces input data that has never been seen before; the present application proposes a hybrid model architecture based on Transformer and graph attention network, and uses multiple independent and parallel decoders to output the geometry and physical parameters of the crack, this model structure can capture the possible nonlinear relationship between all features, and use the graph attention network to make the interaction between features have a clear physical meaning, so that the model can improve the prediction accuracy while ensuring that the learning process meets the physical laws; the present application embeds multiple physical loss terms into the loss function to construct a hybrid loss function with physical meaning, and constantly updates the model parameters through gradient backpropagation during the training process, so that the laws learned by the model have more physical meaning, greatly improving the interpretability of the model. BRIEF DESCRIPTION OF DRAWINGS

[0117] Figure 1 is the overall design flowchart of the data-physical dual-driven crack propagation prediction method of the present application.

[0118] Figure 2 is the model structure diagram of the crack propagation prediction of the present application.

[0119] Figure 3 is the loss diagram in the model training process of the present application.

[0120] Figure 4 is the comparison diagram of the real value and the predicted value of the crack area in the embodiment of the present application.

[0121] Figure 5 is the comparison diagram of the real value and the predicted value of the crack length in the embodiment of the present application.

[0122] Figure 6 is the comparison diagram of the real value and the predicted value of the crack width in the embodiment of the present application.

[0123] Figure 7 is the loss comparison diagram in the training process of the present application and the pure data-driven method.

[0124] Figure 8 is the R 2 comparison diagram of the present application and the pure data-driven method. DETAILED DESCRIPTION

[0125] The application will be described in further detail below in conjunction with the accompanying drawings and specific embodiments:

[0126] The main process of the application is: constructing a sample data set; constructing a hybrid architecture that fuses Transformer and graph neural network, deeply interacting input features and fusing with physical priors to generate physically information-enhanced shared features; constructing three parallel structure decoders to map the shared features to the geometric size and mechanical properties of the cracks; based on the linear elastic fracture mechanics and the principle of material balance, constructing a physical residual term, combining it with the mean square error loss of the data-driven prediction result, establishing a comprehensive loss function with physical meaning and performing model training; after training, the geometric shape and mechanical parameters of the cracks can be quickly predicted through the geomechanical parameters, fluid parameters and construction parameters. As shown in the figure, a data-physical dual-driven crack propagation prediction method specifically comprises the following steps: Figure 1

[0127] Step 1, construct a real observation data set with labels and a sampling data set without labels, and perform normalization processing; the specific process is:

[0128] Step 1.1, obtain the geomechanical parameters of the reservoir rock, fluid parameters, fracturing construction parameters and fracturing crack geometric parameters interpreted by microseismic, and construct a real observation data set with labels after processing the missing values and abnormal values , and then divide the training set and the test set in the ratio of 8.5:1.5;

[0129] Step 1.2, based on the numerical range of each input feature in the real observation data set, determine the minimum and maximum values of each feature to define the boundary of the physical parameter space, then use Latin Hypercube Sampling (LHS) to generate a large number of sampling points in the multi-dimensional parameter space to establish a physical sampling data set without labels , the entire physical sampling data set is not involved in testing as the training set;

[0130] Step 1.3, the geomechanical parameters, fluid parameters, fracturing construction parameters and fracturing crack geometric parameters in the real observation data set and the physical sampling data set are constructed into a sample set and normalized;

[0131] The geomechanical parameters, fluid parameters and fracturing construction parameters in the real observation data set and the physical sampling data set are arranged as samples respectively, and two input feature sample sets are obtained correspondingly​ and ; the fracture geometry parameters in the real observation data set are arranged as output sample sets ; wherein, is the group sample in ; is the group sample in ; is the group sample in ;

[0132] Each group of samples in each sample set is normalized, and the normalization method adopts the following formula:

[0133] ;

[0134] wherein, is the normalized data value; is the group sample in the input feature sample set; is the minimum value in the original data; is the maximum value in the original data;

[0135] The two normalized input feature sample sets are , ; and the normalized output sample set is .

[0136] Step 2, a hybrid architecture combining the Transformer encoder and the graph attention network is constructed to perform deep interaction on the input features and fuse them with the physical prior, thereby generating physical information enhanced shared features; the hybrid architecture first classifies the normalized sample data, then performs feature embedding and fusion, acquires the complex nonlinear relationship between the features through the Transformer encoder, constructs a feature relationship graph based on the physical prior through the graph attention network according to the adjacency matrix, and then fuses the outputs of the Transformer encoder and the graph attention network, as shown in the model structure Figure 2 , and the specific working process is as follows:

[0137] Step 2.1, the two normalized input feature sample sets , The parameters are classified into geomechanical parameters, fluid parameters and fracturing operation parameters according to their physical properties; the classified parameters are linearly projected by using a fully connected neural network to map them to high-dimensional vector spaces to obtain three high-dimensional feature vectors; then the high-dimensional feature vectors are spliced and a learnable position code is introduced to form a high-dimensional embedding feature with a unified dimension, and the specific process is as follows:

[0138] Step 2.1.1, splice the normalized input feature sample set (including , ) According to their physical properties, they are classified into three subsets of geomechanical parameters, fluid parameters and fracturing operation parameters; then an independent linear layer (MLP) is used for each type of parameter for projection:

[0139] ;

[0140] ;

[0141] ;

[0142] Among them, , , are the three subsets of geomechanical parameters, fluid parameters and fracturing operation parameters; , , are the learnable weight matrices corresponding to the three subsets; , , are the high-dimensional feature vectors obtained after linear projection of the three subsets, , , are the high-dimensional embedding dimensions of the three subsets after projection; is the shape factor, which is the sum of the dimensions of each type of feature; is the batch size;

[0143] Step 2.1.2, splice the high-dimensional feature vectors obtained by projection to form a preliminary fusion feature, and the formula is as follows:

[0144] ;

[0145] Among them, is the preliminary fusion feature, is the total embedding dimension;

[0146] Step 2.1.3, introduce a learnable position code to form a high-dimensional embedding feature with a unified dimension; the formula is as follows:

[0147] ;

[0148] wherein, is the high-dimensional embedding feature; is the learnable position encoding of the embedding feature;

[0149] Step 2.2, the Transformer encoder is constructed to capture the complex relationship between features; the specific process is as follows:

[0150] Step 2.2.1, the high-dimensional embedding feature is converted into a query matrix , a key matrix and a value matrix through linear transformation; then the context vector of the attention head is generated through , , The formula is as follows (for the th attention head):

[0151] ;

[0152] wherein, is the context vector of the th attention head; , , are the query matrix, key matrix and value matrix of the th attention head, respectively; is the transpose operation; is the scaling factor; is the softmax function, which is used to normalize the weight matrix;

[0153] Step 2.2.2, the context vectors of multiple attention heads are spliced and the output of the multi-head attention module is obtained through a linear projection matrix The formula is as follows:

[0154] ;

[0155] wherein, is a learnable linear projection matrix; is the splicing operation; is the context vector of the th attention head;

[0156] Step 2.2.3, the output of the multi-head attention module is connected in residual with the high-dimensional embedding feature and the layer normalization operation is performed to obtain the embedding feature :

[0157] ;

[0158] wherein, is a layer normalization operation;

[0159] Step 2.2.4, performing a nonlinear transformation on the to obtain a nonlinearly transformed feature :

[0160] ;

[0161] wherein, is a ReLU activation function; , is a learnable weight; , is a bias;

[0162] Step 2.2.5, performing a second residual connection and a layer normalization operation on and to obtain the output of the first encoder :

[0163] ;

[0164] wherein, is a layer normalization operation, is a nonlinearly transformed feature; is a high-dimensional embedding feature;

[0165] Step 2.2.6, stacking two encoder layers, i.e., taking as again inputting Step 2.2.1 until the final output of the encoder part is obtained .

[0166] Step 2.3, constructing a graph attention network (GAT), including a first graph attention layer, a Dropout layer, a second graph attention layer, and an activation layer; the graph attention network defines the physical connection relationship between feature nodes according to the adjacency matrix constructed in Step 1.4, then adaptively calculates the weight between nodes based on the attention mechanism, and finally performs weighted summation on the features of neighbor nodes to generate an enhanced feature matrix that combines physical prior and data-driven ; the specific process is as follows:

[0167] Step 2.3.1, inputting the high-dimensional embedding feature obtained in Step 2.1 into the first graph attention layer for linear transformation, attention calculation, and multi-head attention output; the specific process is as follows:

[0168] Step 2.3.1.1, first, linear transformation is performed:

[0169] ;

[0170] wherein, is a learnable linear transformation matrix; is the transformed node feature;

[0171] Step 2.3.1.2, perform attention calculation, realize adaptive weight distribution based on physical adjacency; define the connection relationship between nodes according to the prior knowledge in the field of hydraulic fracturing, construct the adjacency matrix based on physical prior, and calculate the attention weight based on the adjacency matrix; when the node and the node have physical connection (i.e. ), the attention coefficient is calculated:

[0172] ;

[0173] wherein, is the adjacency matrix; when the node and the node have strong interaction, physical connection is established, i.e. ;

[0174] ;

[0175] wherein, is the attention weight vector; is the transpose operation; , respectively the transformed node feature of the node and the node , is to splice the node and the node ; is the attention coefficient between the node and the node ; is a linear rectifier function with leakage coefficient;

[0176] Then, normalization operation is performed on all neighbor nodes:

[0177] ;

[0178] wherein, is the attention weight between the node and the node ; is the exponential operation; is the element in the th row and the th column of the adjacency matrix, the row corresponding node , the column corresponding node ; is an element in the adjacency matrix in the row and the column; is an attention coefficient between the node and the node ; is the number of neighbor nodes;

[0179] Step 2.3.1.3, multi-head attention output; the present application uses two attention heads, and first calculates and integrates node neighborhood features respectively to obtain :

[0180] ;

[0181] ;

[0182] wherein, are output features of the node in the first and second attention heads respectively; is a nonlinear activation function; are attention weights between the node and the node in the first and second attention heads respectively;

[0183] Then, the obtained is spliced to obtain the final enhanced feature:

[0184] ;

[0185] is the enhanced feature of the node ;

[0186] Step 2.3.2, in order to prevent overfitting after multi-head feature splicing, the spliced output is input into a Dropout layer for Dropout operation:

[0187] ;

[0188] wherein, is a random mask matrix with the same dimension as ; is an inactivation probability; is a corresponding element multiplication; is a node enhanced matrix after the Dropout operation; is a Dropout layer;

[0189] Step 2.3.3, repeat the operation of the first graph attention layer to obtain higher-order node features: As the input of the second graph attention layer, repeat the operation of the first graph attention layer to obtain higher-order node features:

[0190] ;

[0191] wherein, is the output of the second graph attention layer; is the second graph attention layer.

[0192] Step 2.3.4, finally input the output of the second graph attention layer into the activation layer, apply a nonlinear activation function, and obtain the final enhanced feature matrix:

[0193] ;

[0194] wherein, is the enhanced feature matrix.

[0195] Step 2.4, the final output of the Transformer encoder is concatenated with the enhanced feature matrix output by the graph attention network to obtain shared features enhanced by physical information :

[0196] ;

[0197] Step 3, construct a physical information attention layer, use a fully connected neural network to construct three independent and parallel decoders to map the shared features to the geometric dimensions and mechanical parameters of the crack, the three decoders are respectively a geometric parameter prediction decoder with label (predicting crack area , crack length , crack width ), a geometric parameter prediction decoder without label (predicting crack height ), and a mechanical state prediction decoder (predicting net pressure , stress intensity factor ), the geometric dimensions include crack area, crack length, crack width, and crack height, and the mechanical parameters include net pressure and stress intensity factor; input the shared features obtained in step 2 into the physical information attention layer to obtain shared physical features , and then input into the three independent decoders respectively to output the respective prediction results, using three independent decoders can make the model select more suitable weights and biases for its own output parameters, and the model structure is shown in Figure 2 , and the specific process is as follows:

[0198] Step 3.1: Construct a physical information attention layer, which incorporates the shared features obtained in Step 2.4. Input physical information into the attention layer to obtain shared physical features. The specific process is as follows:

[0199] Step 3.1.1: Generate the query matrix Key matrix Sum matrix :

[0200] ;

[0201] ;

[0202] ;

[0203] in, , , These are three learnable weight coefficients;

[0204] Step 3.1.2: Calculate the query matrix Bond matrix The similarity is then used to obtain a set of attention weights after passing through the softmax function. :

[0205] ;

[0206] in, This is a transpose operation; This is the scaling factor; It is the softmax function;

[0207] Step 3.1.3: Calculate the weights right The weighted sum yields the shared physical characteristics. :

[0208] ;

[0209] Step 3.2, Simultaneously, the inputs are fed into three parallel and independent multilayer perceptron (MLP) decoders, which respectively obtain six target physical quantities that need to be predicted. , , , , , ) is used to calculate the mixing loss, where The area of ​​the crack. The length of the crack. The width of the crack. The height of the crack. Net pressure, This is the stress intensity factor.

[0210] Step 4: Establish a physical loss function based on linear elastic fracture mechanics and the principle of mass balance, and combine it with the data loss function to construct a physically meaningful hybrid loss function. Use the predicted physical quantities output in Step 3 to calculate the hybrid loss, and use the Adam optimizer to minimize the hybrid loss. Then, perform backpropagation to update all learnable parameters in the model, so that the model's prediction results gradually approach and satisfy the physical constraints, thereby completing the model training. The specific process is as follows:

[0211] Step 4.1: Construct the residuals of the elasticity equations based on linear elastic fracture mechanics to calculate physical losses. The specific process is as follows:

[0212] Step 4.1.1: For a medium subjected to uniform internal pressure in an infinitely elastic medium For a two-dimensional (plane strain) crack, the width distribution of the crack surface is calculated using the residuals of the elasticity equation (for a length of...). (cracks)

[0213] ;

[0214] in, The width of the crack; The length of the crack is half its length; It is the length coordinate along the crack, from the center ( ) to the tip ( ); The equivalent elastic modulus is given by the formula:

[0215] ;

[0216] in, It is the elastic modulus; It is Poisson's ratio;

[0217] but maximum value for:

[0218] ;

[0219] Calculate the stress intensity factor (For cracks subjected to internal pressure):

[0220] ;

[0221] United maximum value Formulas and stress intensity factors Formula elimination , we get , , The relationship between them is:

[0222] ;

[0223] Step 4.1.2, the model is derived to a more complex three-dimensional crack, although the analytical solution is different, but the dimension is universal, must be proportional to , only the proportionality coefficient will change with the crack geometry, so the three-dimensional crack model is:

[0224] ;

[0225] Step 4.1.3, according to the above elastic mechanics equation residual, construct the physical loss based on the elastic mechanics equation residual , the specific formula is:

[0226] ) ;

[0227] Where, is the mean square error; is the predicted value of the stress intensity factor; is the predicted value of the crack width; is the predicted value of the net pressure;

[0228] Step 4.2, based on the linear elastic fracture mechanics to construct crack propagation criterion residual, to calculate the physical loss based on the crack propagation criterion residual , the specific process is:

[0229] Step 4.2.1, the propagation criterion: when the stress intensity factor at the crack tip reaches the fracture toughness of the material (i.e.

[0230] Step 4.2.2, in this model, it is assumed that crack propagation occurs, thereby introducing the physical loss based on the flow equation residual This loss ensures that the crack is not randomly expanded, but is strictly limited by the mechanical properties of the rock:

[0231] ;

[0232] Step 4.3, based on the material balance principle to construct the flow equation residual, to calculate the physical loss based on the flow equation residual , the specific process is:

[0233] Step 4.3.1, according to the material balance relationship between the injected fluid volume, the total fracture volume and the filtration volume, the following formula is obtained:

[0234] ;

[0235] Where, is the injected fluid volume; is the total fracture volume; is the filtration volume;

[0236] ;

[0237] ;

[0238] Where, is the displacement; is the time; is the shape factor;

[0239] In order to simplify the model, it can be assumed that is 0, and the simplified relationship is obtained:

[0240] ;

[0241] Step 4.3.2, according to the simplified model, the physical loss based on the flow equation residual is calculated , which ensures that the fracture size predicted by the network meets the basic law of conservation of mass:

[0242] ;

[0243] Where, is the predicted value of the fracture length; is the predicted value of the fracture height;

[0244] Step 4.4, according to the above formula, the total physical loss is :

[0245] ;

[0246] Where, , , is the weight coefficient of different loss terms;

[0247] Step 4.5, calculate the data loss :

[0248] ;

[0249] Where, is the number of samples; the first predicted output for the first sample; the first predicted output for the first sample; the first predicted output for the first sample; the first real value for the first sample; Step 4.6, add the data loss and the total physical loss in a certain proportion to obtain a total hybrid loss function

[0250] Step 4.6, add the data loss and the total physical loss in a certain proportion to obtain a total hybrid loss function

[0251]

[0252] wherein, are weight coefficients of different losses;

[0253] Step 4.7, using the predicted physical quantity (P) output in step 3 to calculate the hybrid loss function (L) and using the Adam optimizer to minimize the hybrid loss, then updating all learnable parameters in the model based on the training set, so that the prediction result of the model gradually approaches and meets the physical constraints, thereby completing the model training. Step 4.7, using the predicted physical quantity (P) output in step 3 to calculate the hybrid loss function (L) and using the Adam optimizer to minimize the hybrid loss, then updating all learnable parameters in the model based on the training set, so that the prediction result of the model gradually approaches and meets the physical constraints, thereby completing the model training. Step 5, using the trained model to predict the geometric size and mechanical parameters of the fracture; inputting the samples in step 1 test set and the adjacency matrix in step 2 into the trained model can quickly predict the fracture propagation, including the fracture area, the fracture length and the fracture width.

[0254] In order to prove the feasibility and superiority of the present application, the following examples are given.

[0255] This embodiment verifies the proposed method and performs subsequent analysis and discussion, and then performs the fracture propagation prediction based on the data and physical double driving according to the following steps. First, collect the real observed data (including feature data and label data) on site, since the real data on site often has missing values and abnormal values, so it is necessary to detect and process the missing values and abnormal values, and finally build a real observed data set with 1000 samples and 4000 data

[0256] Through Latin hypercube sampling, simultaneously perturb within the maximum and minimum value range of each feature to obtain a physical sampling data set (no label data) with 1000 samples, and table 1 is the input and output parameter table of the model. Table 1 Model input and output parameter table

[0257] Table 1 Model input and output parameter table ​​​​​​​​​​

[0258] .

[0259] Then according to the process of step 1, the real observation data set is divided into training set and test set in the ratio of 8.5:1.5, and the real observation data set and the physical sampling data set are normalized. According to the process of step 2, a hybrid architecture based on the Transformer encoder and the graph attention network is constructed, and then the normalized samples are divided into three categories according to their attributes and input into the hybrid model based on the Transformer encoder and the graph attention network. The features are mapped to a high-dimensional space, and then after feature interaction and fusion, a shared feature is obtained. In step 2, the adjacency matrix between the features is manually constructed according to the physical prior knowledge, and the adjacency matrix is as follows:

[0260] ;

[0261] According to the process of step 3, a physical information attention layer and three independent and parallel decoders are constructed, and the shared feature is input into the physical information attention layer to output a shared physical feature with physical meaning , , then the , , and are input into the three parallel decoders, and the high-dimensional data is mapped to six prediction outputs of crack area, crack length, crack width, crack height, net pressure and stress intensity factor by the decoder. According to the process of step 4, a hybrid loss function is constructed, including physical loss and data loss based on the residual error of elastic mechanics equation, crack propagation criterion residual error and flow equation residual error , , , , , , , , , , , , , , , , , , , , ; in this embodiment, , , ; using the Adam optimizer to minimize the hybrid loss, and then performing back propagation to update all the learnable parameters in the model, so that the prediction results of the model gradually approach and meet the physical constraints, thereby completing the model training; in the model training process, the batch size is set to 64, the learning rate is 0.001, the training round is 30 rounds, and the random seed is set to 42. According to the process of step 5, the feature can be input into the trained model to quickly predict the area, length and width of the crack.

[0262] Figure 3 The loss graph of the model in the training process is shown in the figure, and it can be seen from the figure that the loss of the model in the training process decreases rapidly and converges to a very small value, indicating that the model well fits the data in the training set and can learn the rules therein.

[0263] Figure 4 The comparison graph of the predicted value and the true value of the crack area is shown in the figure, and it can be seen from the figure that the overall trend of the true value and the predicted value is consistent, and most of the curves can be overlapped, indicating that the model can correctly capture the overall law of the crack area change and the prediction effect is good.

[0264] Figure 5 The comparison graph of the predicted value and the true value of the crack length is shown in the figure, and it can be seen from the figure that the overall trend of the true value and the predicted value is consistent, and the peak and valley positions of the two curves basically correspond, indicating that the model can correctly capture the overall law of the crack height change and the prediction effect is good.

[0265] Figure 6 The comparison graph of the predicted value and the true value of the crack width is shown in the figure, and it can be seen from the figure that the overall trend of the true value and the predicted value is consistent, and basically coincides, indicating that the model can correctly capture the overall law of the crack width change and the prediction effect is good.

[0266] The present application will be compared with the pure data-driven machine learning prediction method, Figure 7 The loss comparison graph of the pure data-driven method and the data-physical double-driven method proposed in the present application in the training process can be seen that the method proposed in the present application converges faster and has smaller training loss.

[0267] Figure 8 The R 2 of the pure data-driven method and the data-physical double-driven method proposed in the present application is shown in the figure, and it can be seen from the figure that the R 2 of the method proposed in the present application is all above 0.9, and is all greater than the R 2 value of the pure data-driven method.

[0268] Of course, the above description is not a limitation of the present application, and the present application is not limited to the above examples. Changes, modifications, additions or substitutions made by those skilled in the art within the spirit and scope of the present application should also be included in the protection scope of the present application.

Claims

1. A data physically dual-driven crack propagation prediction method, characterized by, Comprising the following steps: Step 1, construct a real observation data set with labels and a physical sampling data set without labels, and perform normalization processing; Step 2, build a hybrid architecture that fuses a Transformer encoder with a graph attention network, perform deep interaction on input features, and fuse with physical priors to generate physically information-enhanced shared features; the specific process is as follows: Step 2.1, normalizing two input feature sample sets , According to their physical properties, the parameters are classified into geomechanical parameters, fluid parameters and fracturing operation parameters; linear projection is performed on the classified parameters by using a fully connected neural network, so as to map the parameters to a high-dimensional vector space and obtain three high-dimensional feature vectors; then, the high-dimensional feature vectors are spliced, and a learnable position code is introduced, so as to form a high-dimensional embedding feature with a unified dimension; Step 2.2, build a Transformer encoder to capture complex relationships between features; Step 2.3, build a graph attention network to generate an enhanced feature matrix that combines physical priors and data-driven features, the graph attention network includes a first graph attention layer, a Dropout layer, a second graph attention layer, and an activation layer; Step 2.4, concatenate the output of the Transformer encoder with the output of the graph attention network to obtain physically information-enhanced shared features; Step 3, build three independent and parallel decoders to map shared features to the geometric size and mechanical parameters of the fracture; Step 4, establish a physical loss function based on linear elastic fracture mechanics and material balance principles, and combine it with a data loss function to build a hybrid loss function with physical meaning for model training; Step 5, use the trained model to predict the geometric size and mechanical parameters of the fracture.

2. The data-physical dual-driven crack propagation prediction method of claim 1, wherein, The specific process of step 1 is as follows: Step 1.1, obtain reservoir rock geomechanics parameters, fluid parameters, fracturing construction parameters, and fracturing fracture geometry parameters interpreted by microseismic, and construct a real observation data set with labels after processing missing values and abnormal values Then the training set and the test set are divided in proportion; Step 1.2, based on the numerical range of each input feature in the real observation dataset, determine the minimum and maximum value of each feature to define the boundary of the physical parameter space, and then use Latin hypercube sampling to generate a large number of sampling points in the multi-dimensional parameter space to establish a no-label physical sampling dataset The entire physical sampling dataset is used as the training set; Step 1.3: Transfer the real observation dataset and physical sampling dataset Geomechanical parameters, fluid parameters, fracturing operation parameters, and fracturing fracture geometric parameters were used as samples to construct a sample set, which was then normalized; specifically: geomechanical parameters, fluid parameters and fracturing operation parameters in the real observation dataset and the physical sampling dataset are respectively arranged as samples, and two input feature sample sets and are obtained correspondingly; the fracturing fracture geometry parameters in the real observation dataset are arranged as samples to obtain an output sample set ; wherein, is the th sample in the real observation dataset ; is the th sample in the physical sampling dataset ; is the th sample in the fracture geometry parameter dataset ; The two input feature sample sets after normalization are respectively , ; and the output sample set after normalization is .

3. The data-physical dual-driven crack propagation prediction method of claim 2, wherein, The specific process of step 2.1 is as follows: Step 2.1.1, normalizing two input feature sample sets , , classified into three subsets of geomechanical parameters, fluid parameters and fracturing operation parameters according to their physical attributes; then project each type of parameters using independent linear layers respectively: ； ； ； wherein, , , are three subsets of geomechanical parameters, fluid parameters and fracturing operation parameters, respectively; , , are the learnable weight matrices corresponding to the three subsets; , , are the high-dimensional feature vectors obtained by linear projection of the three subsets, , , are the high-dimensional embedding dimensions after projection of the three subsets; is the batch size; Step 2.1.2, concatenate the high-dimensional feature vectors obtained by projection to form preliminary fusion features, the formula is as follows: ； wherein, is the preliminary fused feature, is the total embedding dimension; Step 2.1.3, introduce learnable position encoding to form high-dimensional embedding features of uniform dimension; the formula is as follows: ； wherein, are high-dimensional embedding features; are learnable positional encodings for the embedding features.

4. The data-physical dual-driven crack propagation prediction method of claim 3, wherein, The specific process of step 2.2 is as follows: Step 2.2.1, high-dimensional embedding features converted into a query matrix by a linear transformation , a key matrix , and a value matrix ; then generate the context vector of the attention head by , , The formula is as follows: ； wherein, is a context vector for the th attention head; , , are a query matrix, a key matrix, and a value matrix for the th attention head, respectively; is a transpose operation; is a scaling factor; is a softmax function; Step 2.2.2, concatenating the context vectors of the plurality of attention heads and obtaining the output of the multi-head attention module by a linear projection matrix , the formula is: ； wherein, is a learnable linear projection matrix; is a concatenation operation; is a first context vector for the h-th attention head; Step 2.2.3, concatenating the output of the multi-headed attention module with the high-dimensional embedding feature performing a residual connection and performing a layer normalization operation to obtain the embedding feature : ； wherein is a layer normalization operation; Step 2.2.

4. Non-linear transformation of the features based on a feed-forward neural network to obtain non-linearly transformed features : ； wherein, is a ReLU activation function; , is a learnable weight; , is a bias; Step 2.2.5, on and performing a second residual connection and layer normalization operation to obtain an output of the first encoder : ； Step 2.2.6, superimpose the two encoder layers, i.e. add When acting as Again, input step 2.2.1 until the encoder part final output is obtained .

5. The data-physical dual-driven crack propagation prediction method of claim 4, wherein, The specific process of step 2.3 is as follows: Step 2.3.1, high-dimensional embedding features obtained in step 2.1 The first graph attention layer is input for linear transformation, attention calculation and multi-head attention output; the specific process is as follows: First, perform linear transformation: ； wherein, is a learnable linear transformation matrix; is the transformed node feature; Then, attention calculation is performed; the connection relationship between nodes is defined according to prior knowledge in the field of hydraulic fracturing, an adjacency matrix based on physical prior is constructed, and the attention weight is calculated based on the adjacency matrix, when there is a physical connection between nodes and nodes , the attention coefficient thereof is calculated: ； wherein is an adjacency matrix; ； wherein, is an attention weight vector; , respectively node with node transformed node features, is a feature concatenation of node with node ; is an attention coefficient between node and node ; is a linear rectifier function with a leakage coefficient; Perform normalization operation on all neighbor nodes: ； wherein, is a node is an attention weight between a node and a node is an exponential operation; is an element in an adjacency matrix at row and column , where row corresponds to a node and column corresponds to a node ; is an element in an adjacency matrix at row and column ; is an attention coefficient between a node and a node ; is a number of neighbor nodes; Then, perform multi-head attention output; use two attention heads to calculate and integrate node neighborhood features: ； ； wherein, are output features of the 1st and 2nd attention head nodes, respectively; are output features of the 1st and 2nd attention head nodes, respectively; is a nonlinear activation function; are attention weights between the 1st and 2nd attention head nodes, respectively; are attention weights between the 1st and 2nd attention head nodes, respectively; are attention weights between the 1st and 2nd attention head nodes, respectively; The resulting The final augmented features are obtained by stitching ； to the node enhanced features; Step 2.3.2, to Dropout operation with an input Dropout layer: ； wherein, is a random mask matrix consistent with the dimension of the input data; is a random mask matrix consistent with the dimension of the input data; is a deactivation probability; is a corresponding element multiplication; is a node enhancement matrix after the Dropout operation; is a Dropout layer; Step 2.3.3, to As the input of the second graph attention layer, the operation of the first graph attention layer is repeated to obtain the high-order node features: ； wherein, is an output of the second graph attention layer; is the second graph attention layer; Step 2.3.4, input the output of the second graph attention layer into the activation layer to obtain the final enhanced feature matrix: ； wherein, is an enhanced feature matrix.

6. The data-physical dual-driven crack propagation prediction method of claim 5, wherein, In step 2.4, the calculation formula of shared features is as follows: ； wherein are shared features.

7. The data-physical dual-driven crack propagation prediction method of claim 6, wherein, The specific process of step 3 is as follows: Step 3.1, constructing the physical information attention layer, inputting the shared features obtained in step 2.4 into the physical information attention layer to obtain shared physical features ; the specific process is as follows: ; the specific process is as follows: Step 3.1.

1. Generating a query matrix , a key matrix , and a value matrix : ； ； ； wherein, , , are three learnable weight coefficients; Step 3.1.2, compute similarity of query matrix and key matrix and then get a set of attention weights after softmax function : ； Step 3.1.3, calculation of a weighted sum of shared physical features : ； Step 3.2, Simultaneously, the data is input into three parallel and independent multilayer perceptron decoders, yielding six predicted values: crack area, crack length, crack width, crack height, net pressure, and stress intensity factor. Among these, crack area, crack length, crack width, and crack height are geometric dimensions, while net pressure and stress intensity factor are mechanical parameters.

8. The data-physical dual-driven crack propagation prediction method of claim 1, wherein, The specific process of step 4 is as follows: Step 4.1, Constructing physical loss based on the residual of the equation of elasticity The specific formula is: ) ； wherein, is the mean square error; is the predicted value of the stress intensity factor; is the predicted value of the crack width; is the predicted value of the net pressure; is the proportionality coefficient; is the equivalent elastic modulus; Step 4.2, Constructing a physical loss based on crack propagation criterion residual : ； wherein, Gc is the fracture toughness; Step 4.3, Constructing the physical loss based on the flow equation residual : ； wherein, is a predicted value of crack length; is a predicted value of crack height; is a displacement; is a time; is a shape factor; Step 4.4, total physical loss is : ； wherein, , , are weight coefficients for different loss terms; Step 4.5, Calculating Data Loss : ； wherein, is the number of samples; is the first predicted output for the th sample; is the first true value for the th sample; is the first predicted output for the th sample; Step 4.6, add the data loss and the total physical loss proportionally to get the total hybrid loss function : ； wherein , are weight coefficients for different losses; Step 4.7, use the predicted value of step 3 output to calculate the hybrid loss function And use the Adam optimizer to minimize the hybrid loss, and then update all learnable parameters in the model based on the training set, so that the prediction result of the model gradually approaches and meets the physical constraints, thereby completing the model training.

9. The data-physical dual-driven crack propagation prediction method of claim 8, wherein, The specific process of step 5 is as follows: input the samples in step 1 test set and the adjacency matrix in step 2 into the trained model to predict the fracture propagation, including fracture area, fracture length, and fracture width.

Images