Numerical solution of wellbore multiphase flow model and gas-liquid distribution state inversion method and system

Abstract

The present application relates to a wellbore multiphase flow model numerical solution and gas-liquid distribution state inversion method and system, belonging to the technical field of petroleum engineering, comprising: step 1: constructing and training the physical information neural network of drilling wellbore multiphase flow dynamic simulation and overflow gas distribution state inversion; determining the input and output of the physical information neural network; determining the loss function of the physical information neural network; training the physical information neural network; step 2: designing an adaptive optimization algorithm to optimize the final solution accuracy and convergence speed of the physical information neural network, obtaining an adaptive physical information neural network; designing an adaptive activation function; designing a residual-based adaptive sampling mechanism; step 3: realizing the numerical solution of the wellbore multiphase flow model and the gas-liquid distribution state inversion based on the adaptive physical information neural network. The present application effectively solves the problem that the traditional numerical method usually requires high-precision grid division and a large amount of computing resources.

Description

Technical Field

[0001] This invention belongs to the field of petroleum engineering technology, specifically relating to a method and system for numerical solution of multiphase flow model in wellbore and inversion of gas-liquid distribution state. Background Technology

[0002] The annular air-liquid distribution characteristics after gas intrusion are a major basis for subsequent well control design and construction. Although traditional numerical methods can describe multiphase flow in the wellbore, they usually require high-precision mesh generation and a large amount of computational resources. Furthermore, the complex downhole conditions and limited observational data make it difficult for traditional numerical methods to achieve high-precision solutions within a reasonable computational cost. With the increasing demand for energy security and intelligent drilling, developing efficient, stable, and physically consistent wellbore multiphase flow inversion methods has become an inevitable trend in the current industry's technological development.

[0003] Currently, numerical solutions for multiphase flow in wellbores mainly employ traditional numerical methods such as finite difference, finite element, or volumetric methods. These methods can accurately describe the characteristics of gas-liquid flow, but they require extremely high mesh accuracy and computational resources, resulting in high computational costs, slow convergence speeds, and poor stability under strongly nonlinear conditions. While data-driven deep learning methods can accelerate computation, they lack physical constraints, are prone to prediction bias and physical inconsistencies, and are difficult to apply in complex conditions. Existing technologies generally suffer from low computational efficiency, poor generalization, and insufficient physical consistency. Summary of the Invention

[0004] To address the shortcomings of existing technologies, this invention proposes a method and system for numerical solution of wellbore multiphase flow models and inversion of gas-liquid distribution state driven by Physics-Informed Neural Network (PINN).

[0005] To address this, this invention proposes a numerical solution and gas-liquid distribution state inversion method for wellbore multiphase flow models based on a Physics-Informed Neural Network (PINN). This method embeds the wellbore multiphase flow control equations as physical constraints into the neural network, directly solving the residual terms of the partial differential equations using automatic differentiation techniques, achieving high-precision prediction of wellbore pressure and gas-liquid integrals. Simultaneously, an adaptive activation function and a residual-based adaptive sampling mechanism are designed to dynamically adjust the training focus in high-gradient regions, thereby improving the model's accuracy and convergence speed in critical areas. This method effectively reduces the grid dependency and computational overhead of traditional numerical methods, enabling high-precision physical field reconstruction under limited data conditions, providing efficient and reliable technical support for intelligent well control, well control operations, and wellbore safety prediction.

[0006] The technical solution of this invention is as follows:

[0007] Numerical solution methods for multiphase flow models in wellbores and methods for inverting gas-liquid distribution states include:

[0008] Step 1: Construct and train a physical information neural network for dynamic simulation of multiphase flow in drilling wellbore and inversion of overflow gas distribution state;

[0009] Determine the input and output of the physical information neural network;

[0010] Determine the loss function of the physical information neural network;

[0011] Training a physical information neural network;

[0012] Step 2: Design an adaptive optimization algorithm to optimize the final solution accuracy and convergence speed of the physical information neural network, and obtain the adaptive physical information neural network;

[0013] Design an adaptive activation function;

[0014] Design an adaptive sampling mechanism based on residuals;

[0015] Step 3: Implement numerical solution of wellbore multiphase flow model and gas-liquid distribution state inversion based on adaptive physical information neural network.

[0016] According to a preferred embodiment of the present invention, step 3 involves numerically solving the wellbore multiphase flow model and inverting the gas-liquid distribution state based on an adaptive physical information neural network; including:

[0017] Step 3.1: Input Preparation Stage;

[0018] The input variable is: spatial coordinates x With time coordinates t ;

[0019] The output variables are: predicting 5 physical quantities, including: wellbore pressure. P Gas phase velocity u g Liquid phase velocity u l Gas phase volume fraction E g and liquid volume fraction E l ;

[0020] Step 3.2: Establish physical equation constraints

[0021] Based on the law of conservation of mass, a continuity equation is established along the flow direction in the wellbore annulus. The coordinate along the flow direction in the wellbore annulus is set as z. A small element segment dz is selected for study. The cross-sectional area of ​​the annulus flow is A. The following partial differential equations are used as physical constraints:

[0022] Continuity equation for gas-phase producing sections:

[0023] (1);

[0024] The continuity equation for the liquid phase:

[0025] (2);

[0026] According to the law of conservation of momentum: the rate of change of momentum over time within a unit cell is equal to the sum of all external forces acting on that unit cell, therefore:

[0027]

[0028] (3);

[0029] Slip velocity relationship:

[0030] (4);

[0031] Volume fraction constraint equation:

[0032] (5);

[0033] In the above formula, The density of free gas at the temperature and pressure within the annulus, in kg·m³. -3 ; Density of drilling fluid, kg·m -3 ; , The upward return velocities of free gas and drilling fluid, respectively, are in m / s. -1 ; Let be the gas slip velocity, in m·s -1 ; , These are the integrals of free gas and drilling fluid, respectively, dimensionless; A is the annular cross-sectional area, in meters. 2 A= Where D is the inner diameter of the wellbore, in meters; d is the outer diameter of the drill string in different well sections, in meters; The mass of gas produced per unit time per unit thickness of reservoir, in kg·s -1 ·m -1 Rs represents the solubility of the gas in the drilling fluid, m 3 ·m -3 ; The density of the gas under standard conditions is expressed in kg·m³. -3 ; The volume coefficient of the drilling fluid in the local area is dimensionless.

[0034] Step 3.3: Design a physical information neural network; the physical information neural network includes an input layer, 6 hidden layers, and an output layer;

[0035] Step 3.4: Residual-based adaptive sampling; initially, training points are uniformly sampled from the spatial-temporal domain; during the middle of training, the residuals of each equation at different locations are calculated. r ( x Construct the residual probability density function. p ( x Automatically add sampling points in areas with large residuals;

[0036] Step 3.5: Training the physical information neural network;

[0037] Step 3.6: Gas-liquid distribution inversion and visualization: After the physical information neural network is trained, the gas-liquid distribution and pressure field output by the neural network can be used to realize the gas-liquid distribution inversion inside the wellbore.

[0038] According to a preferred embodiment of the present invention, the physical information neural network includes a multi-layer feedforward neural network, which includes an input layer, a hidden layer and an output layer.

[0039] The input layer is the spatiotemporal solution domain of the physical problem. The hidden layer consists of several layers of neurons and nonlinear activation functions to approximate the physical model. By optimizing the weights and biases, the input is gradually mapped to the output. The output layer obtains the approximate distribution of the physical field in the spatiotemporal solution domain through the nonlinear mapping of multiple hidden layers.

[0040] The output parameters of the physical information neural network are five physical quantities, including wellbore pressure. P Gas phase true velocity u g True velocity of liquid phase u l Gas phase volume fraction E g and liquid phase volume fraction E l ;

[0041] The physical information neural network is shown below:

[0042] (6);

[0043] (7);

[0044] In the formula, α 0 As the input layer of a physical information neural network, it consists of space x and time t Coordinate composition; α kFor the first k Layer output; For the first k-1 Layer output; σ It is a non-linear activation function; w k and b k The first k Layer weights and bias coefficients; L This represents the number of layers in the physical information neural network.

[0045] According to a preferred embodiment of the present invention, the loss function Loss of the physical information neural network includes data loss. Loss Data Residual loss Loss Res Initial loss Loss IC and boundary loss Loss BC As shown in equations (8) to (12):

[0046] (8);

[0047] (9);

[0048] (10);

[0049] (11);

[0050] (12);

[0051] In equations (8) to (12), , , , The weights for the loss function are set to 1 if both the data and the equations have been dimensionless; n, m, p, and k represent the nth observation data point, the mth internal sampling point, the pth initial condition point, and the kth boundary point, respectively. N , M , P , K This represents the total number of data points, interior points, initial points, and boundary points. y (·) represents the physical quantity predicted by the physical information neural network at (x, t). P , u g , u l , E g , E l The output of )I (·) represents the initial physical quantity predicted by the physical information neural network; B (·) represents the boundary physical quantity output by the physical information neural network at the boundary. y n , I p , B k These represent the observed values ​​corresponding to the data points, the initial condition values ​​corresponding to the initial points, and the boundary condition values ​​corresponding to the boundary points. ， The data is known. e i (·) is the output value of the physical information neural network, which is then input into the first... i The result after the physical equation q This represents the total number of physical equations.

[0052] The five physical equations for residual loss are shown in equations (13) to (17):

[0053] (13);

[0054] (14);

[0055] (15);

[0056] (16);

[0057] (17);

[0058] In the formula, ρ g The density of free gas at the temperature and pressure within the annulus, in kg·m³. -3 ; ρ l Density of drilling fluid, kg·m -3 ; u g , u l The upward return velocities of free gas and drilling fluid, respectively, are in m / s. -1 ; E g , E l These are the integrals of free gas and drilling fluid, respectively, and are dimensionless. A Let m be the cross-sectional area of ​​the annulus. 2 , ,in D The inner diameter of the well shaft is in meters (m). d The outer diameter of the drill string for different well sections, in meters; q gThe mass of gas produced per unit time per unit thickness of reservoir, in kg·s -1 ·m -1 ; R s m represents the solubility of the gas in the drilling fluid. 3 ·m -3 ; ρ gs The density of the gas under standard conditions is expressed in kg·m³. -3 ; B l The volume coefficient of the drilling fluid in the local area is dimensionless. u m The apparent flow rate of the mixed phase is in m·s. -1 , u m = u sl + u sg C0 is the gas phase distribution coefficient, which is dimensionless. u gr Let be the gas phase slip velocity, in m·s -1 ; P Pressure, Pa.

[0059] According to a preferred embodiment of the present invention, training a physical information neural network includes:

[0060] This paper utilizes the built-in objects of the TensorFlow framework to build the network, design physical constraints, calculate loss, optimize and update, and perform graph execution. The weight matrix of the physical information neural network is generated using the Xavier initialization method. The parameters of the physical information neural network are optimized using the Adam-LBFGS hybrid optimizer. Initially, a fixed sampling point strategy is used in conjunction with the Adam optimizer for preliminary training. Subsequently, the optimization is performed by switching to the L-BFGS-B optimizer. An adaptive sampling mechanism based on residuals is introduced to dynamically increase training points in high residual regions, thereby improving the resolution of key regions.

[0061] According to a preferred embodiment of the present invention, an adaptive activation function is designed, comprising:

[0062] The training process of the physical information neural network is as follows: find suitable weights w and bias terms b for each neuron in each layer of the physical information neural network so that the loss function gradually decreases and eventually reaches the global minimum; the update process of weights and bias terms is shown in equations (18) and (19):

[0063] (18);

[0064] (19);

[0065] In the formula, ηFor learning rate, η >0, adopt dynamic learning rate, that is, after setting the initial learning rate, it decays once every certain number of training steps, so that the learning rate decreases in a step-like manner; J m For the first m The loss function at the next iteration; For the first m The weight vector of the physical information neural network at the next iteration. For the first m+ The weight vector of the physical information neural network in the first iteration; loss function Weights The gradient; For the first m The bias vector of the physical information neural network at the next iteration. For the first m The bias vector at +1 iterations loss function For bias b The gradient;

[0066] In activation function Introducing hyperparameters c To obtain the output of the k-th layer As shown in equation (20):

[0067] (20);

[0068] in, This refers to the output of the (k-1)th layer; It refers to the first k Layer bias vector; It refers to the first k Layer weight vector; c It is an adjustable parameter used to control the slope of the activation function, thereby changing the topology of the loss function during neural network training; c The gradient descent method is used to optimize the weights and biases together, as shown in equation (21):

[0069] (twenty one);

[0070] in, For the first m The slope of the activation function is controlled in the next iteration. For the first m+ In the first iteration, control the slope of the activation function; loss function For hyperparameters c The gradient.

[0071] Further preferred, the adaptive activation function is the Tanh function.

[0072] According to a preferred embodiment of the present invention, an adaptive sampling mechanism based on residuals is designed, comprising:

[0073] A probability density function based on residuals is used to increase the number of sampling points. As shown in equation (22):

[0074] (twenty two);

[0075] in, X For the candidate point set, by space x and time t Coordinate composition; r ( X ) represents the residual value corresponding to the candidate point set; This is a threshold used to filter out small residuals; when the residual value is less than... When the logarithmic function returns a negative value, the max function is used to convert the residual value to 0; Ω refers to the entire space-time domain.

[0076] For the five physical equations shown in equations (13) to (17) e 1 - e 5 The corresponding residual loss is used to calculate its probability density function, thus giving each physical equation a different sampling point distribution. Specifically, in each round of sampling, new sampling points at different locations are added to these five physical equations to capture the residual distribution characteristics of each equation. The implementation process of the residual-based adaptive sampling mechanism is as follows:

[0077] 1) Initialization: Determine the initial set of sampling points and set the relevant parameters for adaptive sampling, including the candidate point set, the number of iterations to start sampling, the number of new points added in each round of sampling, the total number of samplings, and the corresponding loop period. e 1 - e 5 The sampling priority and residual threshold parameter ε are determined; at the same time, a neural network for solving the PDE is constructed, and an appropriate number of hidden layers and neuron size are set.

[0078] 2) Calculate the residuals and construct the probability density function: Combine automatic differentiation technology to calculate the residual loss at the candidate point, that is, the residuals corresponding to the 5 physical equations; at the same time, construct the probability density function corresponding to each physical equation according to equation (22);

[0079] 3) Determine new sampling points: Using the five constructed probability density functions, extract five new sampling points from the candidate point set and add the five new sampling points to the training set;

[0080] 4) Network training and optimization: On the updated training set, use the optimization algorithm to train the physical information neural network. After reaching the set number of training iterations, start a new round of sampling and repeat steps 2)-4) until the total loss is lower than the set threshold or the maximum number of iterations is reached, and obtain the adaptive physical information neural network.

[0081] According to a preferred embodiment of the present invention, adaptive physical information neural network structure optimization refers to:

[0082] In the adaptive physical information neural network, there are 6 hidden layers and 60 neurons.

[0083] A computer device includes a memory and a processor, the memory storing a computer program, and the processor executing the computer program to implement the steps of numerical solution of a wellbore multiphase flow model and gas-liquid distribution state inversion method.

[0084] A computer-readable storage medium having a computer program stored thereon, wherein the computer program, when executed by a processor, implements the steps of a method for numerically solving a multiphase flow model in a wellbore and inverting the gas-liquid distribution state.

[0085] Numerical solution system for multiphase flow model in wellbore and gas-liquid distribution state inversion system, including:

[0086] The physical information neural network construction module is configured to: construct and train a physical information neural network for dynamic simulation of multiphase flow in drilling wellbore and inversion of overflow gas distribution; determine the input and output of the physical information neural network; determine the loss function of the physical information neural network; and train the physical information neural network.

[0087] The adaptive physical information neural network design module is configured to: design an adaptive optimization algorithm to optimize the final solution accuracy and convergence speed of the physical information neural network, thereby obtaining the adaptive physical information neural network; design an adaptive activation function; and design an adaptive sampling mechanism based on residuals.

[0088] The module for numerical solution of multiphase flow model in wellbore and inversion of gas-liquid distribution state is configured to: realize numerical solution of multiphase flow model in wellbore and inversion of gas-liquid distribution state based on adaptive physical information neural network.

[0089] The beneficial effects of this invention are as follows:

[0090] The adaptive activation function in this invention can effectively accelerate the optimization process, enabling the network to approach a better solution more quickly. Final application examples show that the adaptive PINN performs well in gas flow processes, effectively capturing key features during gas intrusion and recirculation exhaust, especially exhibiting small errors in porosity prediction. This invention effectively solves the problem that traditional numerical methods typically require high-precision mesh generation and large computational resources. Attached Figure Description

[0091] Figure 1 This is a schematic diagram of the PINN network structure;

[0092] Figure 2 For different hyperparameters c Schematic diagram of activation functions at various values;

[0093] Figure 3 This is a schematic diagram of a residual-based three-stage adaptive sampling process (taking the residuals of the gas phase continuity equation and momentum equation as an example);

[0094] Figure 4 Here is a flowchart of the annular air-liquid distribution state inversion based on adaptive PINN;

[0095] Figure 5 This is a schematic diagram illustrating the impact of the adaptive activation function on the loss.

[0096] Figure 6 The simulation results of the adaptive PINN for the gas intrusion process are shown in the figure. Detailed Implementation

[0097] The present invention will be further defined below with reference to the accompanying drawings and embodiments, but is not limited thereto.

[0098] Example 1

[0099] Numerical solution and gas-liquid distribution state inversion method for wellbore multiphase flow model driven by physical information neural network, including:

[0100] Step 1: Construct and train a physical information neural network for dynamic simulation of multiphase flow in drilling wellbore and inversion of overflow gas distribution state;

[0101] Determine the input and output of the physical information neural network;

[0102] Determine the loss function of the physical information neural network;

[0103] Training a physical information neural network;

[0104] Step 2: Design an adaptive optimization algorithm to optimize the final solution accuracy and convergence speed of the physical information neural network, and obtain the adaptive physical information neural network;

[0105] Design an adaptive activation function;

[0106] Design an adaptive sampling mechanism based on residuals;

[0107] Step 3: Implement numerical solution of wellbore multiphase flow model and gas-liquid distribution state inversion based on adaptive physical information neural network.

[0108] Example 2

[0109] The difference between the numerical solution and gas-liquid distribution state inversion method of wellbore multiphase flow model driven by physical information neural network as described in Example 1 and the method described in Example 1 is as follows:

[0110] In step 3, the numerical solution of the wellbore multiphase flow model and the inversion of the gas-liquid distribution state are realized based on the adaptive physical information neural network; including:

[0111] Step 3.1: Input Preparation Stage;

[0112] The input variable is: spatial coordinates x With time coordinates t ;

[0113] The output variables are: five physical quantities predicted by the network in this invention, including: wellbore pressure. P Gas phase velocity u g Liquid phase velocity u l Gas phase volume fraction E g and liquid volume fraction E l ;

[0114] Step 3.2: Establish physical equation constraints

[0115] Based on the law of conservation of mass, a continuity equation is established along the flow direction in the wellbore annulus. The coordinate along the flow direction in the wellbore annulus is set as z. A small element segment dz is selected for study. The cross-sectional area of ​​the annulus flow is A. The following partial differential equations are used as physical constraints:

[0116] Continuity equation for gas-phase producing sections:

[0117] (1);

[0118] The continuity equation for the liquid phase:

[0119] (2);

[0120] According to the law of conservation of momentum: the rate of change of momentum over time within a unit cell is equal to the sum of all external forces acting on that unit cell, therefore:

[0121]

[0122] (3);

[0123] Slip velocity relationship:

[0124] (4);

[0125] Volume fraction constraint equation:

[0126] (5);

[0127] In the above formula, The density of free gas at the temperature and pressure within the annulus, in kg·m³. -3 ; Density of drilling fluid, kg·m -3 ; , The upward return velocities of free gas and drilling fluid, respectively, are in m / s. -1 ; Let be the gas slip velocity, in m·s -1 ; , These are the integrals of free gas and drilling fluid, respectively, dimensionless; A is the annular cross-sectional area, in meters. 2 A= Where D is the inner diameter of the wellbore, in meters; d is the outer diameter of the drill string in different well sections, in meters; The mass of gas produced per unit time per unit thickness of reservoir, in kg·s -1 ·m -1 Rs represents the solubility of the gas in the drilling fluid, m 3 ·m -3 ; The density of the gas under standard conditions is expressed in kg·m³. -3 ; The volume coefficient of the drilling fluid in the local area is dimensionless.

[0128] Step 3.3: Design the physical information neural network; the physical information neural network includes an input layer, 6 hidden layers, and an output layer; the activation function is the adaptive activation function designed in Step 2; the weight parameters are initialized using Xavier; the training optimizer is the Adam-LBFGS hybrid optimizer;

[0129] Step 3.4: Residual-based adaptive sampling; initially, training points are uniformly sampled from the spatial-temporal domain; during the middle of training, the residuals of each equation at different locations are calculated. r ( x Construct the residual probability density function. p ( xAutomatically add sampling points in areas with large residuals;

[0130] Step 3.5: Training the physical information neural network; use the Adam optimizer for initial training (4000 steps), then switch to the L-BFGS-B optimizer for training. After reaching the set maximum number of iterations in the final stage, the network completes training.

[0131] Step 3.6: Gas-Liquid Distribution Inversion and Visualization: After the physical information neural network is trained, the gas-liquid distribution output by the neural network is visualized. and pressure field P ( x , t This allows for the inversion of gas-liquid distribution inside the wellbore.

[0132] Adaptive physical information neural network structure optimization;

[0133] Adaptive activation function;

[0134] PINN network parameters;

[0135] The annular air-liquid distribution and wellbore pressure characteristics during gas intrusion and circulating exhaust are simulated using an adaptive physical information neural network. The output results are compared with numerical solutions under high-density discrete grid conditions to verify the feasibility and accuracy of the adaptive PINN.

[0136] In this invention, the PINN network structure design achieves efficient solution of partial differential equations (PDEs) describing physical systems under physical constraints. Compared to traditional numerical methods, PINN directly guides network learning through the physical model, avoiding truncation errors caused by equation discretization. It also enables rapid inference of physical parameters and completes the fitting and inference of physical models even with limited or unlabeled data. The PINN network in this invention is as follows: Figure 1 As shown.

[0137] The physical information neural network includes a multi-layer feedforward neural network, which consists of an input layer, a hidden layer, and an output layer.

[0138] The input layer is the spatiotemporal solution domain of the physical problem, ensuring that the network can cover the entire range of changes in the physical field; the hidden layer is the core of the neural network, consisting of several layers of neurons and nonlinear activation functions, used to approximate the physical model, and gradually maps the input to the output by optimizing the weights and biases; the output layer obtains the approximate distribution of the physical field in the spatiotemporal solution domain through the nonlinear mapping of multiple hidden layers.

[0139] The input layer is shown in equation (6) below, where x For spatial coordinates, tUsing time as the coordinate, in actual calculations, the original input data... x After normalization to the range [-1, 1], the output is sent to the hidden layer. The hidden layer is shown in Equation (7) below. In actual operation, each layer will combine the input of the previous layer with the input of the previous layer. w Perform a linear transformation and add a bias coefficient. b The result is then passed to the next layer via the activation function designed in step 2. Unlike the hidden layers, the output layer does not have an activation function; it utilizes the weights of the last layer in the network. w and bias coefficient b Performing linear transformations and biasing is a typical approach for the output layer of a regression task, with the aim of outputting continuous values.

[0140] The output parameters of the physical information neural network are five physical quantities, including wellbore pressure. P Gas phase true velocity u g True velocity of liquid phase u l Gas phase volume fraction E g and liquid phase volume fraction E l ;

[0141] The physical information neural network is shown below:

[0142] (6);

[0143] (7);

[0144] In the formula, α 0 As the input layer of a physical information neural network, it consists of space x and time t Coordinate composition; α k For the first k Layer output; For the first k-1 Layer output; σ It is a non-linear activation function; w k and b k The first k Layer weights and bias coefficients; L This represents the number of layers in the physical information neural network.

[0145] The loss function of a physical information neural network includes data loss. Loss Data Residual loss Loss Res Initial loss Loss IC and boundary loss Loss BC As shown in equations (8) to (12):

[0146] Among them, data loss Loss Data Similar to the loss function of traditional machine learning algorithms, it represents the data error between the network output value and the actual observed value; while the residual loss... Loss Res This indicates the degree of agreement between the output value and the physical equations within the solution domain. It is used to ensure that the network prediction results follow existing physical laws. It is mainly achieved by calculating the partial derivatives of the output variables using automatic differentiation techniques and substituting them into the partial differential equations of the physical model to obtain the corresponding residuals; initial loss. Loss IC and boundary loss Loss BC This compares the network's output values ​​at the initial and boundary points with the determined initial and boundary conditions.

[0147] (8);

[0148] (9);

[0149] (10);

[0150] (11);

[0151] (12);

[0152] In equations (8) to (12), , , , The weights for the loss function are set to 1 if both the data and the equations have been dimensionless; n, m, p, and k represent the nth observation data point, the mth internal sampling point, the pth initial condition point, and the kth boundary point, respectively. N , M , P , K This represents the total number of data points, interior points, initial points, and boundary points. y (·) I (·) B (·) is actually the output of the same physical information at different sampling positions of the neural network. y (·) represents the physical quantity predicted by the physical information neural network at (x, t). P , u g ,u l , E g , E l The output of ) I (·) represents the initial physical quantities predicted by the physical information neural network, such as initial pressure, initial gas volume fraction, initial liquid volume fraction, initial gas flow rate, and initial liquid flow rate. B (·) represents the boundary physical quantity output by the physical information neural network at the boundary; such as wellhead pressure, bottom hole velocity, etc. y n , I p , B k These represent the observed values ​​corresponding to the data points, the initial condition values ​​corresponding to the initial points, and the boundary condition values ​​corresponding to the boundary points. ， The data is known. e i (·) is the output value of the physical information neural network, which is then input into the first... i The result after the physical equation q This represents the total number of physical equations.

[0153] The loss associated with 'n' is the data loss, representing the nth observation data point. It originates from the location of measured physical quantities in the field or numerical simulation and serves as guiding data for the network fitting. The loss associated with 'm' is the residual loss, representing the mth internal sampling point. It originates from within the wellbore solution domain and is used to force the network to satisfy the multiphase flow physical equations (Equations (8)-(12)) where no observation data is available. The loss associated with 'p' is the initial loss, representing the pth initial condition point. It corresponds to the known state at time t=0, ensuring that the network solution satisfies the initial conditions. The loss associated with 'k' is the boundary loss, representing the kth boundary point. It corresponds to the conditions at boundary locations such as the wellhead and bottom, and is used to ensure that the physical constraints of the solution domain boundary satisfy the boundary conditions.

[0154] The five physical equations for residual loss are shown in equations (13) to (17):

[0155] (13);

[0156] (14);

[0157] (15);

[0158] (16);

[0159] (17);

[0160] In the formula, ρ g The density of free gas at the temperature and pressure within the annulus, in kg·m³. -3 ; ρ l Density of drilling fluid, kg·m -3 ; u g , u l The upward return velocities of free gas and drilling fluid, respectively, are in m / s. -1 ; E g , E l These are the integrals of free gas and drilling fluid, respectively, and are dimensionless. A Let m be the cross-sectional area of ​​the annulus. 2 , ,in D The inner diameter of the well shaft is in meters (m). d The outer diameter of the drill string for different well sections, in meters; q g The mass of gas produced per unit time per unit thickness of reservoir, in kg·s -1 ·m -1 ; R s m represents the solubility of the gas in the drilling fluid. 3 ·m -3 ; ρ gs The density of the gas under standard conditions is expressed in kg·m³. -3 ; B l The volume coefficient of the drilling fluid in the local area is dimensionless. u m The apparent flow rate of the mixed phase is in m·s. -1 , u m = u sl + u sg C0 is the gas phase distribution coefficient, which is dimensionless. u gr Let be the gas phase slip velocity, in m·s -1 ; P Pressure, Pa.

[0161] Training a physical information neural network; including:

[0162] This paper utilizes TensorFlow's built-in objects to construct the network (using TensorFlow.Variable, TensorFlow.matmul, and TensorFlow.tanh to build the forward computation graph), design physical constraints (using TensorFlow.gradients to automatically calculate partial derivatives and construct PDE residuals), calculate loss (using TensorFlow.reduce_mean and TensorFlow.square to combine data error and physical error), optimize and update (using TensorFlow.train.AdamOptimizer and TensorFlow.contrib.opt.ScipyOptimizerInterface to automatically backpropagate and update weights), and execute the graph (using TensorFlow.Session.run to run the computation graph). The Xavier initialization method, proposed by Xavier Glorot and Yoshua Bengio in 2010, is used to generate the weight matrix of the physical information neural network. Xavier is a method for initializing weights in deep neural networks, aiming to maintain consistency in signal variance during forward and backward propagation to address the vanishing and exploding gradient problems. The weights are specified in the Xavier initialization method. The variance should be set to be equal to the reciprocal of the average of the number of input nodes and the number of output nodes, which effectively alleviates the gradient vanishing or gradient exploding problem during training. Furthermore, by comparing commonly used optimizers and their advantages and disadvantages, the Adam-LBFGS hybrid optimizer is ultimately adopted to optimize the parameters of the physical information neural network. Commonly used optimizers and their advantages and disadvantages are shown in Table 1. Initial training is performed using a fixed sampling point strategy combined with the Adam optimizer, followed by switching to the L-BFGS-B optimizer for efficient optimization. The residual-based adaptive sampling mechanism designed in step 2 is introduced to dynamically increase training points in high residual regions, improving the resolution of key regions.

[0163] Table 1. Commonly used optimizers and their advantages and disadvantages;

[0164]

[0165] Design an adaptive activation function; including:

[0166] The training process of the physical information neural network is as follows: find suitable weights w and bias terms b for each neuron in each layer of the physical information neural network so that the loss function gradually decreases and eventually reaches the global minimum; the update process of weights and bias terms is shown in equations (18) and (19):

[0167] (18);

[0168] (19);

[0169] In the formula, η For learning rate, η >0. During training, a large learning rate may exceed the global minimum, while a small learning rate, although it will slowly move towards the global minimum, will increase the computational cost. Therefore, this invention adopts a dynamic learning rate, that is, after setting the initial learning rate, it decays once every certain number of training steps, so that the learning rate decreases in a step-like manner, thereby improving training efficiency while ensuring convergence. J m For the first m The loss function at the next iteration; For the first m The weight vector of the physical information neural network at the next iteration. For the first m+ The weight vector of the physical information neural network in the first iteration; loss function Weights The gradient; For the first m The bias vector of the physical information neural network at the next iteration. For the first m The bias vector at +1 iterations loss function For bias b The gradient;

[0170] To optimize the physical information neural network for best performance, the activation function... Introducing hyperparameters c To obtain the output of the k-th layer As shown in equation (20):

[0171] (20);

[0172] in, This refers to the output of the (k-1)th layer; It refers to the first k Layer bias vector; It refers to the first k Layer weight vector; c It is an adjustable parameter used to control the slope of the activation function, thereby changing the topological structure of the loss function during neural network training; similarly, this parameter also needs to be optimized, and this invention incorporates it into the training process. c The gradient descent method is used to optimize the weights and biases together, as shown in equation (21):

[0173] (twenty one);

[0174] in, For the first m The slope of the activation function is controlled in the next iteration. For the first m+ In the first iteration, control the slope of the activation function; loss function For hyperparameters c The gradient.

[0175] Figure 2 Demonstrates different hyperparameters c The changing trends of various activation functions under different values. From Figure 2 As can be seen from this, the slope of the activation function will change with... c The slope changes with the adjustment of the value, which has a significant impact on gradient propagation and model optimization during training. A larger slope can accelerate the optimization process and help to quickly explore the solution space, but may lead to gradient explosion; a smaller slope can make training smoother and is suitable for fine-tuning when approaching the optimal solution, but may cause gradient vanishing.

[0176] Due to the high complexity of neural networks, selecting a suitable activation function often requires manual experience, which necessitates repeated trials and increases training costs. Therefore, this invention aims to optimize the neural network to achieve the best performance. The adaptive activation function of this invention introduces hyperparameters based on the traditional Tanh function. c And it is incorporated into the training process for optimization, which is the origin of "adaptive". As shown in Equation (20), it is the adaptive activation function of this invention.

[0177] The adaptive activation function is the Tanh function.

[0178] Design an adaptive sampling mechanism based on residuals; including:

[0179] The adaptive sampling mechanism uses a residual-based probability density function to increase the number of sampling points. More points are added to regions with high probability (high residual regions), while fewer points are added to regions with low probability (low residual regions). These points are then added to the training set and trained simultaneously with the initial fixed sampling points to optimize the network parameters, thereby improving the PINN network's ability to solve problems in regions with drastic local feature changes. The residual-based probability density function... As shown in equation (22):

[0180] (twenty two);

[0181] in, X For the candidate point set, by space x and time t Coordinate composition; r( X ) represents the residual value corresponding to the candidate point set; This is a threshold used to filter out small residuals; when the residual value is less than... When the logarithmic function returns a negative value, the max function is used to turn the residual value to 0, that is, the sampling probability of the region with particularly small residual is 0, thus avoiding the waste of computing resources; Ω refers to the entire space-time domain.

[0182] Taking the logarithm of the residuals compresses the dynamic range of large residual values, thereby avoiding excessive clustering of new sampling points.

[0183] For the five physical equations shown in equations (13) to (17) e 1 - e 5 The corresponding residual loss is used to calculate its probability density function, thus giving each physical equation a different sampling point distribution. Specifically, in each round of sampling, new sampling points at different positions are added to these five physical equations to more accurately capture the residual distribution characteristics of each equation, further improving the training efficiency and solution accuracy of the PINN network. The implementation process of the residual-based adaptive sampling mechanism is as follows:

[0184] 1) Initialization: Determine the initial set of sampling points and set the relevant parameters for adaptive sampling, including the candidate point set, the number of iterations to start sampling, the number of new points added in each round of sampling, the total number of samplings, and the corresponding loop period. e 1 - e 5 The sampling priority (which can be handled by allocating different numbers of sampling points) and the residual threshold parameter ε (which divides the residual into two regions, high residual and low residual, sampling in the high residual region and ignoring the low residual region) are determined. At the same time, a neural network for solving the PDE is constructed, and an appropriate number of hidden layers and neuron size are set.

[0185] The PDE is a set of physical constraint partial differential equations describing multiphase flow in a wellbore, specifically corresponding to equations (1) to (5) in this invention. These equations are written as a residual term, and the residual is calculated using the automatic differentiation function of the neural network and used as part of the loss function, so that the network training and solution results conform to physical laws. The specific architecture data of the neural network for solving the PDE is shown in Table 3. The specific settings for the number of hidden layers and the size of neurons are as follows: the number of hidden layers is 6, and the number of neurons in each layer is 60.

[0186] 2) Calculate the residuals and construct the probability density function: Combine automatic differentiation technology to calculate the residual loss at the candidate point, that is, the residuals corresponding to the 5 physical equations; at the same time, construct the probability density function corresponding to each physical equation according to equation (22);

[0187] In Physical Information Neural Networks (PINNs), residuals refer to the imbalance on both sides of each physical constraint equation after substituting the neural network's predicted solution into the original equation. Residual calculation in PINNs can be achieved using automatic differentiation techniques. e 1 - e 5 The five physical equations correspond to the residuals of the gas phase mass conservation equation, the liquid phase mass conservation equation, the momentum conservation equation, the drift flow model equation, and the volume fraction constraint equation, respectively. The specific implementation steps are as follows: ① Use TensorFlow.GradientTape to start recording the computation operation; ② Use TensorFlow.GradientTape.watch to ensure that the spatiotemporal coordinates are tracked; ③ Use TensorFlow.GradientTape.gradient to calculate the derivative; ④ Substitute the calculated derivative into the physical equation to obtain the residual; ⑤ For each candidate point X, calculate the ratio of the absolute value of the residual to the threshold. ⑥ Take the logarithm of the residual ratio to make the differences between residuals of different sizes relatively smooth; ⑦ Apply the max function for threshold filtering, that is, do not sample regions with particularly small residuals, saving computational resources; ⑧ Calculate the integral over the entire space-time domain Ω, and divide the value of each point by this integral value to obtain the standard probability distribution.

[0188] 3) Determine new sampling points: Using the five constructed probability density functions, extract five new sampling points from the candidate point set and add the five new sampling points to the training set;

[0189] 4) Network training and optimization: On the updated training set, the physical information neural network is trained using optimization algorithms (mainly including two parts: adaptive activation function and residual-based adaptive sampling mechanism, which correspond to the two parts in step 2). After reaching the set number of training iterations, a new round of sampling is started, and steps 2)-4) are repeated until the total loss is lower than the set threshold or the maximum number of iterations is reached, thus obtaining the adaptive physical information neural network.

[0190] To intuitively demonstrate the residual-based adaptive sampling process, uniform sampling is first used to select internal points ( x × t=100×100), and similarly, two optimizers were used in the training process, one for the early stage and one for the later stage. After training with the Adam optimizer for 4000 steps, it was switched to the L-BFGS-B optimizer. The corresponding probability density function was calculated using equation (22) based on the residuals corresponding to the five physical constraint equations, thereby increasing the sampling points in the region with large residuals to improve the accuracy of the model. This invention performed three adaptive samplings. Taking the residuals of the gas phase continuity equation and momentum equation as examples, 400×5 sampling points were added each time. The adaptive sampling results are as follows. Figure 3 As shown, the black dots represent the initial sampling points, and the red crosses represent newly added sampling points. Figure 3 In the middle, the vertical axis x This represents spatial coordinates, specifically well depth, in meters (m); the horizontal axis... t This represents the time coordinate, i.e., the time when the air intrusion occurs, measured in seconds. Through adaptive sampling, the model can better capture the dynamic changes of the flow front. It can not only identify high residual regions but also optimize the distribution of sampling points through multiple iterations, avoiding over-concentration in a small area, thereby improving the model's generalization ability.

[0191] Figure 3 (a) shows the first adaptive sampling result based on the residuals of the gas-phase continuity equation; from Figure 3 As can be seen in (a), the model can better capture the dynamic changes of the flow front through adaptive sampling. The newly added points are mainly concentrated in the lower part of the wellbore and some wellhead locations, indicating that the residual of the gas phase continuity equation in this region is large, that is, the calculation error caused by the multiphase flow region and the boundary location is high, requiring more computing resources. Figure 3 (b) shows the first adaptive sampling result based on the momentum equation residuals; comparison Figure 3 The sampling results shown in (b) reveal that the newly added points are relatively dispersed, indicating that the momentum equation residuals are relatively evenly distributed in the spatiotemporal domain, eliminating the need for additional attention to specific regions. Furthermore, the overall trend shows that the distribution of sampling points gradually becomes more uniform, covering more areas with large residuals. This demonstrates that the adaptive algorithm can not only identify high residual regions but also optimize the distribution of sampling points through multiple iterations, avoiding over-concentration in a small area and thus improving the model's generalization ability.

[0192] Figure 3 (c) shows the second adaptive sampling result based on the residual of the gas phase continuity equation; Figure 3 (d) represents the second adaptive sampling result based on the momentum equation residual; Figure 3 (e) represents the third adaptive sampling result based on the residuals of the gas-phase continuity equation; Figure 3 In the middle (f), the result of the third adaptive sampling is based on the residual of the momentum equation.

[0193] The inversion process of the annular air-liquid distribution state based on adaptive PINN is as follows: Figure 4 As shown.

[0194] To find a suitable depth and width for the adaptive physical information neural network, this invention tested the computational error under different hidden layers and numbers of neurons. First, the number of internal points was set to 10,000 and the number of boundary points to 180. Then, different combinations of hidden layers (2, 4, 6, and 8 layers) and different numbers of neurons (20, 40, 60, 80, and 100) were used for training. The computational error was also calculated using relative L2 error. The results are shown in Table 2.

[0195] Table 2. Error table for different hidden layers and different numbers of neurons;

[0196]

[0197] As can be seen from Table 2, the computational error generally decreases as the number of hidden layers and the number of neurons per layer increase, indicating that deeper networks can better fit complex problems and improve the solution accuracy. When the number of hidden layers increases to more than 6 and the number of neurons exceeds 60, the trend of error reduction slows down and may even lead to an increase in error. This indicates that when the number of network layers and neurons exceeds a certain range, the improvement in computational power is limited, while the computational cost increases and overfitting may occur, affecting the generalization ability.

[0198] Adaptive physical information neural network structure optimization refers to:

[0199] In the adaptive physical information neural network, there are 6 hidden layers and 60 neurons. Network structure optimization mainly explores how to achieve a balance between computational efficiency and accuracy when the number of hidden layers and neurons per layer are reached. Further increasing the number of layers or neurons results in longer computation time and may even increase the error.

[0200] To investigate the impact of adaptive activation functions on the training loss of PINN under different training strategies, the optimization method went through two different stages during training: the Adam optimizer was used in the early and middle stages of training, while the L-BFGS-B optimizer was used in the later stages of training. The trend of loss change during training is shown in the figure. Figure 5 As shown, Figure 5 In the diagram, the green curve represents the dynamic changes of the hyperparameter c during the training process.

[0201] Looking at the loss curves, in the early stages of training (before approximately 1500 steps), the losses for both using the adaptive activation function (red curve) and not using the adaptive activation function (blue curve) decrease rapidly, exhibiting the same downward trend, consistent with the hyperparameters. cThe trend of change is consistent (remaining around 1); in the middle of training (approximately 1500 to 4000 steps), the loss curve shows significant fluctuations. This is mainly due to the Adam optimizer used, which combines momentum and adaptive learning rate, enabling self-adjustment during optimization, but may cause oscillations in the loss function. Furthermore, the slope of the activation function varies with hyperparameters... c The changes in parameters lead to significant variations in the gradient update magnitude, making the fluctuations in the red curve more pronounced. In the later stages of training (after 4000 steps), the optimizer switches to L-BFGS-B, at which point the loss curve tends to stabilize and decreases further at a faster rate. During this phase, the red curve consistently maintains a lower loss value, indicating that the adaptive activation function effectively accelerates the optimization process, allowing the network to approach a better solution more quickly. The blue curve, however, decreases relatively slowly, indicating that the traditional activation function converges more slowly under the same optimization conditions. Meanwhile, the hyperparameters... c During this phase, the gradient increases rapidly at first and then stabilizes. This dynamic adjustment mechanism helps to provide appropriate gradient updates at different training stages.

[0202] In summary, employing an adaptive activation function combined with the Adam+L-BFGS-B optimization strategy can effectively improve the training efficiency of neural networks. This is achieved by dynamically adjusting hyperparameters. c This method can accelerate optimization, enhance stability, and allow the loss to converge faster, while avoiding getting trapped in local optima. Compared to the traditional approach with a fixed activation function, this method has stronger adaptability and optimization capabilities, and can significantly improve the overall performance of neural networks.

[0203] The optimized adaptive PINN network parameter settings are shown in Table 3.

[0204] Table 3. Optimized Adaptive PINN Network Parameters;

[0205]

[0206] Adaptive PINN was used to simulate the annular air-liquid distribution and wellbore pressure characteristics during gas intrusion and circulating exhaust. The output results were compared with the numerical solution under high-density discrete grid conditions to verify the feasibility and accuracy of adaptive PINN.

[0207] After three rounds of adaptive sampling and model training optimization, the inversion results of adaptive PINN with respect to annular porosity and wellbore pressure during gas intrusion were obtained, such as... Figure 6 As shown. Figure 6 (a) shows the distribution characteristics of the numerical solution of porosity; Figure 6 (b) shows the distribution characteristics of the numerical solution of wellbore pressure; Figure 6 (c) shows the porosity inversion results; Figure 6The middle (d) diagram shows the wellbore pressure inversion results; Figure 6 (e) is the absolute error diagram for porosity inversion; Figure 6 (f) is the absolute error diagram of wellbore pressure inversion; the results show the process of continuous gas intrusion and gradual diffusion and migration from the bottom of the well upwards, in which the inversion results of porosity ( Figure 6 The solution in (c) is very close to the original numerical solution, accurately reflecting the upward diffusion trend of gas along the wellbore. This indicates that PINN successfully captured the dynamic changes in porosity, combined with the absolute error distribution ( Figure 6 As can be seen from (e)), the inversion error of porosity is relatively small, mainly concentrated at the gas front; the inversion results of wellbore pressure ( Figure 6 The variation of the solution in (d) is basically consistent with that of the numerical solution, but compared with the porosity, it shows some deviations in local areas, especially in the multiphase flow region when the time is close to 600 s, with a maximum error of about 0.5 MPa.

[0208] Example 3

[0209] A computer device includes a memory and a processor, wherein the memory stores a computer program, and the processor executes the computer program to implement the steps of the method for numerical solution of wellbore multiphase flow model and gas-liquid distribution state inversion driven by physical information neural network as described in Embodiment 1 or 2.

[0210] Example 4

[0211] A computer-readable storage medium having a computer program stored thereon, wherein the computer program, when executed by a processor, implements the steps of the method for numerical solution of wellbore multiphase flow model and gas-liquid distribution state inversion driven by physical information neural network as described in Embodiment 1 or 2.

[0212] Example 5

[0213] A physical information neural network-driven numerical solution system for wellbore multiphase flow models and a gas-liquid distribution state inversion system, including:

[0214] The physical information neural network construction module is configured to: construct and train a physical information neural network for dynamic simulation of multiphase flow in drilling wellbore and inversion of overflow gas distribution; determine the input and output of the physical information neural network; determine the loss function of the physical information neural network; and train the physical information neural network.

[0215] The adaptive physical information neural network design module is configured to: design an adaptive optimization algorithm to optimize the final solution accuracy and convergence speed of the physical information neural network, thereby obtaining the adaptive physical information neural network; design an adaptive activation function; and design an adaptive sampling mechanism based on residuals.

[0216] The module for numerical solution of multiphase flow model in wellbore and inversion of gas-liquid distribution state is configured to: realize numerical solution of multiphase flow model in wellbore and inversion of gas-liquid distribution state based on adaptive physical information neural network.

Claims

1. A numerical solution method for multiphase flow model in wellbore and a method for inverting gas-liquid distribution state, characterized in that, include: Step 1: Construct and train a physical information neural network for dynamic simulation of multiphase flow in drilling wellbore and inversion of overflow gas distribution state; Determine the input and output of the physical information neural network; Determine the loss function of the physical information neural network; Training a physical information neural network; Step 2: Design an adaptive optimization algorithm to optimize the final solution accuracy and convergence speed of the physical information neural network, and obtain the adaptive physical information neural network; Design an adaptive activation function; Design an adaptive sampling mechanism based on residuals; Step 3: Numerical solution of wellbore multiphase flow model and inversion of gas-liquid distribution state based on adaptive physical information neural network; In step 3, the numerical solution of the wellbore multiphase flow model and the inversion of the gas-liquid distribution state are realized based on the adaptive physical information neural network; include: Step 3.1: Input Preparation Stage; The input variable is: spatial coordinates x With time coordinates t ; The output variables are: predicting 5 physical quantities, including: wellbore pressure. P Gas phase velocity u g Liquid phase velocity u l Gas phase volume fraction E g and liquid volume fraction E l ; Step 3.2: Establish physical equation constraints Based on the law of conservation of mass, a continuity equation is established along the flow direction in the wellbore annulus. The coordinate along the flow direction in the wellbore annulus is set as z. A small element segment dz is selected for study. The cross-sectional area of ​​the annulus flow is A. The following partial differential equations are used as physical constraints: Continuity equation for gas-phase producing sections: (1)； The continuity equation for the liquid phase: (2)； According to the law of conservation of momentum: the rate of change of momentum over time within a unit cell is equal to the sum of all external forces acting on that unit cell, therefore: (3)； Slip velocity relationship: (4)； Volume fraction constraint equation: (5)； In the above formula, The density of free gas at the temperature and pressure within the annulus, in kg·m³. -3 ; Density of drilling fluid, kg·m -3 ; , The upward return velocities of free gas and drilling fluid, respectively, are in m / s. -1 ; Let be the gas slip velocity, in m·s -1 ; , These are the integrals of free gas and drilling fluid, respectively, dimensionless; A is the annular cross-sectional area, in meters. 2 A= Where D is the inner diameter of the wellbore, in meters; d is the outer diameter of the drill string in different well sections, in meters; The mass of gas produced per unit time per unit thickness of reservoir, in kg·s -1 ·m -1 Rs represents the solubility of the gas in the drilling fluid, m 3 ·m -3 ; The density of the gas under standard conditions is expressed in kg·m³. -3 ; C0 is the volume coefficient of the drilling fluid in the local area, dimensionless; C0 is the gas phase distribution coefficient, dimensionless. Step 3.3: Design a physical information neural network; the physical information neural network includes an input layer, 6 hidden layers, and an output layer; Step 3.4: Residual-based adaptive sampling; initially, training points are uniformly sampled from the spatial-temporal domain; during the middle of training, the residuals of each equation at different locations are calculated. r ( x Construct the residual probability density function. p ( x Automatically add sampling points in areas with large residuals; Step 3.5: Training the physical information neural network; Step 3.6: Gas-liquid distribution inversion and visualization: After the physical information neural network is trained, the gas-liquid distribution and pressure field output by the neural network are used to realize the gas-liquid distribution inversion inside the wellbore; The physical information neural network includes a multi-layer feedforward neural network, which consists of an input layer, a hidden layer, and an output layer. The input layer is the spatiotemporal solution domain of the physical problem. The hidden layer consists of several layers of neurons and nonlinear activation functions to approximate the physical model. By optimizing the weights and biases, the input is gradually mapped to the output. The output layer obtains the approximate distribution of the physical field in the spatiotemporal solution domain through the nonlinear mapping of multiple hidden layers. The output parameters of the physical information neural network are five physical quantities, including wellbore pressure. P Gas phase true velocity u g True velocity of liquid phase u l Gas phase volume fraction E g and liquid phase volume fraction E l ; The physical information neural network is shown below: (6)； (7)； In the formula, α 0 As the input layer of a physical information neural network, it consists of space x and time t Coordinate composition; α k For the first k Layer output; For the first k-1 Layer output; σ It is a non-linear activation function; w k and b k The first k Layer weights and bias coefficients; L This represents the number of layers in the physical information neural network.

2. The method for numerical solution of multiphase flow model in wellbore and inversion of gas-liquid distribution state according to claim 1, characterized in that, The loss function of a physical information neural network includes data loss. Loss Data Residual loss Loss Res Initial loss Loss IC and boundary loss Loss BC As shown in equations (8) to (12): (8)； (9)； (10)； (11)； (12)； In equations (8) to (12), , , , This is the weight of the loss function; if both the data and the equation have been dimensionless, then it is set to 1. n, m, p, and k represent the nth observation data point, the mth internal sampling point, the pth initial condition point, and the kth boundary point, respectively. N , M , P , K This represents the total number of data points, interior points, initial points, and boundary points. y (·) represents the output of the physical quantity predicted by the physical information neural network at (x, t); I (·) represents the initial physical quantity predicted by the physical information neural network; B (·) represents the boundary physical quantity output by the physical information neural network at the boundary. y n , I p , B k These represent the observed values ​​corresponding to the data points, the initial condition values ​​corresponding to the initial points, and the boundary condition values ​​corresponding to the boundary points. ， The data is known. e i (·) is the output value of the physical information neural network, which is then input into the first... i The result after the physical equation q This represents the total number of physical equations. The five physical equations for residual loss are shown in equations (13) to (17): (13)； (14)； (15)； (16)； (17)； In the formula, ρ g The density of free gas at the temperature and pressure within the annulus, in kg·m³. -3 ; ρ l Density of drilling fluid, kg·m -3 ; u g , u l The upward return velocities of free gas and drilling fluid, respectively, are in m / s. -1 ; E g , E l These are the integrals of free gas and drilling fluid, respectively, and are dimensionless. A Let m be the cross-sectional area of ​​the annulus. 2 , ,in D The inner diameter of the well shaft is in meters (m). d The outer diameter of the drill string for different well sections, in meters; q g The mass of gas produced per unit time per unit thickness of reservoir, in kg·s -1 ·m -1 ; R s m represents the solubility of the gas in the drilling fluid. 3 ·m -3 ; ρ gs The density of the gas under standard conditions is expressed in kg·m³. -3 ; B l The volume coefficient of the drilling fluid in the local area is dimensionless. u m The apparent flow rate of the mixed phase is in m·s. -1 , u m = u sl + u sg C0 is the gas phase distribution coefficient, which is dimensionless. u gr Let be the gas phase slip velocity, in m·s -1 ; P Pressure, Pa.

3. The method for numerical solution of multiphase flow model in wellbore and inversion of gas-liquid distribution state according to claim 1, characterized in that, Training a physical information neural network; including: This paper utilizes the built-in objects of the TensorFlow framework to build the network, design physical constraints, calculate loss, optimize and update, and perform graph execution. The weight matrix of the physical information neural network is generated using the Xavier initialization method. The parameters of the physical information neural network are optimized using the Adam-LBFGS hybrid optimizer. Initially, a fixed sampling point strategy is used in conjunction with the Adam optimizer for preliminary training. Subsequently, the optimization is performed by switching to the L-BFGS-B optimizer. An adaptive sampling mechanism based on residuals is introduced to dynamically increase training points in high residual regions, thereby improving the resolution of key regions.

4. The method for numerical solution of multiphase flow model in wellbore and inversion of gas-liquid distribution state according to claim 1, characterized in that, Design an adaptive activation function; including: The training process of the physical information neural network is as follows: find suitable weights w and bias terms b for each neuron in each layer of the physical information neural network so that the loss function gradually decreases and eventually reaches the global minimum; the update process of weights and bias terms is shown in equations (18) and (19): (18)； (19)； In the formula, η For learning rate, η >0, adopt dynamic learning rate, that is, after setting the initial learning rate, it decays once every certain number of training steps, so that the learning rate decreases in a step-like manner; J m For the first m The loss function at the next iteration; For the first m The weight vector of the physical information neural network at the next iteration. For the first m+ The weight vector of the physical information neural network in the first iteration; loss function Weights The gradient; For the first m The bias vector of the physical information neural network at the next iteration. For the first m The bias vector at +1 iterations loss function For bias b The gradient; In activation function Introducing hyperparameters c To obtain the output of the k-th layer As shown in equation (20): (20)； in, This refers to the output of the (k-1)th layer; It refers to the first k Layer bias vector; It refers to the first k Layer weight vector; c It is an adjustable parameter used to control the slope of the activation function, thereby changing the topology of the loss function during neural network training; c The gradient descent method is used to optimize the weights and biases together, as shown in equation (21): (21)； in, For the first m The slope of the activation function is controlled in the next iteration. For the first m+ In the first iteration, control the slope of the activation function; loss function For hyperparameters c The gradient; The adaptive activation function is the Tanh function.

5. The method for numerical solution of multiphase flow model in wellbore and inversion of gas-liquid distribution state according to claim 2, characterized in that, Design an adaptive sampling mechanism based on residuals; including: A probability density function based on residuals is used to increase the number of sampling points. As shown in equation (22): (22)； in, X For the candidate point set, by space x and time t Coordinate composition; r ( X ) represents the residual value corresponding to the candidate point set; This is a threshold used to filter out small residuals; when the residual value is less than... When the logarithmic function returns a negative value, the max function is used to convert the residual value to 0; Ω refers to the entire space-time domain. For the five physical equations shown in equations (13) to (17) e 1 - e 5 The corresponding residual loss is used to calculate its probability density function, thus giving each physical equation a different sampling point distribution. Specifically, in each round of sampling, new sampling points at different locations are added to these five physical equations to capture the residual distribution characteristics of each equation. The implementation process of the residual-based adaptive sampling mechanism is as follows: 1) Initialization: Determine the initial set of sampling points and set the relevant parameters for adaptive sampling, including the candidate point set, the number of iterations to start sampling, the number of new points added in each round of sampling, the total number of samplings, and the corresponding loop period. e 1 - e 5 The sampling priority and residual threshold parameter ε are determined; at the same time, a neural network for solving the PDE is constructed, and an appropriate number of hidden layers and neuron size are set. 2) Calculate the residuals and construct the probability density function: Combine automatic differentiation technology to calculate the residual loss at the candidate point, that is, the residuals corresponding to the 5 physical equations; at the same time, construct the probability density function corresponding to each physical equation according to equation (22); 3) Determine new sampling points: Using the five constructed probability density functions, extract five new sampling points from the candidate point set and add the five new sampling points to the training set; 4) Network training and optimization: On the updated training set, use the optimization algorithm to train the physical information neural network. After reaching the set number of training iterations, start a new round of sampling and repeat steps 2)-4) until the total loss is lower than the set threshold or the maximum number of iterations is reached, and obtain the adaptive physical information neural network. Adaptive physical information neural network structure optimization refers to: In the adaptive physical information neural network, there are 6 hidden layers and 60 neurons.

6. A numerical solution system for a wellbore multiphase flow model and a gas-liquid distribution state inversion system, used to implement the numerical solution system for a wellbore multiphase flow model and a gas-liquid distribution state inversion method as described in any one of claims 1-5, characterized in that, include: The physical information neural network construction module is configured to: construct and train a physical information neural network for dynamic simulation of multiphase flow in drilling wellbore and inversion of overflow gas distribution; determine the input and output of the physical information neural network; determine the loss function of the physical information neural network; and train the physical information neural network. The adaptive physical information neural network design module is configured to: design an adaptive optimization algorithm to optimize the final solution accuracy and convergence speed of the physical information neural network, and obtain the adaptive physical information neural network. Design an adaptive activation function; design an adaptive sampling mechanism based on residuals; The module for numerical solution of multiphase flow model in wellbore and inversion of gas-liquid distribution state is configured to: realize numerical solution of multiphase flow model in wellbore and inversion of gas-liquid distribution state based on adaptive physical information neural network.

Images