A gas hydrate phase equilibrium prediction method based on physical thermodynamic model and machine learning fusion

Abstract

This invention discloses a method for predicting the phase equilibrium of gas hydrates based on the fusion of physical thermodynamic models and machine learning, comprising: S1: collecting mole fractions and temperature data of components in a gas mixture to generate an experimental dataset; S2: based on the experimental dataset, using a thermodynamic model to calculate the equilibrium pressure and determine the crystal structure type, generating a physical calculation result set; S3: based on the physical calculation result set, calculating the logarithmic difference between the experimental pressure and the predicted physical pressure, generating a logarithmic residual vector; S4: based on the physical calculation result set, constructing a feature matrix containing nonlinear and interactive features; S5: based on the logarithmic residual vector and the feature matrix, performing data partitioning and standardization to generate a training dataset; S6: based on the training dataset, training a gradient boosting model and an ensemble learning model and performing a weighted combination to generate a residual prediction model; S7: based on the residual prediction model and the predicted physical pressure, calculating the corrected equilibrium pressure.

Description

Technical Field

[0001] This invention relates to the field of phase equilibrium thermodynamics, and in particular to a method for predicting the phase equilibrium of gas hydrates based on the fusion of physical thermodynamic models and machine learning. Background Technology

[0002] Natural gas hydrates are crystalline inclusions formed by gas and water molecules under high pressure and low temperature conditions. They are characterized by large reserves, wide distribution, high energy density, and clean combustion, and are considered one of the most promising alternative energy sources of the 21st century. Natural gas hydrates are widely found in seabed sediments and permafrost zones, and their formation and stability are influenced by various thermodynamic factors such as temperature, pressure, and gas composition. Hydrates can exist in various phases, including gaseous, liquid, and solid states, under different environments. The complex three-phase dynamic transformation process plays a decisive role in the formation depth, occurrence state, and extraction stability of hydrates. Accurately predicting the phase equilibrium conditions of gas hydrates is a core issue in hydrate resource assessment and extraction design.

[0003] Existing techniques for predicting the thermodynamic behavior of gas hydrates mainly fall into two categories, both of which have significant limitations. One category is based on physical mechanism-based thermodynamic models, such as the van der Waals–Platteeuw (vdW-P) model, which combines the Kihara potential function with the Peng-Robinson equation of state to estimate equilibrium pressure at a given temperature. While these methods have clear physical meaning, they rely on idealized assumptions such as cavity symmetry, single-molecule occupancy, and no molecular interactions, leading to significant prediction errors when dealing with multi-component gas mixtures, complex boundary conditions, or highly nonlinear systems. Furthermore, the computational complexity and difficulty in parameter tuning hinder their widespread adoption in large-scale data analysis and rapid engineering applications. The other category is based on data-driven machine learning methods, which construct nonlinear regression models by fitting experimental data, improving prediction efficiency and accuracy to some extent. However, purely data-driven models lack physical constraints and are considered "black box" algorithms. They are prone to overfitting and prediction distortion under boundary extrapolation or sparse data conditions, exhibiting insufficient generalization ability and failing to meet high-reliability engineering requirements.

[0004] Currently, no single method can fully leverage data-driven capabilities to achieve high-precision, highly generalizable, and physically interpretable predictions of hydrate phase equilibrium while maintaining physical consistency. Especially in engineering applications involving multi-component gases, extreme conditions, and the identification of three-phase stability boundaries, existing methods exhibit significant shortcomings in accuracy, stability, and scalability. Therefore, there is an urgent need to develop a novel prediction scheme that integrates the physical mechanisms of thermodynamic models with the advantages of machine learning and data mining. Summary of the Invention

[0005] This invention provides a gas hydrate phase equilibrium prediction method based on the fusion of physical thermodynamic models and machine learning, in order to solve the above-mentioned problems existing in the prior art.

[0006] To achieve the above objectives, the present invention provides the following technical solution: A method for predicting gas hydrate phase equilibrium based on the fusion of physical thermodynamic models and machine learning includes: S1: Collect the mole fraction and temperature data of the gas mixture to generate an experimental dataset; S2: Based on the experimental dataset, the equilibrium pressure is calculated using a thermodynamic model and the crystal structure type is determined, generating a set of physical calculation results. S3: Based on the physical calculation result set, calculate the logarithmic difference between the experimental pressure and the physical predicted pressure, and generate a logarithmic residual vector; S4: Based on the physical calculation result set, construct a feature matrix containing nonlinear and interactive features; S5: Based on the logarithmic residual vector and feature matrix, perform data partitioning and standardization to generate a training dataset; S6: Based on the training dataset, train the gradient boosting model and the ensemble learning model, and combine them with weights to generate a residual prediction model; S7: Calculate the corrected equilibrium pressure based on the residual prediction model and physical prediction pressure.

[0007] Furthermore, S2 includes: S21: Traverse the experimental dataset and load potential energy parameters, state equation parameters and binary interaction parameters according to the component information; S22: Determine the crystal structure type of hydrates based on the structural formation tendency of the components; S23: Based on the selected crystal structure, the equilibrium pressure is calculated by coupling the van der Waals-Platteeuw model with the Peng-Robinson equation of state and the Kihara potential model; S24: Filter the data points that converge to generate a physical calculation result set containing experimental data, physical predicted pressure, and structure type.

[0008] Furthermore, in S22, the determination of the hydrate crystal structure type includes: When the system is a single component, the preset preferred structure of that component is selected; When the components of a system tend to form the same structure, that common structure is chosen; When the system contains different structure-forming agents, the equilibrium pressure of each structure is calculated, and the structure with the smaller deviation between the pressure value and the experimental value is selected.

[0009] Furthermore, S4 includes: S41: Extract temperature, component mole fraction, logarithmic value of physical predicted pressure, and structure type as basic features; S42: Construct nonlinear temperature characteristics, including: power terms, reciprocal terms, and radical terms of temperature; S43: Calculate the sum of the mole fractions of a specific component and the ratio of the mole fractions between components as characteristics of component polymerization; S44: Construct interactive features, including: the product of component mole fraction and temperature; S45: Construct physical predictions as nonlinear characteristics of pressure, including: power terms of pressure; S46: Combine various features to generate a feature matrix.

[0010] Furthermore, S5 includes: S51: Use a clustering algorithm to group the logarithmic residual vector; S52: Based on the clustering results, perform stratified sampling on the feature matrix and log residual vector to divide the training set and the test set; S53: Calculate statistical parameters for the training set and standardize the training and test sets.

[0011] Furthermore, the gradient boosting model in S6 is the LightGBM model, and the ensemble learning model is the Random Forest model.

[0012] Furthermore, S6 includes: S61: Optimize the hyperparameters of the LightGBM model through grid search cross-validation, and train the model by assigning sample weights based on residual size and stress value; S62: Train a random forest model using the same sample weights; S63: Perform a linear weighted combination of the prediction results from the two models to generate a residual prediction model.

[0013] Furthermore, S7 includes: S71: Input the test data into the residual prediction model to obtain the predicted logarithmic residual; S72: Calculate the exponential function value of the sum of the logarithmic value of the physical predicted pressure and the logarithmic residual of the prediction to obtain the corrected equilibrium pressure.

[0014] Furthermore, the formula for calculating the logarithmic residual is: in, This represents the corresponding actual observed target value. This represents the equilibrium pressure predicted by the physical model; The revised formula for calculating equilibrium pressure is as follows: in, This indicates the corrected equilibrium pressure. This represents the ensemble model.

[0015] Furthermore, it also includes: S8: Persistently store the trained gradient boosting model, ensemble learning model, standardized parameters, and weighting coefficients.

[0016] Compared with the prior art, the present invention has the following advantages: This invention organically combines the theoretical rigor of physical models with the data-driven capabilities of machine learning models, correcting systematic biases in complex, non-ideal systems through residual learning. Compared to traditional physical models, this method significantly improves prediction accuracy while maintaining good extrapolation performance. Compared to purely data-driven models, this method, with its physical mechanism as a "skeleton," is less prone to producing absurd predictions that violate the laws of thermodynamics, exhibiting stronger robustness in conditions where training data is sparse or exceeds the training range. This method can be reliably applied to a wide range of operating conditions, from pure components to multi-component systems, and from low-temperature, low-pressure systems to high-temperature, high-pressure systems, demonstrating strong practicality.

[0017] Other features and advantages of the invention will be set forth in the description which follows, and will be apparent in part from the description, or may be learned by practicing the invention.

[0018] The technical solution of the present invention will be further described in detail below with reference to the accompanying drawings and embodiments. Attached Figure Description

[0019] The accompanying drawings are provided to further illustrate the invention and form part of the specification. They are used in conjunction with embodiments of the invention to explain the invention and do not constitute a limitation thereof. In the drawings: Figure 1 This is a flowchart of a gas hydrate phase equilibrium prediction method based on the fusion of physical thermodynamic models and machine learning in an embodiment of the present invention; Figure 2 This is a schematic diagram of the gas hydrate phase equilibrium prediction method based on the fusion of physical thermodynamic model and machine learning in an embodiment of the present invention. Detailed Implementation

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

[0021] The embodiments of the present invention provide, as follows Figure 1 As shown, a method for predicting the phase equilibrium of gas hydrates based on the fusion of physical thermodynamic models and machine learning includes: S1: Collect the mole fraction and temperature data of the gas mixture to generate an experimental dataset; S2: Based on the experimental dataset, the equilibrium pressure is calculated using a thermodynamic model and the crystal structure type is determined, generating a set of physical calculation results. S3: Based on the physical calculation result set, calculate the logarithmic difference between the experimental pressure and the physical predicted pressure, and generate a logarithmic residual vector; S4: Based on the physical calculation result set, construct a feature matrix containing nonlinear and interactive features; S5: Based on the logarithmic residual vector and feature matrix, perform data partitioning and standardization to generate a training dataset; S6: Based on the training dataset, train the gradient boosting model and the ensemble learning model, and combine them with weights to generate a residual prediction model; S7: Calculate the corrected equilibrium pressure based on the residual prediction model and physical prediction pressure.

[0022] The following is a detailed description with reference to specific embodiments.

[0023] like Figure 2As shown, input the original experimental data, including temperature (T), mole fraction of gas components (X), and experimental equilibrium pressure (Pexp). Based on the input data, a preliminary prediction is made using a physical mechanism model (vdW-P model). This stage includes an intelligent crystal structure selection module that automatically determines the optimal hydrate crystal structure (sI or sII) based on thermodynamic stability principles and outputs the physical model's predicted pressure (Pphys). Residual calculation: The predicted pressure (Pphys) of the physical model is compared with the experimentally measured pressure (Pexp), and the logarithmic residual (ln(Pexp) - ln(Pphys)) is calculated. This residual serves as the target for the machine learning model to learn and correct. Feature engineering: Before inputting the data into the machine learning model, advanced feature engineering is performed to construct a high-dimensional feature matrix containing basic features (such as temperature and composition), nonlinear features (such as T², 1 / T), and interaction features to fully exploit the data's potential. The feature data is standardized, and K-Means clustering is performed based on the residuals. Then, a stratified sampling method is used to divide the dataset into training and test sets. Two machine learning models are trained in parallel using the training set: the LightGBM model and Random Forest. The model employs grid search cross-validation for hyperparameter optimization during training. Finally, the predictions from the two models are weighted and ensembled to obtain an ensemble residual prediction model. Test set data is input into the trained ensemble model to obtain residual prediction values. These residual predictions are then combined with the predictions from the physical model to generate the final prediction pressure (Pfinal). The model performance is comprehensively evaluated using metrics such as RMSE, MAE, and R² by comparing Pfinal with the true value Pexp. The trained machine learning model, normalizer, and other components are persistently saved (e.g., as a .joblib file), forming a deployment model that can be directly used for rapid prediction of new data, and outputs the final prediction results.

[0024] This scheme is based on predicting a methane-containing ( ), ethane ( ), propane ( ) and carbon dioxide ( Taking the equilibrium pressure of the hydrate phase of a four-component gas mixture at a specific temperature as an example, the implementation process of this method is explained in detail.

[0025] S1: Collect the mole fraction and temperature data of the gas mixture to generate an experimental dataset; Software environment setup: Configure a Python 3.8 or higher environment on a computer and install the necessary scientific computing libraries, including Pandas for data processing, NumPy for numerical computation, SciPy for scientific computing and optimization, and machine learning libraries Scikit-learn, LightGBM, and Joblib.

[0026] Physical parameter library preparation: Create an Excel file named hydrate_parameters.xlsx. This file contains three worksheets: The Components worksheet stores the components (e.g., ...). , Critical temperature (Tc), critical pressure (Pc), acentric factor, and preferred hydrate structure type (e.g., For sI, (for sII).

[0027] Kihara worksheets: store the Kihara spherical core potential model parameters for each component, including the core radius ( ), potential well depth parameter ( ) and potential well depth ( ).

[0028] BIPs worksheet: Stores binary interaction parameters (kᵢⱼ) between different component pairs.

[0029] Raw experimental data collection: Collect and organize experimental data on gas hydrate phase equilibrium containing the above four components from publicly available academic literature or authoritative databases. Data collection covers multiple gases (such as...). , , , Experimental data on the phase equilibrium of gas hydrates in the system (etc.). The data were compiled into a CSV file named date.csv. Each line in this file represents an independent experimental point, containing the mole fraction of the components and temperature data, including at least the following: T(K) (experimental temperature Texp). (molar fraction of each component xexp), and p (MPa) (equilibrium pressure Pexp measured experimentally).

[0030] Generate the experimental dataset: Organize the collected data into a structured experimental dataset. Each record in this dataset contains three core pieces of information: an experimental temperature Texp (unit: K), a gas component mole fraction vector xexp, and a corresponding experimental equilibrium pressure Pexp (unit: MPa).

[0031] S2: Based on the experimental dataset, a thermodynamic model is used to calculate the equilibrium pressure and determine the crystal structure type, generating a set of physical calculation results. In this step, a Python script is written to process each data point in date.csv using thermodynamic model functions.

[0032] S21: Traverse the experimental dataset and load potential energy parameters, state equation parameters, and binary interaction parameters based on the component information. Load the experimental dataset: Load the date.csv file and read each experimental condition point (Texp, xexp, Pexp) line by line.

[0033] Physical parameter loading: For each operating point, based on the component information at the operating point, the physical property parameters of each component are loaded from the parameter library (including Kihara potential energy parameters, Peng-Robinson equation of state parameters, and binary interaction parameters BIPs).

[0034] S22: Determine the crystal structure type of hydrates based on the structural formation tendency of the components. Intelligent crystal structure selection: Based on prior physical knowledge and the principle of thermodynamic stability, the optimal hydrate crystal structure (sI or sII) is automatically determined. The decision logic for determining the hydrate crystal structure type is as follows: a) When the system is a single component, select the preset preferred structure of that component: If the system is a single component, select the preset preferred generation structure of that component in the parameter library.

[0035] b) When the components of the system tend to form the same structure, select the common structure: if all components in the system tend to form the same structure, then directly select that structure.

[0036] c) When the system contains different structure-forming agents, calculate the equilibrium pressure of different structures and select the structure with the smaller deviation from the experimental value: If the system contains different types of structure-forming agents, calculate the equilibrium pressure of the two structures sI and sII respectively, and select the structure that makes the equilibrium pressure closer to the true value as the optimal structure, because it is more thermodynamically stable.

[0037] Example: Assume the current operating condition is T = 278.15 K, and the composition is { 0.85, 0.10, 0.03, : 0.02}. This function first detects the presence of an existing Si-forming agent in the mixture ( , , ), and also sII forming agents ( This triggers the thermodynamic competition rule.

[0038] S23: Based on the selected crystal structure, the equilibrium pressure is calculated by coupling the van der Waals-Platteeuw model with the Peng-Robinson equation of state and the Kihara potential model; Equilibrium pressure calculation: Based on the selected crystal structure, the core calculation is performed using the van der Waals–Platteeuw (vdW-P) thermodynamic model. This process internally couples the Peng-Robinson (PR) equation of state for calculating gas fugacity and the Kihara spherical potential model for calculating the Langmuir constant. Finally, the predicted equilibrium pressure Pphys of the physical model is obtained by solving the problem using numerical optimization methods (such as minimize_scalar).

[0039] S24: Filter the data points that converge to the calculation and generate a set of physical calculation results containing experimental data, physical predicted pressure, and structure type; Data filtering: Remove operating points that the physical model cannot successfully calculate or converge from the dataset to ensure the quality of subsequent machine learning training data.

[0040] Generate a set of physical calculation results: Store the calculated Pphys and Ssel (selected crystal structure types, 1 for sI, 2 for sII) along with the original operating condition data. If a calculation fails at a certain operating point (e.g., numerical optimization does not converge), it is marked and removed from subsequent processing. Finally, a filtered and enhanced set of physical calculation results is generated. Each record in this set contains: the original (Texp, xexp, Pexp), and two newly added calculated values: a physical predicted pressure Pphys and a selected crystal structure type Ssel.

[0041] S3: Based on the physical calculation result set, calculate the logarithmic difference between the experimental pressure and the physical predicted pressure, and generate a logarithmic residual vector; Constructing the logarithmic residual: Using the physical calculation result set generated in the previous step, calculate the logarithmic difference between the experimental pressure and the physical predicted pressure. The formula for calculating the logarithmic residual is: ΔY = ln(Y) - ln(Yphys) Where Y represents the corresponding actual observed target value (i.e., experimental equilibrium pressure Pexp), and Yphys represents the equilibrium pressure Pphys predicted by the physical model.

[0042] Generate a logarithmic residual vector: This logarithmic residual vector, ytarget, will serve as the target for subsequent machine learning model training. The machine learning part aims to capture the non-linear relationship between this logarithmic difference and the input features.

[0043] S4: Based on the physical calculation results set, construct a feature matrix containing nonlinear and interactive features. To fully explore the potential correlation between input variables and logarithmic residuals, a high-dimensional feature matrix is ​​constructed. This process not only includes basic variables but also creatively introduces physically meaningful nonlinear and interactive features.

[0044] S41: Extraction temperature, component mole fraction, logarithm of physical predicted pressure, and structure type as basic features. Basic characteristics: experimental temperature Texp, mole fraction of each gas component xexp, logarithmic value of pressure predicted by the physical model log(Pphys), and selected crystal structure type Ssel.

[0045] S42: Construct power terms, reciprocal terms, and radical terms of temperature as nonlinear characteristics of temperature; Temperature nonlinearity characteristics: Introducing the square term (T²), reciprocal term (1 / T), and square root term of temperature (T²). This is used to capture the nonlinear changes in heat and mass transfer effects within different temperature ranges.

[0046] S43: Calculate the sum of the mole fractions of a specific component and the ratio of the mole fractions between components as characteristics of component polymerization; Component polymerization characteristics: Calculation of heavy components (e.g.) The total mole fraction of methane and the ratio of methane to other gaseous components are used to characterize the overall properties of the gas mixture.

[0047] S44: Construct the product of component mole fraction and temperature as an interaction feature; Interactive features: Constructing key components (such as...) The product of mole fraction and temperature is used to describe the weighted changes in the influence of a specific component at different temperatures.

[0048] S45: Construct power terms for physical prediction of pressure as nonlinear characteristics of pressure; Nonlinear characteristics of pressure: Introducing the square root term of physical pressure prediction ( The sum and square term (Pphys²) are used to capture the complex effects of pressure on system non-idealities.

[0049] S46: Combine various features to generate a feature matrix; Generating the feature matrix: For each data point, the `create_enhanced_features` function is used to combine various features to construct a multi-dimensional feature vector. The feature vectors of all data points are then stacked to form a high-dimensional, information-rich feature matrix, `Xfeatures`.

[0050] S5: Based on the logarithmic residual vector and feature matrix, perform data partitioning and standardization to generate a training dataset; S51: Use a clustering algorithm to group the logarithmic residual vector; Residual Clustering: The K-Means clustering algorithm is used to perform unsupervised clustering of the logarithmic residual vector ytarget. sklearn.cluster.KMeans is used to cluster the logarithmic residual vector, with the number of clusters n_clusters=4. This will divide the samples into four groups with different error characteristics. This step aims to group sample points with similar prediction errors of the physical model (e.g., overestimated, underestimated, and accurate regions) into the same cluster.

[0051] S52: Based on the clustering results, perform stratified sampling on the feature matrix and log residual vector to divide the training set and the test set. Stratified Sampling: Based on the clustering results described above, stratified sampling is performed on the feature matrix and log residual vector to divide the dataset into training and test sets. The `sklearn.model_selection.train_test_split` function is used to split the feature matrix `Xfeatures` and the log residual vector `ytarget` into 80% training and 20% test sets. Crucially, the clustering labels of `KMeans` are passed as the `stratify` parameter to ensure that the proportions of the four error classes are consistent between the training and test sets. This method ensures that the error distributions of the training and test sets are consistent, avoiding model bias towards certain error patterns during training and thus improving the model's generalization ability.

[0052] S53: Calculate statistical parameters for the training set and standardize the training and test sets; Feature standardization: The feature matrix Xfeatures is standardized using StandardScaler. An sklearn.preprocessing.StandardScaler object is instantiated. The mean and standard deviation (fit) are calculated on the training set. The .fit_transform() method is called on the training set to obtain the standardized training features Xtrain_scaled. This transformation is then applied to the test set (transform). The .transform() method is called on the test set, using the mean and variance learned from the training set to obtain the standardized test features Xtest_scaled. This eliminates the dimensional differences between different features, accelerating the convergence speed of the machine learning model and improving performance.

[0053] Generate training datasets: split and standardized training datasets (Xtrain_scaled, ytrain) and test datasets (Xtest_scaled, ytest).

[0054] S6: Based on the training dataset, train the gradient boosting model and the ensemble learning model, and combine them with weights to generate a residual prediction model; This application employs two heterogeneous machine learning models for parallel training and combines their advantages through weighted ensemble learning to achieve more robust residual prediction capabilities. The gradient boosting model is the LightGBM model, and the ensemble learning model is the Random Forest model.

[0055] S61: Optimize the hyperparameters of the LightGBM model through grid search cross-validation, and train the model by assigning sample weights based on residual size and stress value; LightGBM model training: a) The LightGBM (Lightweight Gradient Boosting Machine) model is used, which is an efficient gradient boosting decision tree algorithm. A lightgbm.LGBMRegressor model is instantiated.

[0056] b) Define a hyperparameter grid, for example: {'n_estimators': [200, 400], 'learning_rate': [0.05, 0.1], 'num_leaves': [31, 50]}. Utilize GridSearchCV (Grid Search Cross-Validation) to automatically find the optimal hyperparameter combination within a predefined hyperparameter space (e.g., learning rate, number of trees, number of leaf nodes) to maximize model performance on the validation set. Use `sklearn.model_selection.GridSearchCV`, passing in the model and parameter grid, and perform 3-fold cross-validation optimization using `neg_mean_squared_error` as the scoring criterion.

[0057] c) During training, an intelligent sample weighting strategy is introduced, assigning higher weights to samples with large logarithmic residuals or those in low-pressure regions, guiding the model to focus on these regions that are difficult to predict or have significant engineering implications. When calling the .fit() method, pass in Xtrain_scaled, ytrain, and a sample_weight array calculated based on the residual size and pressure value of the training samples.

[0058] d) After training, gscv.best_estimator_ is the optimal LightGBM model lgb_model.

[0059] S62: Train a random forest model using the same sample weights; Random Forest model training: a) In parallel, a random forest model is used, which is an ensemble learning algorithm based on the Bagging idea, known for its high accuracy and resistance to overfitting. Instantiate an sklearn.ensemble.RandomForestRegressor model and set appropriate parameters (e.g., n_estimators=300, max_features='sqrt').

[0060] b) Train using the same sample weights defined in S61 to complement LightGBM. Train directly on Xtrain_scaled, ytrain, and the same sample_weights using the .fit() method to obtain rf_model.

[0061] S63: Perform a linear weighted combination of the prediction results from the two models to generate a residual prediction model; Weighted combination: The trained LightGBM model and the Random Forest model are linearly weighted and combined. The log residual predictions of the two models for the same sample are summed according to preset weights (e.g., LightGBM weight is 0.65, and Random Forest weight is 0.35) to obtain the final ensemble log residual prediction value Δypred.

[0062] Generate a residual prediction model: an optimized and integrated residual prediction model fusemble, in which the physical model undertakes the main trend fitting, while machine learning focuses on correcting the multiplicative nonlinear biases that the physical model struggles to capture through formulas.

[0063] S7: Calculate the corrected equilibrium pressure based on the residual prediction model and the physical predicted pressure; Now, for a completely new and unseen operating condition, prediction is needed: temperature Tnew = 280.0 K, composed of { 0.90, 0.05, 0.05, : 0.0}.

[0064] Calling the physical model: First, call the thermodynamic model function to obtain the physical predicted pressure Pphys_new and structural parameters Ssel_new for the new operating condition. Assume the results are Pphys_new = 5.0 MPa and Ssel_new = sI.

[0065] Constructing a feature vector: Use the methods in S4 to construct a multi-dimensional feature vector for this new operating point.

[0066] Feature standardization: Using the StandardScaler object already trained in S5, the feature vectors of the new working condition are transformed using .transform().

[0067] S71: Input the test data into the residual prediction model to obtain the predicted logarithmic residual; Residual prediction: The standardized feature vectors are input into the trained lgb_model and rf_model respectively to obtain two log residual prediction values: lgb_pred_residual and rf_pred_residual.

[0068] Based on the preset weighting weights in S63 (e.g., 0.65 and 0.35), calculate the ensemble predicted log residual: predicted_log_residual = 0.65 lgb_pred_residual + 0.35 rf_pred_residual.

[0069] S72: Calculate the exponential function value of the sum of the logarithmic value of the physical predicted pressure and the logarithmic residual of the prediction to obtain the corrected equilibrium pressure; Let the input features be X (temperature, gas composition), and the equilibrium pressure predicted by the physical model be: The corresponding actual observed target value is Y. The residual is defined as: The machine learning part aims to capture the non-linear relationship between the residual and the input features, i.e., to train the following mapping: in, This is an ensemble model consisting of LightGBM and Random Forest, with ε representing Gaussian white noise. The final corrected prediction is: Where Yfinal represents the corrected equilibrium pressure, and fML(X) represents the prediction result of the ensemble model (i.e., the residual prediction model fensemble).

[0070] Temperature Tnew = 280.0 K, composition is { 0.90, 0.05, 0.05, : 0.0}.

[0071] Invoking the physical model: First, call the `predict_structure_and_pressure` function to obtain the physical predicted pressure P_phys_new and structure S_sel_new for this new operating condition. Assume the results are Pphys_new = 5.0 MPa, Ssel_new = sI.

[0072] Constructing a feature vector: Use create_enhanced_features to construct a multi-dimensional feature vector for this new operating condition.

[0073] Feature standardization: Using the StandardScaler object trained in step five, transform the feature vectors of the new working condition using the .transform() function.

[0074] Residual prediction: The standardized feature vectors are input into the trained lgb_model and rf_model respectively to obtain two log residual prediction values: lgb_pred_residual and rf_pred_residual.

[0075] Calculate the logarithmic residual of the ensemble based on the preset ensemble weights (e.g., 0.65 and 0.35): predicted_log_residual = 0.65 lgb_pred_residual + 0.35 rf_pred_residual.

[0076] Generate the final predicted stress: The final corrected pressure is calculated using the formula P_final = exp(log(P_phys) + Δy_pred): P_final_new = np.exp(np.log(P_phys_new) + predicted_log_residual).

[0077] Assuming predicted_log_residual is -0.1, then P_final_new = np.exp(np.log(5.0) -0.1) ≈ 4.52 MPa.

[0078] Using various statistical indicators such as root mean square error (RMSE), mean absolute error (MAE), and coefficient of determination (R²), the prediction performance of the original physical model (Pphys) and the hybrid model of this invention (Pfinal) on the test set is comprehensively evaluated and compared, and the superiority of this method is quantitatively verified.

[0079] S8: Persistently store the trained gradient boosting model, ensemble learning model, standardized parameters, and weighting coefficients; Model persistence: Use the joblib.dump function to save the trained lgb_model, rf_model, StandardScaler objects and the configuration dictionary containing weighted weights and other information as .joblib files respectively.

[0080] Deploying the application: Package these model files with prediction scripts (such as the `predict_pressure_for_new_conditions` function). This allows new prediction tasks to be executed quickly and accurately by loading these files in any Python-enabled environment, without repeating the entire training process. This enables fast prediction by directly loading the model when faced with a new prediction task, without having to repeat the above training process.

[0081] The technical effects of this application are as follows: This method organically combines the theoretical rigor of physical models with the data-driven capabilities of machine learning models, using residual learning to correct systematic biases in complex, non-ideal systems. Compared to traditional physical models, this method significantly improves prediction accuracy (e.g., higher R² value, lower RMSE) while maintaining good extrapolation performance. Compared to purely data-driven models, this method, with its physical mechanism as a "skeleton," is less prone to absurd predictions that violate the laws of thermodynamics, exhibiting stronger robustness in conditions where training data is sparse or exceeds the training range. This method can be reliably applied to a wide range of operating conditions, from pure components to multi-component systems, and from low-temperature, low-pressure systems to high-temperature, high-pressure systems, demonstrating strong practicality.

[0082] Traditional pure physics model solving processes, especially when iterative optimization is needed to find the equilibrium point, are computationally intensive and time-consuming. This invention completes the prediction through a one-time physics model calculation and a one-time, efficient machine learning forward propagation. Because the trained machine learning model (such as LightGBM and Random Forest) performs predictions extremely quickly, the overall prediction time is far lower than that of pure physics methods involving repeated iterations. This high efficiency enables this method to meet the high computational speed requirements of applications such as online monitoring, dynamic simulation, and real-time optimization.

[0083] This invention provides more reliable data support for safe production and process optimization, offering significant economic and safety value. In industrial fields such as natural gas extraction, transportation, and storage, accurately predicting the hydrate formation boundary (i.e., the phase equilibrium pressure-temperature curve) is crucial. The high-precision prediction results provided by this invention can assist engineers in: In the design and operation of oil and gas pipelines, it is possible to more accurately determine the critical conditions for the addition of inhibitors or the safe operating temperature and pressure window to maintain, thereby effectively avoiding pipeline blockage caused by hydrate formation, reducing production interruptions and huge economic losses.

[0084] In deepwater oil and gas drilling and other scenarios, accurately predicting the stable region of hydrates in the formation helps to formulate safe drilling fluid density and temperature plans, and prevent major safety accidents such as wellbore instability and blowouts caused by accidental decomposition of hydrates.

[0085] In gas separation or storage processes based on hydrate technology, accurate phase equilibrium data is fundamental for optimizing operating parameters (such as pressure and temperature), improving separation efficiency, and reducing energy consumption. This invention provides crucial and reliable basic data for the development and optimization of such processes.

[0086] Obviously, those skilled in the art can make various modifications and variations to this invention without departing from the spirit and scope of this invention.

Claims

1. A method for predicting the phase equilibrium of gas hydrates based on the fusion of physical thermodynamic models and machine learning, characterized in that, include: S1: Collect the mole fraction and temperature data of the gas mixture to generate an experimental dataset; S2: Based on the experimental dataset, the equilibrium pressure is calculated using a thermodynamic model and the crystal structure type is determined, generating a set of physical calculation results. S3: Based on the physical calculation result set, calculate the logarithmic difference between the experimental pressure and the physical predicted pressure, and generate a logarithmic residual vector; S4: Based on the physical calculation result set, construct a feature matrix containing nonlinear and interactive features; S5: Based on the logarithmic residual vector and feature matrix, perform data partitioning and standardization to generate a training dataset; S6: Based on the training dataset, train the gradient boosting model and the ensemble learning model, and combine them with weights to generate a residual prediction model; S7: Calculate the corrected equilibrium pressure based on the residual prediction model and physical prediction pressure.

2. The gas hydrate phase equilibrium prediction method based on the fusion of physical thermodynamic model and machine learning according to claim 1, characterized in that, S2 include: S21: Traverse the experimental dataset and load potential energy parameters, state equation parameters and binary interaction parameters according to the component information; S22: Determine the crystal structure type of hydrates based on the structural formation tendency of the components; S23: Based on the selected crystal structure, the equilibrium pressure is calculated by coupling the van der Waals-Platteeuw model with the Peng-Robinson equation of state and the Kihara potential model; S24: Filter the data points that converge to generate a physical calculation result set containing experimental data, physical predicted pressure, and structure type.

3. The gas hydrate phase equilibrium prediction method based on the fusion of physical thermodynamic models and machine learning according to claim 2, characterized in that, S22 defines the crystal structure type of hydrates as follows: When the system is a single component, the preset preferred structure of that component is selected; When the components of a system tend to form the same structure, that common structure is chosen; When the system contains different structure-forming agents, the equilibrium pressure of different structures is calculated separately, and the structure whose pressure value deviates from the experimental value is less than the preset value is selected.

4. The gas hydrate phase equilibrium prediction method based on the fusion of physical thermodynamic model and machine learning according to claim 1, characterized in that, S4 include: S41: Extract temperature, component mole fraction, logarithmic value of physical predicted pressure, and structure type as basic features; S42: Construct nonlinear temperature characteristics, including: power terms, reciprocal terms, and radical terms of temperature; S43: Calculate the sum of the mole fractions of a specific component and the ratio of the mole fractions between components as characteristics of component polymerization; S44: Construct interactive features, including: the product of component mole fraction and temperature; S45: Construct nonlinear characteristics of physical prediction pressure, including: power terms of pressure; S46: Combine various features to generate a feature matrix.

5. The gas hydrate phase equilibrium prediction method based on the fusion of physical thermodynamic model and machine learning according to claim 1, characterized in that, S5 include: S51: Use a clustering algorithm to group the logarithmic residual vector; S52: Based on the clustering results, perform stratified sampling on the feature matrix and log residual vector to divide the training set and the test set; S53: Calculate statistical parameters for the training set and standardize the training and test sets.

6. The gas hydrate phase equilibrium prediction method based on the fusion of physical thermodynamic models and machine learning according to claim 1, characterized in that, The gradient boosting model in S6 is the LightGBM model, and the ensemble learning model is the Random Forest model.

7. The gas hydrate phase equilibrium prediction method based on the fusion of physical thermodynamic model and machine learning according to claim 6, characterized in that, S6 include: S61: Optimize the hyperparameters of the LightGBM model through grid search cross-validation, and train the model by assigning sample weights based on residual size and stress value; S62: Train a random forest model using the same sample weights; S63: Perform a linear weighted combination of the prediction results from the two models to generate a residual prediction model.

8. The gas hydrate phase equilibrium prediction method based on the fusion of physical thermodynamic model and machine learning according to claim 1, characterized in that, S7 includes: S71: Input the test data into the residual prediction model to obtain the predicted logarithmic residual; S72: Calculate the exponential function value of the sum of the logarithmic value of the physical predicted pressure and the logarithmic residual of the prediction to obtain the corrected equilibrium pressure.

9. The gas hydrate phase equilibrium prediction method based on the fusion of physical thermodynamic models and machine learning according to claim 1, characterized in that, The formula for calculating the logarithmic residual is: ； in, This represents the corresponding actual observed target value. This represents the equilibrium pressure predicted by the physical model; The revised formula for calculating equilibrium pressure is as follows: ； in, This indicates the corrected equilibrium pressure. This represents the ensemble model.

10. The gas hydrate phase equilibrium prediction method based on the fusion of physical thermodynamic models and machine learning according to claim 1, characterized in that, Also includes: S8: Persistently store the trained gradient boosting model, ensemble learning model, standardized parameters, and weighting coefficients.

Images