Pipeline pollution source tracing method based on coupling of pinn and isotope split flow field

Abstract

The application discloses a pipeline pollution source tracing method based on PINN and isotope diversion field coupling, belongs to the technical field of pipeline pollution source tracing based on data processing, and aims to solve the technical problem that existing pipeline micro-leakage points are difficult to trace, and the technical solution points are as follows: S1, constructing a multidimensional space-time training data set; S2, establishing an anisotropic reactive solute transport control equation group; including a concentration conservation equation and an isotope mass conservation equation, and introducing an isotope mass density intermediate variable in the isotope mass conservation equation; S3, constructing a full-connection deep neural network; S4, constructing a hybrid loss function containing physical residual and sparse regularization; S5, parameter iteration and source term reconstruction based on two-stage optimization; after training, the distribution of the output layer is extracted, a heat map is drawn based on the distribution, and the pipeline leakage point is determined based on the heat map.

Description

Technical Field

[0001] This invention belongs to the field of pipeline pollution source tracing technology based on data processing, specifically involving a pipeline pollution source tracing method based on the coupling of PINN and isotope flow field. Background Technology

[0002] Long-distance oil pipelines, as core infrastructure for petroleum energy transportation, are mostly buried several meters underground, traversing complex terrains such as farmland, forests, and wetlands. They are constantly exposed to complex environments including soil erosion, groundwater soaking, and ground subsidence. Furthermore, they are highly susceptible to hidden micro-leakage (such as pinhole leaks and micro-crack leaks) due to factors such as weld defects, pipe aging, and damage to anti-corrosion layers during pipeline laying, as well as disturbances from third-party construction and media erosion during operation. These leaks are characterized by small leakage volumes (typically milliliters / minute to liters / hour), extremely low flow rates, slow diffusion, and no significant instantaneous pressure fluctuations. They are also covered by soil and vegetation, making them difficult to detect visually. Long-term leakage not only wastes petroleum resources but also pollutes groundwater and soil, causing ecological and environmental accidents, and even posing safety hazards such as fires and explosions. However, current detection technologies for such hidden micro-leakage still have many significant shortcomings, making it difficult to meet the needs of accurate detection, tracing, and prevention in practical engineering projects.

[0003] Existing physical detection methods (such as negative pressure wave method, fiber optic sensing method, and acoustic wave detection method) all take "significant physical signal changes" as the core of detection. They have a good response effect to instantaneous large-flow leakage (such as pipe rupture and large-area damage), but they have obvious detection blind spots for chronic micro-leakage (such as pinhole leakage) with extremely low flow rate and insignificant pressure changes, and cannot be effectively identified.

[0004] Among them, the negative pressure wave method is a commonly used method for detecting leaks in long-distance oil pipelines. Its core relies on the propagation of negative pressure waves generated by the instantaneous pressure drop at the leak point. However, when a small leak (such as a pinhole leak) occurs, the amount of medium leakage is extremely small and cannot form a significant pressure change in the pipeline. The negative pressure wave signal is weak and dispersed, far below the detection threshold of existing pressure sensors. Even with signal amplification and filtering, it is difficult to extract effective signals from pressure fluctuations during normal pipeline operation (such as changes in medium flow rate, slight adjustments of pumps and valves, and formation pressure disturbances), resulting in leaks going undetected for a long time.

[0005] While fiber optic sensing methods (such as distributed fiber optic temperature and vibration measurements) can achieve long-distance, blind-zone-free coverage, they are limited by detection accuracy and environmental interference. The soil temperature changes and vibration signals caused by minor leaks are extremely weak and easily masked by interference signals such as ambient temperature fluctuations, soil settlement, and plant root growth. Especially in complex environments such as cold regions and wetlands, the environmental background signal superimposed on the leak signal further reduces the detection sensitivity. At the same time, fiber optic sensing systems have high requirements for installation technology. Construction disturbances around the pipeline and uneven soil compaction can all lead to damage to the sensing fiber or signal distortion, further expanding the detection blind zone.

[0006] In addition, other physical detection methods (such as ultrasonic testing and cathodic protection potential testing) also have similar problems: ultrasonic testing is affected by soil attenuation and cannot penetrate thick soil layers to capture minute leak signals; cathodic protection potential testing can only reflect the damage to the pipeline's anti-corrosion layer, cannot determine whether there is an actual leak at the damaged point, and is difficult to locate the specific location of minute damage points, thus failing to achieve accurate leak identification.

[0007] When a minor leak occurs in a long-distance oil pipeline, the leaked medium (oil and associated contaminants) gradually seeps into the ground and spreads through groundwater flow. In such cases, groundwater source tracing models are needed to pinpoint the leak location and provide a basis for leak response. However, existing groundwater source tracing models are mostly based on the ideal assumption of "homogeneous and isotropic media," which differs significantly from actual pipeline laying scenarios. This leads to source tracing results that deviate severely from the actual leak point, often biased towards the downstream area, and thus unable to provide precise guidance for leak response.

[0008] The core issue lies in the fact that the "dominant flow channel" formed during pipeline laying was not considered in the model. When laying long-distance oil pipelines, trenches 1.5-3 meters deep need to be excavated. This excavation process disrupts the original soil structure, resulting in backfill soil (mostly disturbed soil or mixed soil) having significantly higher porosity and permeability than the surrounding undisturbed soil (especially in areas with cohesive or silty soil). This creates a "dominant flow channel" extending longitudinally along the trench. When a minor leak occurs, the leaking medium does not diffuse laterally to the surrounding soil as the model assumes. Instead, it preferentially "flows" longitudinally along the more permeable backfill soil in the trench. Only after flowing a certain distance does it gradually diffuse laterally to the surrounding undisturbed soil.

[0009] Traditional hydrological models neglect this crucial effect, still calculating pollutant diffusion paths and velocities as if the medium were homogeneous. This leads to significant discrepancies between model-predicted pollutant diffusion range and direction and actual conditions. Firstly, the models fail to capture the longitudinal flow of pollutants, mistakenly identifying the starting point of diffusion after longitudinal flow as the leak point, resulting in source tracing results biased towards the downstream (flow direction). Secondly, the models' calculations of pollutant diffusion rates have significant errors, making it impossible to accurately determine the duration and extent of the leak, further increasing the difficulty of leak control and potentially leading to a continued expansion of the pollution area and more severe ecological damage. Furthermore, the permeability of the backfill soil in the trench is affected by backfilling techniques, soil type, and compaction degree, exhibiting significant spatial heterogeneity. Traditional models cannot quantify this heterogeneity, further exacerbating source tracing errors.

[0010] Traditional numerical simulation methods (such as MODFLOW and MT3DMS) are core tools for tracing groundwater pollution sources and predicting leakage and diffusion. The accuracy of their calculation results highly depends on precise boundary conditions (such as groundwater level, flow direction, and recharge rate) and hydrogeological parameters (such as soil permeability, porosity, and dispersion coefficient). However, in actual long-distance oil pipeline projects, due to the long pipeline length (usually hundreds of kilometers), complex terrain crossings, and extremely limited number of monitoring wells with incomplete monitoring coverage, key parameters cannot be accurately obtained, resulting in high uncertainty in parameter inversion and thus affecting the reliability of numerical simulation results.

[0011] On the one hand, boundary conditions are difficult to determine accurately: long-distance oil pipelines often pass through remote areas with weak groundwater monitoring infrastructure, a small number of monitoring wells, and sparse distribution. It is impossible to obtain boundary conditions such as groundwater level, flow direction, and recharge / discharge along the pipeline in real time and comprehensively. Only limited monitoring data can be used for interpolation estimation, which results in large estimation errors. At the same time, groundwater level and flow direction are affected by seasonal changes, precipitation, and human activities (such as agricultural irrigation and groundwater extraction), and have strong dynamics. Traditional numerical simulations cannot update boundary conditions in real time, which further reduces the simulation accuracy.

[0012] On the other hand, hydrogeological parameters exhibit significant spatial heterogeneity and uncertainty: soil types and geological structures vary considerably along the pipeline route, with soil permeability coefficients and porosity differing by one to two orders of magnitude in different regions. The limited number of monitoring wells cannot fully cover this spatial heterogeneity, and calculations can only be performed using "average parameters" or "zonal parameters," which cannot accurately reflect the hydrogeological characteristics of local areas. Furthermore, some hydrogeological parameters (such as dispersion coefficients) cannot be directly obtained through on-site monitoring and can only be estimated using empirical formulas. The estimated results deviate significantly from the actual values, leading to strong dispersion in parameter inversion results. Numerical simulations cannot accurately predict the diffusion path of pollutants and the location of leak points, making it difficult to meet the precise source tracing requirements of actual engineering projects.

[0013] In recent years, data-driven methods (especially deep learning) have been applied to some extent in the field of leak detection. Their core is to train models using a large amount of historical leak data to achieve automatic identification of leak signals and precise location of leak points. However, in the field of detecting hidden, minute leaks in long-distance oil pipelines, purely data-driven methods face significant application bottlenecks and are difficult to implement and promote. The core reason is the "data scarcity" in the field of environmental accidents, which cannot meet the training requirements of deep learning models.

[0014] First, hidden micro-leaks are low-probability events, and historical leakage data is extremely scarce: the safety management requirements for long-distance oil pipelines are extremely high, and pipeline operators will take a variety of prevention and control measures to avoid leakage accidents, especially micro-leaks. Due to their strong concealment and slow manifestation of harm, the probability of detection is extremely low, resulting in very little amount of historical leakage data (including leakage signals, leakage locations, leakage amounts, environmental parameters, etc.) that can be used for model training. This makes it impossible to cover leakage scenarios with different soil types, different leakage amounts, and different environmental conditions, resulting in insufficient model training and extremely poor generalization ability.

[0015] Secondly, labeling leakage data is difficult and costly: the signals from minute leaks are weak and easily interfered with, requiring professional technicians to combine on-site detection, laboratory analysis, and other methods to accurately label the characteristics of the leakage signal, the location of the leak point, and other key information. The labeling process is time-consuming, labor-intensive, and extremely costly. At the same time, the leakage data formats under different detection equipment and different detection scenarios are not uniform, requiring a lot of data preprocessing (such as data cleaning, standardization, and feature extraction), which further increases the difficulty of data preparation.

[0016] Furthermore, purely data-driven models lack physical mechanism support and are not robust enough: the leakage process of long-distance oil pipelines is affected by a variety of factors such as soil environment, groundwater flow, and pipeline operating parameters, and has complex physical and chemical mechanisms. However, simple deep learning models rely only on statistical data patterns and cannot incorporate these physical mechanisms. When encountering scenarios that do not appear in the training data (such as minor leaks under special soil types or extreme weather conditions), the detection accuracy and positioning accuracy of the model will drop sharply, failing to meet the reliability requirements of actual engineering.

[0017] In summary, existing detection technologies have insurmountable shortcomings in the identification, tracing, and simulation of concealed micro-leaks, failing to meet the actual needs of safety control in long-distance oil pipelines. There is an urgent need to develop new detection technologies that integrate physical mechanisms and data-driven approaches, balancing detection accuracy and practicality, to address the pain points of "difficult to detect, difficult to locate, and difficult to trace" micro-leaks. Summary of the Invention

[0018] The purpose of this invention is to provide: A pipeline pollution source tracing method based on PINN and isotope fractionation field coupling is proposed to solve the technical problem of difficulty in tracing the source of small pipeline leaks. This method does not rely on a large amount of training data, but embeds the isotope Rayleigh fractionation law and the hydrodynamic equation of the dominant flow in the trench as physical constraints into the loss function of the neural network. By minimizing the physical residual, the location of the leak source is intelligently inverted.

[0019] Terminology Explanation: Unless otherwise defined, all technical terms used herein have the same meanings as commonly understood by one of ordinary skill in the art to which this subject matter pertains. Unless otherwise stated, all patents, patent inventions, and disclosures cited throughout this document are incorporated herein by reference in their entirety. Where multiple definitions exist for terms herein, the definitions provided in this chapter shall prevail.

[0020] It should be understood that the above brief description and the following detailed description are exemplary and for illustrative purposes only, and do not limit the subject matter of the invention in any way. In this invention, the singular is used in conjunction with the plural unless otherwise specifically stated. It should also be noted that, unless otherwise stated, the use of “or” or “or” means “and / or”. Furthermore, the use of the term “comprising” and other forms such as “including,” “containing,” and “contains” are not limiting.

[0021] Unless specifically defined herein, the use of various commercially available products herein employs standard techniques. For example, they may be implemented using the manufacturer's instructions for use, or in accordance with methods known in the art or the description of this invention. The techniques and methods described herein can generally be implemented according to conventional methods well known in the art, based on the descriptions in the various general and more specific documents cited and discussed in this specification.

[0022] The terms “optional / arbitrary” or “optionally / arbitrarily” mean that the event or situation described below may or may not occur, including both the occurrence and non-occurrence of the event or situation.

[0023] The term "monomer isotope ratio" as used in this article refers to the determination of the elemental isotope abundance ratio at the molecular level of a single compound, and expressed as... This indicates a difference relative to the standard.

[0024] In a first aspect, the present invention provides: A pipeline pollution source tracing method based on PINN and isotopic flow field coupling includes: S1. Construct a multidimensional spatiotemporal training dataset; Collect data from sparse monitoring wells around the pipeline, including pollutant concentrations. isotope ratio with monomer It also records the coordinates of the monitoring wells and the sampling time to form a spatiotemporal coordinate system. As a training set ; S2. Establish a set of governing equations for anisotropic reactive solute transport; It includes the concentration conservation equation and the isotope mass conservation equation, and introduces the isotope mass density as an intermediate variable into the isotope mass conservation equation. ; S3. Construct a fully connected deep neural network; Constructing based on spatiotemporal coordinates As input, with The output deep neural network, where, This indicates the output pollutant concentration. This represents an intermediate variable representing the output isotopic mass density. The leakage source phase field to be inverted; S4. Construct a hybrid loss function that includes physical residuals and sparse regularization; Constructing the loss function loss function Due to observation data error Physical equation residuals Source-phase sparsity constraint Composition, including the physical equation residuals It was calculated using the concentration conservation equation and the isotope mass conservation equation; S5. Parameter iteration and source term reconstruction based on two-stage optimization; By minimizing Update network weights and extract the output layer after training. Distribution, based on Draw a heat map; determine the pipeline leak point based on the heat map.

[0025] Based on further solutions to the technical problems of the present invention, or simultaneous solutions to multiple technical problems, the preferred solution in the technical solution provided in the first aspect of the present invention includes: First preferred option: Step S1, when constructing the training set, also includes normalizing the monitoring well data.

[0026] By normalizing the data, deep neural networks can perform iterative training with simplified data while maintaining accuracy, thus improving training speed.

[0027] Second preferred option: The concentration conservation equation is:

[0028] The isotopic mass conservation equation is:

[0029] Where C represents the pollutant in spatial coordinates and molar concentration at time t; I represents the intermediate variable of isotopic mass density; ; Indicates the isotope ratio of a single monomer; , This represents the components of the velocity vector v in the x and y directions; , These represent the hydrodynamic dispersion coefficients along the x-axis and y-axis, respectively; The rate constant representing the first-order biodegradation rate of a pollutant; Represents the Rayleigh fractionation coefficient; S represents the intensity of the leakage source phase to be inverted; This indicates the initial isotope ratio of the crude oil in the pipeline.

[0030] By introducing intermediate variables, numerical oscillations caused by the neural network directly predicting isotope ratios (which are negative and not conserved) are avoided.

[0031] Third preferred option: The components of the velocity vector are calculated based on Darcy's law as follows:

[0032]

[0033] in, Indicates water head. Indicates effective porosity. This represents the permeability coefficient at various points in space, expressed in m / d.

[0034] Fourth preferred option: Permeability coefficient at various points in space The calculation is based on the following formula: Define the permeability field Assume the centerline of the pipeline is located at The half-width of the trench area is Then the permeability coefficient distribution function is:

[0035] in: This represents the background permeability coefficient of the undisturbed soil surrounding the trench; This indicates the permeability coefficient of the backfill soil inside the trench; This indicates the coordinate position of the oil pipeline centerline on the y-axis; Indicates the smoothing factor; This represents the hyperbolic tangent function.

[0036] Fifth preferred option: The activation function of deep neural networks is the hyperbolic tangent function; The output layer is constrained, specifically including: concentration output. and source term output Applying the Softplus activation function to isotope intermediate variables No non-negativity constraints are applied.

[0037] By imposing constraints on the output layer, it can be ensured that the source terms of the output are not negative, and the output results conform to physical reality.

[0038] Sixth preferred option: loss function for:

[0039] in, , and These are the weights for each sub-item.

[0040] Among them, the error of the observation data for:

[0041] in, This indicates the total number of observation data points (monitoring well sampling points); This represents the predicted concentration output by the neural network at the i-th observation point; This represents the predicted isotope quantity output by the neural network at the i-th observation point; This represents the measured pollutant concentration at the i-th observation point; This represents the measured isotope ratio at the i-th observation point; This represents the L2 norm (mean squared error), used to measure the distance between the predicted value and the true value.

[0042] Physical equation residuals for:

[0043] in, This represents the total number of configuration points randomly sampled within the computational domain; This represents the algebraic residual calculated after substituting the network predictions into the concentration conservation equation; This represents the algebraic residual value calculated after substituting the network prediction values ​​into the isotopic mass conservation equation.

[0044] Source phase sparsity constraint for:

[0045] This indicates the number of grid points used to evaluate the distribution of source terms; This represents the predicted source term strength output by the network at grid point k; This represents the L1 norm.

[0046] Seventh preferred option: Step S5 specifically includes: S51, Global Optimization: Using the Adam optimizer, setting the initial learning rate and momentum parameters, and iteratively training. Second, reduce data loss ; S52. Physical Consistency Convergence: Switch to the L-BFGS-B optimizer and perform second-order optimization until the physical equation residuals are found. < , The preset convergence threshold is used; S53. Source Term Reconstruction: Extract the global source term output of the trained grid at time t=0. ,based on Draw a heat map; set a peak intensity threshold, select bright spots in the heat map whose peak intensity is greater than the peak intensity threshold, obtain the centroid coordinates of the bright spots, and then determine the leak point.

[0047] The present invention has at least the following beneficial effects: This invention introduces an isotope fractionation field coupling mechanism, which incorporates the residuals obtained from the concentration conservation equation and the isotope mass conservation equation when calculating the loss function. This enables the precise identification of the specific leak point and significantly improves the accuracy of leak point tracing.

[0048] This invention also imposes physical constraints on the deep neural network to ensure that the source term of the output is not negative, thereby solving the ill-posedness of the inversion problem and further improving the accuracy of leak point location. Attached Figure Description

[0049] Figure 1 This is a flowchart of the present invention. Detailed Implementation

[0050] The following non-limiting embodiments are intended to enable those skilled in the art to gain a more comprehensive understanding of the present invention, but do not limit the invention in any way. The following content is merely an exemplary description of the scope of protection claimed by the present invention, and those skilled in the art can make various changes and modifications to the present invention based on the disclosed content, and such changes should also fall within the scope of protection claimed by the present invention.

[0051] The present invention will be further described below by way of specific embodiments. Unless otherwise specified, all instruments, devices, equipment and other items used in the embodiments of the present invention are obtained through conventional commercial means.

[0052] Example 1 This embodiment provides a pipeline pollution source tracing method based on the coupling of PINN and isotopic flow field, such as... Figure 1 As shown, the specific steps include: S1. Construct a multidimensional spatiotemporal training dataset; Collect data from sparse monitoring wells around the pipeline, including pollutant concentrations. isotope ratio with monomer For input spatiotemporal coordinates The corresponding monitoring well data were subjected to dimensionless normalization, where x is the axial coordinate along the pipeline, y is the lateral coordinate perpendicular to the pipeline direction, and t represents the sampling time. The spatiotemporal coordinates and the corresponding normalized monitoring well data were correlated to form a training set. Normalization simplifies the data and increases training speed.

[0053] S2. Establish a set of governing equations for anisotropic reactive solute transport; Based on Darcy's law and Rayleigh fractionation principle, a system of partial differential equations (PDEs) describing the transport of pollutants in heterogeneous media in pipe trenches is constructed, and intermediate variables are introduced. To eliminate the computational singularity of the isotope ratio nonlinear term; The governing equations for anisotropic reactive solute transport are specifically expressed in the following coupled partial differential equation form: Concentration conservation equation:

[0054] Isotope mass conservation equation:

[0055] Where C represents the pollutant in spatial coordinates and the molar concentration at time t, i.e. Unit: mol / L; I represents the intermediate variable of isotopic mass density; This represents the isotope ratio of the monomer, i.e. The unit is ‰; this variable is introduced to transform the highly nonlinear isotope fractionation equation into a linear transport equation. t represents the sampling time, in days; Represents spatial Cartesian coordinates, where x is the axial coordinate along the pipe and y is the lateral coordinate perpendicular to the pipe's direction, both in meters (m). , The components of the velocity vector v in the x and y directions represent the dominant flow in the trench, and the units are m / d. , These represent the hydrodynamic dispersion coefficients along the x-axis and y-axis, respectively; For hydrodynamic dispersion coefficient and The hydrodynamic dispersion coefficient characterizes the macroscopic migration ability of solutes in groundwater along the coordinate axis under hydrodynamic action, and its value is closely related to geological conditions. The hydrodynamic dispersion coefficient is obtained through the following method: Based on the on-site hydrogeological survey, the longitudinal dispersion along the pipeline axis (x direction) was determined. Lateral dispersion perpendicular to the pipe orientation (y-direction) ; typically longitudinal dispersion The value ranges from 0.1m to 100m, and the lateral dispersion... It is approximately one-tenth to one-fifth of the longitudinal dispersion αL.

[0056] Based on the acquired velocity vector components and The hydrodynamic dispersion coefficients in each direction are calculated according to the following formula:

[0057]

[0058] Longitudinal diffusion and lateral dispersion The value can be obtained from field tracer tests, or from empirical values ​​in regional hydrogeological manuals, or from recommended parameter values ​​in publicly available literature for similar soil types. When field testing is lacking, longitudinal dispersion... Empirical formulas can be used to estimate based on the observation scale: =0.1×L, where L is the pollutant migration distance in meters. The typical geological conditions applicable to this invention are: the pipeline is buried in an alluvial plain or lacustrine sedimentary strata; the backfill soil of the trench is mainly silty sand or sandy silt; and the undisturbed soil around the trench is mainly silty clay or clay. In this embodiment, the longitudinal dispersion... The value is 5.0m, representing the lateral dispersion. The value is 0.5m.

[0059] The rate constant for the first-order biodegradation of pollutants is expressed in units of 1 / d, which characterizes the rate at which pollutants naturally decay. α represents the Rayleigh fractionation coefficient, a dimensionless parameter ranging from 0.950 to 1.000. It is used to quantify the degree of separation between light and heavy isotopes during biodegradation and is a core parameter of the coupling equation. When α is less than 1, it indicates that microorganisms preferentially utilize light isotopes (carbon-12), leading to the gradual enrichment of heavy isotopes (carbon-13) in residual pollutants; when α approaches 1, it indicates that the isotope fractionation effect is extremely weak.

[0060] The Rayleigh fractionation coefficient is obtained by at least one of the following methods: The first method involves conducting batch experiments on microbial degradation in the laboratory. Specifically, groundwater or soil samples are collected from the pipeline leak area. The in-situ biodegradation process is simulated under controlled laboratory conditions. The concentration of residual pollutants and the monomer isotope ratio are monitored regularly. The Rayleigh equation is used for linear regression fitting, and the slope of the fitted line is calculated to obtain the α value. The second method involves consulting publicly available databases of isotopic fractionation factors for similar pollutants or literature reports. For typical petroleum hydrocarbon pollutants, including crude oil, diesel, and gasoline, the Rayleigh fractionation coefficient α generally ranges from 0.990 to 0.999, with the α value corresponding to crude oil biodegradation typically ranging from 0.998 to 0.999. For chlorinated hydrocarbon pollutants, such as trichloroethylene and tetrachloroethylene, the α value generally ranges from 0.950 to 0.990. For benzene compounds, such as benzene and toluene, the α value generally ranges from 0.997 to 0.999. The corresponding typical value range can be selected based on the actual pollutant type, or the fractionation coefficient for a specific site can be accurately obtained through the aforementioned indoor batch experiments.

[0061] S is the intensity of the leakage source phase to be inverted, in mol / (L·d), which is positive at the leakage point and 0 in the non-leakage region; This represents the initial isotope ratio of the crude oil in the pipeline, expressed in per mille (‰), and is a known measured constant.

[0062] The components of the velocity vector are calculated based on Darcy's law as follows:

[0063]

[0064] in, Indicates water head. Indicates effective porosity. This represents the permeability coefficient at various points in space, expressed in m / d.

[0065] And the permeability coefficient at various points in space The calculation is based on the following formula: Define the permeability field Assume the centerline of the pipeline is located at The half width of the trench is Then the permeability coefficient distribution function is:

[0066] This function introduces a "smooth step function" to characterize the abrupt change in permeability between the trench and the original soil. In this formula: This represents the background permeability coefficient of the undisturbed soil surrounding the trench, in m / d. This represents the permeability coefficient of the backfill soil inside the trench, expressed in m / d, and typically meets the following requirements. > ; This indicates the coordinate position of the oil pipeline centerline on the y-axis, in meters (m). This represents half the width of the trench area, in meters (m). This represents the smoothing factor, a dimensionless factor used to control the steepness of the change in permeability coefficient from the trench to the original soil. The smoothing factor... The principle that ensures the differentiability of the permeability coefficient distribution function K(x,y) lies in the hyperbolic tangent function. It has continuous derivatives of arbitrary order (i.e., it belongs to the infinite smooth function category of C), and will Introducing the independent variable tanh (i.e., forming the absolute value of the difference of the ordinates within the parentheses of tanh minus half the width and then divided by σ) does not disrupt the smoothness of the function, but only changes its rate of change; when When K(x,y) takes a positive value, it achieves a smooth transition rather than an abrupt change near the trench boundary, ensuring the continuous differentiability of the permeability coefficient field throughout the entire computational domain, thus satisfying the mathematical requirement for the differentiability of the coefficient field when solving partial differential equations based on the Physical Information Neural Network (PINN).

[0067] The smoothing factor σ ranges from 0.1 to 2.0, preferably from 0.2 to 1.0. The smaller the value of σ, the steeper the change in permeability coefficient at the trench boundary, approaching a step function but still remaining differentiable; the larger the value of σ, the smoother the transition. In practical applications, the value of σ is determined based on the actual transition between the backfill and undisturbed soil in the trench: when the trench boundary is clear and the boundary between the backfill and undisturbed soil is obvious, σ takes a smaller value between 0.1 and 0.4; when the trench boundary is blurred and there is a wide disturbance transition zone, σ takes a larger value between 0.8 and 2.0. The above value range and selection principles also apply when using other families of smoothing functions to construct the permeability coefficient transition. As long as the family of functions used itself is continuous and differentiable, a positive value for σ can guarantee the differentiability of the transition function. This represents the hyperbolic tangent function, used to construct continuously differentiable step changes.

[0068] S3. Construct a fully connected deep neural network (PINN). Build with As input, with The output deep neural network, where, The input is the phase field of the leakage source to be inverted. Then, based on the correspondence, the deep neural network can obtain the pollutant concentration and monomer isotope ratio corresponding to the spatiotemporal coordinate, and then perform data prediction based on the deep neural network.

[0069] The specific architecture of a deep neural network includes the design of specific activation functions and output constraints: Hidden layer activation function: Hyperbolic tangent function. Because it possesses C∞ smoothness, it ensures the second derivative It exists and is non-zero; Output layer constraint: concentration output and source term output Apply the Softplus activation function: In the Softplus activation function, the original model output is transformed by exponential and natural logarithms to obtain an output that is always positive, ensuring the non-negativity of physical quantities. This is also relevant for intermediate isotope variables. No nonnegativity constraints are imposed, allowing it to adapt to numerical fluctuations.

[0070] S4. Construct a hybrid loss function that includes physical residuals and sparse regularization; The loss function is constructed by calculating the partial derivatives of the network output with respect to the input parameters using automatic differentiation techniques. This function is determined by the error in the observed data. Physical equation residuals Source-phase sparsity constraint composition.

[0071] Specifically, the loss function for:

[0072] in, Right now Θ is the set of trainable parameters of the neural network. , and These are the weights of each sub-item; In this embodiment, regarding the weights: It is 0.2. It is 0.3. The value is 0.5; this value is suitable for tracing the source of leakage in typical long-distance oil pipelines, where the number of monitoring wells is sparse (2 to 4), and the leakage is a single point of minute seepage.

[0073] When monitoring wells are densely deployed, sampling frequency is high, and the reliability of measured data is strong, one can choose... Larger value Smaller values, for example It is 0.5. It is 0.3. It is 0.2.

[0074] If the number of monitoring wells is extremely small, such as less than two, or if the quality of the monitoring well data is poor (hydrogeological parameters such as sparse permeability and dispersion coefficients are obtained through exploration, which yields higher accuracy), then it is advisable to... The value of is increased, for example It is 0.1. It is 0.6. It is 0.3.

[0075] When a single strong leak point is known and precise geometric location is required, or when multiple sources of interference need to be eliminated, it can improve... The possible values, for example It is 0.2. It is 0.2. It is 0.6.

[0076] Observational data error for:

[0077] in, This indicates the total number of observation data points (monitoring well sampling points); This represents the predicted concentration output by the neural network at the i-th observation point; This represents the predicted isotope quantity output by the neural network at the i-th observation point; This represents the measured pollutant concentration at the i-th observation point; This represents the measured isotope ratio at the i-th observation point; This represents the L2 norm (mean squared error), used to measure the distance between the predicted value and the true value.

[0078] Physical equation residuals for:

[0079] in, This represents the total number of configuration points randomly sampled within the computational domain; This represents the algebraic residual calculated by substituting the network predictions into the left-hand side of the concentration conservation equation; This represents the algebraic residual calculated by substituting the network predictions into the left-hand side of the isotope mass conservation equation.

[0080] Source phase sparsity constraint for:

[0081] This indicates the number of grid points used to evaluate the distribution of source terms; This represents the predicted source term strength output by the network at grid point k; This represents the L1 norm (sum of absolute values), a mathematical operator that induces sparsity, forcing most of the non-leaking regions to become sparsity. The value is 0.

[0082] S5. Parameter iteration and source term reconstruction based on two-stage optimization; By minimizing Update network weights and extract the output layer after training. Distribution, using thresholds to determine and locate pipeline leak points.

[0083] Specifically, it includes: S51, Global Optimization: Using the Adam optimizer, set the initial learning rate. Momentum parameters , Iterative training Secondly, based on the Adam optimizer's own adaptive moment estimation optimization algorithm, data loss is quickly reduced. ; where the initial learning rate Used to control the step size of parameter updates.

[0084] S52, Physical Consistency Convergence: Switch to the L-BFGS-B (Finite Memory Quasi-Newton) optimizer, and use the inverse matrix approximation information of the Hessian matrix for second-order optimization until the physical equation residuals are found. < , To preset the convergence threshold, we need to ensure that the inversion results strictly satisfy the physical and chemical conservation laws.

[0085] Preset convergence threshold Values to .

[0086] S53. Source Term Reconstruction: Extract the global source term output of the trained grid at time t=0. ,based on Draw a heat map; based on the heat map, identify the light spot with a large peak intensity, and then determine the leak point.

[0087] Example This embodiment uses the source tracing and inversion of a minor leak in a buried oil pipeline in a plain area of ​​China as an example. The pipeline is buried at a depth of 1.5m and has a diameter of 720mm. There is a significant difference in permeability between the backfill soil of the trench and the surrounding undisturbed soil. This embodiment aims to accurately invert the location and intensity of the unknown upstream leak point using the method described in this invention, based on the pollutant concentration and monomer isotope data collected from sparse downstream monitoring wells.

[0088] The computing environment in this embodiment is: Intel Xeon Gold 6248R CPU, NVIDIA RTX 3090 GPU (24GB VRAM), developed based on the TensorFlow 2.4.0 deep learning framework.

[0089] This invention adopts a "physical constraint-driven data flow" model, and its core process is as follows: Input layer: Sparse monitoring data , +Physical parameters , , , ; Preprocessing layer: Spatiotemporal normalization + variable decoupling ; Network layer: PINN deep network (8 layers × 64 neurons, hyperbolic tangent activation); Output prediction field: concentration ,isotope , source item ; Physical constraint layer: Calculate the first / second order reciprocals; substitute them into the physical equations for the residuals. Calculate the residuals; Optimization layer: Hybrid loss function ; Decision-making level: Extraction Generate a heatmap of leakage probability and output the coordinates of the leakage point.

[0090] The specific steps include: S1. Construct a multidimensional spatiotemporal training dataset Computational Domain Definition and Meshing: The two-dimensional simulation region Ω is defined as a rectangular region with a length of 500m along the pipe flow direction (x-axis) and a width of 100m perpendicular to the pipe direction (y-axis). The pipeline centerline is located at Although this method is a mesh-free method, a background virtual mesh Nx×Ny=256×64 is constructed for visualization and partial integral calculation.

[0091] Data Acquisition and Variable Transformation: Two monitoring wells were deployed downstream at x=350m and x=450m respectively. Samples were collected at t=50, 100, and 150 days. Sample Data Example: Monitoring well A (x=350, y=0.5) measured the total petroleum hydrocarbon concentration at t=100 days. =15.2 mg / L, isotopic value = 28.5‰.

[0092] Additionally, for each spatiotemporal coordinate The monitoring well data at the location were subjected to dimensionless normalization, including pollutant concentration and monomer isotope ratios. The normalization formula is as follows:

[0093]

[0094] in, This represents the normalized value of the pollutant concentration. This represents the minimum pollutant concentration value. This indicates the highest pollutant concentration value; This represents the normalized value of the monomer isotope ratio. This represents the smallest monomeric isotope ratio. This represents the maximum monomeric isotope ratio. This allows mapping the physical domain to... Finally, the training set was obtained. .

[0095] S2. Establish the governing equations for anisotropic reactive solute transport. First, a smooth anisotropic permeability field is constructed. To address the issue of discontinuous derivatives caused by abrupt changes in the medium in traditional finite element methods, this embodiment employs a hyperbolic tangent smoothing operator to construct the permeability field. .

[0096] Parameter settings: Permeability coefficient of backfill soil inside the trench =40m / d, background permeability coefficient of undisturbed soil =2m / d, half width of the trench =1.2m, smoothing factor It is 0.4.

[0097] For any point within the computational domain Its permeability coefficient is:

[0098] This formula generates a value in y∈[ The value within the interval [1.2, 1.2] is approximately 40, and it is a continuous function that smoothly decays to 2 on both sides.

[0099] Flow field and It can be calculated using the following formula:

[0100]

[0101] Among them, effective porosity The hydraulic gradient is 0.3. It is 0.005.

[0102] Then, variable decoupling is performed. To avoid numerical oscillations caused by the neural network directly predicting the isotope ratio (negative and non-conservative), this embodiment introduces an intermediate variable I of isotope mass density. Example: =15.2×( 28.5)= 433.2, will ( , ) serves as the ground truth for neural network training.

[0103] Then, a set of governing equations for anisotropic reactive solute transport was constructed: Concentration conservation equation:

[0104] Isotope mass conservation equation:

[0105] Where C represents the pollutant in spatial coordinates and the molar concentration at time t, i.e. Unit: mol / L; Indicates an intermediate variable for isotopic mass density; This represents the isotope ratio of the monomer, i.e. The unit is ‰; this variable is introduced to transform the highly nonlinear isotope fractionation equation into a linear transport equation. t represents the sampling time, in days; Represents spatial Cartesian coordinates, where x is the axial coordinate along the pipe and y is the lateral coordinate perpendicular to the pipe's direction, both in meters (m). , The components of the velocity vector v in the x and y directions represent the dominant flow in the trench, and the units are m / d. , These represent the hydrodynamic dispersion coefficients along the x-axis and y-axis, respectively; The rate constant for the first-order biodegradation of pollutants is expressed in units of 1 / d, which characterizes the rate at which pollutants naturally decay. α represents the Rayleigh fractionation coefficient, a dimensionless parameter ranging from 0.950 to 1.000. It is used to quantify the degree of separation between light and heavy isotopes during biodegradation and is a core parameter of the coupling equation. When α is less than 1, it indicates that microorganisms preferentially utilize light isotopes (carbon-12), leading to the gradual enrichment of heavy isotopes (carbon-13) in residual pollutants; when α approaches 1, it indicates that the isotope fractionation effect is extremely weak.

[0106] The Rayleigh fractionation coefficient is obtained by at least one of the following methods: The first method involves conducting batch experiments on microbial degradation in the laboratory. Specifically, groundwater or soil samples are collected from the pipeline leak area. The in-situ biodegradation process is simulated under controlled laboratory conditions. The concentration of residual pollutants and the monomer isotope ratio are monitored regularly. The Rayleigh equation is used for linear regression fitting, and the slope of the fitted line is calculated to obtain the α value. The second method involves consulting publicly available databases of isotopic fractionation factors for similar pollutants or literature reports. For typical petroleum hydrocarbon pollutants, including crude oil, diesel, and gasoline, the Rayleigh fractionation coefficient α generally ranges from 0.990 to 0.999, with the α value corresponding to crude oil biodegradation typically ranging from 0.998 to 0.999. For chlorinated hydrocarbon pollutants, such as trichloroethylene and tetrachloroethylene, the α value generally ranges from 0.950 to 0.990. For benzene compounds, such as benzene and toluene, the α value generally ranges from 0.997 to 0.999. The corresponding typical value range can be selected based on the actual pollutant type, or the fractionation coefficient for a specific site can be accurately obtained through the aforementioned indoor batch experiments.

[0107] S is the intensity of the leakage source phase to be inverted, in mol / (L·d), which is positive at the leakage point and 0 in the non-leakage region; This represents the initial isotope ratio of the crude oil in the pipeline, expressed in per mille (‰), and is a known measured constant.

[0108] S3. Construct a fully connected deep neural network (PINN). Network topology design: The input layer consists of 3 neurons, each corresponding to a different spatiotemporal coordinate. The hidden layers consist of eight fully connected layers, each containing 100 neurons. The activation function chosen is the hyperbolic tangent function, which has C∞ smoothness and ensures the existence and non-zero nature of the second derivative, a necessary condition for solving the diffusion equation.

[0109] The output layer consists of three neurons, each outputting a concentration. ,isotope , source item The source term nonnegativity constraint on... Apply the Softplus activation function to the output, forcing ≥0, which is consistent with physical reality (there is no "negative leakage").

[0110] S4. Construct a hybrid loss function that includes physical residuals and sparsity regularization. (1) Automatic differentiation Using TensorFlow's tf.GradientTape to record the forward propagation computation graph, we can calculate the partial derivatives of the network output with respect to the input, including the first derivative. , , , , , Second derivative: , , , .

[0111] Indicates pollutant concentration The first derivative with respect to time t, Indicates pollutant concentration The first derivative with respect to the x-direction, Indicates pollutant concentration The first derivative with respect to the y-direction; Indicate intermediate variables The first derivative with respect to time t, Indicate intermediate variables The first derivative with respect to the direction x, Indicate intermediate variables The first derivative with respect to the y-direction; Indicates pollutant concentration The second derivative with respect to the x-direction, Indicates pollutant concentration The second derivative with respect to the y-direction, Indicate intermediate variables The second derivative in the x-direction, This represents the second derivative of the intermediate variable with respect to the y-direction. Data prediction is then performed based on the concentration conservation equation and the isotope mass conservation equation.

[0112] (2) Calculate the residuals of the physical equations Define Rayleigh fractionation coefficient =0.999, this value was determined by linear regression fitting of the Rayleigh equation in an indoor microbial degradation batch experiment, corresponding to the crude oil biodegradation scenario, and falls within the typical crude oil fractionation coefficient range of 0.998 to 0.999; degradation rate =0.005 Initial isotope ratios The value is -30‰, and then the concentration conservation equation is used to calculate... The isotope mass conservation equation was used to calculate .

[0113] (3) Weighted summation of the total loss function based on Calculate the loss function .

[0114] S5. Parameter Iteration and Source Term Reconstruction Based on Two-Stage Optimization S51. Fast Convergence: Learning rate set to 1×10 3 Training for 30,000 epochs. When The data drops rapidly; the network learns to fit the observation data, but the physical residuals... It is still quite large.

[0115] S52, Physical Consistency Convergence: Switching to the L-BFGS-B optimizer (quasi-Newton method), after approximately 5000 iterations, Reduced to 1.5×10 5 This indicates that the network prediction results are highly consistent with the laws of physical conservation.

[0116] S53. Source Term Reconstruction: Extract the global source term output of the trained grid at time t=0. ,based on Draw a heat map. The values ​​in the background (non-leakage area) of the heat map are all less than 10. 4 (Suppressed by L1 regularization). A bright spot appears near coordinates (150, 0.2), with a peak intensity of 5.8 mol / L / d; this can be used to determine the leak point.

[0117] In this embodiment, the actual leak point is set at x=150m, y=0m. The centroid coordinates of the bright spot determined by the inversion method of this invention are x=151.2m, y=0.15m. The longitudinal positioning error is only 0.15m, and the axial error is 1.2m (relative error <0.3%). If the isotope coupling residual in step three is removed... Using only concentration data for inversion, the resulting heat map shows a long, trailing shape in the x∈[120,200] range, making it impossible to accurately locate the specific leak point.

[0118] Therefore, this invention introduces an isotope fractionation field coupling mechanism, which solves the ill-posedness of the inversion problem by increasing the physical constraint dimension and significantly improves the positioning accuracy.

[0119] Finally, it should be noted that the above content is only used to illustrate the technical solution of the present invention, and is not intended to limit the scope of protection of the present invention. Simple modifications or equivalent substitutions made by those skilled in the art to the technical solution of the present invention do not depart from the essence and scope of the technical solution of the present invention.

Claims

1. A pipeline pollution source tracing method based on PINN and isotopic flow field coupling, characterized in that, include: S1. Construct a multidimensional spatiotemporal training dataset; Collect data from sparse monitoring wells around the pipeline, including pollutant concentrations. isotope ratio with monomer It also records the coordinates of the monitoring wells and the sampling time to form a spatiotemporal coordinate system. As a training set ; S2. Establish a set of governing equations for anisotropic reactive solute transport; It includes the concentration conservation equation and the isotope mass conservation equation, and introduces the isotope mass density as an intermediate variable into the isotope mass conservation equation. ; S3. Construct a fully connected deep neural network; Building based on spatiotemporal coordinates As input, with The output deep neural network, where, This indicates the output pollutant concentration. This represents an intermediate variable representing the output isotopic mass density. The leakage source phase field to be inverted; S4. Construct a hybrid loss function that includes physical residuals and sparse regularization; Constructing the loss function loss function Due to observation data error Physical equation residuals Source-phase sparsity constraint Composition, including the physical equation residuals It was calculated using the concentration conservation equation and the isotope mass conservation equation; S5. Parameter iteration and source term reconstruction based on two-stage optimization; By minimizing Update network weights and extract the output layer after training. Distribution, based on Draw a heat map; determine the pipeline leak point based on the heat map; Source phase sparsity constraint for: This indicates the number of grid points used to evaluate the distribution of the source terms; This represents the predicted source term strength output by the network at grid point k; This represents the L1 norm.

2. The pipeline pollution source tracing method according to claim 1, characterized in that, Step S1, when constructing the training set, also includes normalizing the monitoring well data.

3. The pipeline pollution source tracing method according to claim 1, characterized in that, The concentration conservation equation is: The isotopic mass conservation equation is: Where C represents the pollutant in spatial coordinates and molar concentration at time t; Indicates an intermediate variable in isotopic mass density; ; Indicates the isotope ratio of a single monomer; , This represents the components of the velocity vector v in the x and y directions; , These represent the hydrodynamic dispersion coefficients along the x-axis and y-axis, respectively; The rate constant representing the first-order biodegradation rate of a pollutant; Represents the Rayleigh fractionation coefficient; S represents the intensity of the leakage source phase to be inverted; This indicates the initial isotope ratio of the crude oil in the pipeline.

4. The pipeline pollution source tracing method according to claim 3, characterized in that, The components of the velocity vector are calculated based on Darcy's law as follows: in, Indicates water head. Indicates effective porosity. This represents the permeability coefficient at various points in space, expressed in m / d.

5. The pipeline pollution source tracing method according to claim 4, characterized in that, Permeability coefficient at various points in space The calculation is based on the following formula: Define the permeability field Assume the centerline of the pipeline is located at The half-width of the trench area is Then the permeability coefficient distribution function is: in: This represents the background permeability coefficient of the undisturbed soil surrounding the trench; This indicates the permeability coefficient of the backfill soil inside the trench; This indicates the coordinate position of the oil pipeline centerline on the y-axis; Indicates the smoothing factor; This represents the hyperbolic tangent function.

6. The pipeline pollution source tracing method according to claim 1, characterized in that, Deep neural networks constrain the output layer, specifically including: concentration output. and source term output Applying the Softplus activation function to isotope intermediate variables No non-negativity constraints are applied.

7. The pipeline pollution source tracing method according to claim 1, characterized in that, loss function for: in, , and These are the weights for each sub-item.

8. The pipeline pollution source tracing method according to claim 7, characterized in that, Physical equation residuals for: in, This represents the total number of configuration points randomly sampled within the computational domain; This represents the algebraic residual calculated after substituting the network predictions into the concentration conservation equation; This represents the algebraic residual value calculated after substituting the network prediction values ​​into the isotopic mass conservation equation.

9. The pipeline pollution source tracing method according to claim 1, characterized in that, Step S5 specifically includes: S51, Global Optimization: Using the Adam optimizer, setting the initial learning rate and momentum parameters, and iteratively training. Second, reduce data loss ; S52. Physical Consistency Convergence: Switch to the L-BFGS-B optimizer and perform second-order optimization until the physical equation residuals are found. < , The preset convergence threshold is used; S53. Source Term Reconstruction: Extract the global source term output of the trained grid at time t=0. ,based on Draw a heat map; set a peak intensity threshold, select bright spots in the heat map whose peak intensity is greater than the peak intensity threshold, obtain the centroid coordinates of the bright spots, and then determine the leak point.

Images