Intelligent management and control system for drill pipe heat treatment based on digital twinning

Abstract

This invention relates to the field of simulation technology, specifically to an intelligent control system for heat treatment of drill pipes based on digital twins. The system collects environmental parameters, medium parameters, and workpiece state parameters from the heat treatment production line and outputs a time-series input vector. It reconstructs heat exchange boundary conditions using the medium parameters contained in the time-series input vector. Combining these boundary conditions, it drives a strongly coupled calculation based on a phase transformation-induced strain mechanism, involving a temperature field submodule, a microstructure evolution submodule, and a stress field submodule. This outputs a virtual state matrix containing the temperature distribution across the entire cross-section, grain size distribution, phase transformation volume fraction distribution, and stress distribution. Based on this virtual state matrix, a multi-objective optimization algorithm is executed in a multi-dimensional process parameter space to find the optimal balance point of strength and toughness indicators while satisfying microstructure and stress constraints. This generates an optimal heating power curve and dynamic quenching parameters. Finally, it optimizes production through a process instruction set.

Description

Technical Field

[0001] This invention relates to the field of simulation technology, specifically to an intelligent control system for heat treatment of drill pipes based on digital twins. Background Technology

[0002] With the development of industrial digitalization, using computer-aided engineering technology to numerically simulate the process has become an important means of optimizing process parameters and predicting product performance. Existing heat treatment simulation technology mainly relies on finite element analysis or finite difference method to solve coupled partial differential equations.

[0003] When performing high-fidelity full-process simulations of large components like drill pipes, static or idealized preset parameters, such as a constant heat transfer coefficient, are often used. However, in actual physical conditions, the rheological properties of the cooling medium, nozzle pressure fluctuations, and uneven furnace temperature distribution cause the real boundary conditions to exhibit highly time-varying and nonlinear characteristics. Current simulation systems lack an effective algorithmic mechanism that can transform discrete, multi-source, and noisy sensor data collected from the industrial field into continuously changing boundary condition parameters in the physical model.

[0004] The lack of this mapping mechanism leads to biases in the input of the simulation model, resulting in a serious disconnect between the calculated results on microstructure evolution and macroscopic mechanical response and the actual physical state.

[0005] To address this, a digital twin-based intelligent control system for heat treatment of drill pipes is proposed. Summary of the Invention

[0006] The purpose of this invention is to provide an intelligent control system for heat treatment of drill pipes based on digital twins. This system generates a curve containing the optimal heating power and dynamic spraying parameters by executing a multi-objective optimization algorithm in a multi-dimensional process parameter space. Production optimization is performed through a process instruction set.

[0007] To achieve the above objectives, the present invention provides the following technical solution:

[0008] A digital twin-based intelligent control system for heat treatment of drill pipes includes:

[0009] The multi-source sensing module collects environmental parameters, medium parameters, and workpiece state parameters from the heat treatment production line. After spatiotemporal alignment and noise reduction preprocessing, it outputs a time-series input vector that reflects the physical working conditions.

[0010] The coupled simulation module receives the time-series input vector and reconstructs the heat exchange boundary conditions using the included medium parameters. Combining the heat exchange boundary conditions, it drives the integrated temperature field submodule, microstructure evolution submodule, and stress field submodule to perform strongly coupled calculations based on the phase transformation induced strain mechanism, and outputs a virtual state matrix containing the full-section temperature distribution, grain size distribution, phase transformation volume fraction distribution, and stress distribution.

[0011] The decision optimization module, based on the virtual state matrix, executes a multi-objective optimization algorithm in the multi-dimensional process parameter space to find the optimal balance point of strength and toughness index under the premise of satisfying microstructure and stress constraints, and generates a process instruction set containing the optimal heating power curve and dynamic spraying parameters.

[0012] The execution control module receives and parses the process instruction set, converts it into equipment control signals to drive the heating and spraying mechanisms, and optimizes the production of the next batch.

[0013] The environmental parameters include furnace temperature; the medium parameters include cooling medium flow rate and spray quenching pressure; the workpiece state parameters include workpiece surface temperature and timestamp.

[0014] The spatiotemporal alignment and noise reduction preprocessing includes: for electromagnetic interference present in industrial sites, a multi-scale wavelet transform algorithm is applied to decompose the collected furnace temperature, flow rate and pressure time series data, and high-frequency noise components are removed by setting an adaptive threshold while retaining low-frequency trend terms; for data asynchrony caused by inconsistent sampling frequencies of different types of sensors, a unified time reference axis is constructed, and a cubic spline interpolation algorithm is used to reconstruct and align the sparse sampled data, so that the time steps of all state variables are aligned with the calculation cycle of the coupled simulation module.

[0015] The heat exchange boundary condition reconstruction process includes:

[0016] The coupled simulation module has a built-in boundary mapping unit based on fluid dynamics, which is used to transform the medium parameters in the time-series input vector into microscopic physical boundary conditions. The boundary mapping unit establishes a nonlinear mapping function, which takes the cooling medium flow rate and quenching pressure as independent variables and the comprehensive heat transfer coefficient of the workpiece surface as the dependent variable. The nonlinear mapping function not only includes the basic natural convection term, but also introduces a forced convection weight term dominated by flow rate and an impact kinetic energy weight term dominated by pressure.

[0017] At each simulation time step, the boundary mapping unit reads the current cooling medium flow rate and spray quenching pressure values, calculates the dynamically changing comprehensive heat transfer coefficient, and assigns the comprehensive heat transfer coefficient to the surface mesh nodes of the drill pipe geometry model to simulate the time-varying process of the cooling intensity of the pipe surface caused by hydraulic fluctuations in the production line and changes in nozzle status.

[0018] The temperature field submodule runs a transient heat conduction solver, using the comprehensive heat transfer coefficient as the boundary drive, to calculate the instantaneous temperature of each discrete node in the radial direction of the drill pipe from the surface to the core;

[0019] During the calculation process, a parameter observer based on the extended Kalman filter algorithm is run synchronously. The parameter observer sets the comprehensive heat transfer coefficient as the state variable to be estimated, reads the workpiece surface temperature in the time-series input vector, calculates the residual between the predicted surface temperature and the calculated surface temperature, and corrects the comprehensive heat transfer coefficient through the gain matrix.

[0020] The calibrated transient heat conduction solver outputs a temperature distribution vector containing the nodal temperature values ​​of the entire cross section.

[0021] The microstructure evolution submodule performs tissue calculations to drive the heating and cooling stages;

[0022] Using the temperature distribution vector as input, the calculation logic of the high-temperature region and the phase transformation region is activated in parallel; in the high-temperature region, the average equivalent diameter of the austenite grains is accumulated and calculated using a grain growth kinetic model based on thermal activation energy, and a grain size distribution vector is generated.

[0023] In the phase transition region, the non-diffusion-type phase transition kinetic equation is used to quantify the exponential relationship between the martensite transformation amount and the supercooling, and generate a phase transition volume fraction distribution vector.

[0024] The stress field submodule constructs a mechanism for the transmission of microstructure to macroscopic mechanics, decomposing the total strain rate into elastic strain rate, plastic strain rate, thermal strain rate, phase transformation volume expansion strain rate, and phase transformation induced strain rate.

[0025] Among them, the phase transformation volume expansion strain rate is proportional to the rate of change of the phase transformation volume fraction and the volume expansion coefficient of the material; the phase transformation induced plastic strain rate is proportional to the rate of change of the phase transformation volume fraction and the current deviatoric stress tensor.

[0026] When solving the thermo-elastic-plastic constitutive equation, the volume expansion caused by the change in crystal structure is superimposed on the geometric equation as an intrinsic strain source term to simulate the internal force confrontation and stress relaxation process caused by the asynchronous phase transformation between the surface and the core. The stress distribution vector containing the axial, tangential and radial stress values ​​of the entire cross section is calculated. A virtual state matrix is ​​constructed based on the temperature distribution vector, grain size distribution vector, phase transformation volume fraction distribution vector and stress distribution vector.

[0027] The decision optimization module performs the following operations during the generation of the process instruction set:

[0028] A macro-micro performance mapping model is constructed, using the phase transformation volume fraction distribution, grain size distribution and stress distribution in the virtual state matrix as input features to predict the tensile strength and impact energy of the drill pipe after heat treatment.

[0029] Subsequently, a multi-objective evolutionary algorithm was applied to iteratively search within the feasible region of process parameters, with optimization objectives covering maximizing strength and toughness and minimizing residual stress.

[0030] After the algorithm generates the Pareto optimal solution set, it selects the best process scheme according to the preset engineering weights and converts it into heating temperature curve setting data and spray flow rate and pressure time series data to form the process instruction set.

[0031] Compared with the prior art, the beneficial effects of the present invention are as follows:

[0032] 1. This invention utilizes the window period between batch production runs for in-depth computation, allowing the system ample time to run complex, strongly coupled physical models. This accurately reproduces the entire physical mechanism, including austenite grain growth, martensitic phase transformation volume expansion, and internal force counteraction. By converting the high-precision inversion results of the previous batch into feedforward control commands for the next batch, the computational bottleneck of online calculations is successfully avoided without sacrificing the fidelity of the physical model. This achieves true closed-loop control in industrial production, significantly improving the scientific rigor and effectiveness of process optimization.

[0033] 2. This invention establishes a two-layer adaptive correction mechanism. During a single production process, the extended Kalman filter algorithm is used to fuse sensor data, eliminate measurement noise, and correct state variables online to ensure the accuracy of the current observation. After the batch production is completed, a regression task is constructed using the full-process data, and the weight coefficients in the boundary mapping function are updated in reverse using the least squares method. This endows the system with the ability to self-evolve, enabling it to automatically capture and adapt to environmental changes such as nozzle aging and medium deterioration. It can maintain high prediction accuracy of the model throughout its entire life cycle without frequent manual calibration, greatly reducing the system's maintenance costs and application threshold.

[0034] 3. This invention constructs a high-precision macro-micro performance mapping model, transforming invisible microstructural features into visible macroscopic mechanical property predictions. Based on this, combined with a multi-objective evolutionary algorithm, it can perform automated iterative searches within the parameter space covering the entire process of heating, heat preservation, and quenching. The algorithm can accurately identify the Pareto optimal front among strength, toughness, and residual stress, automatically balance the contradiction between surface hardness and core toughness, and generate the optimal cooling path that minimizes residual stress. This significantly improves the batch stability and yield rate of drill pipe materials, providing core technical support for the stable mass production of high-end oil well pipe materials. Attached Figure Description

[0035] Figure 1 This is a schematic diagram of the intelligent control system for heat treatment of drill pipes based on digital twins according to the present invention.

[0036] Figure 2 This is a logical schematic diagram of the intelligent control system for heat treatment of drill pipes based on digital twins according to the present invention.

[0037] Figure 3 This is a schematic diagram of the coupled simulation module of the present invention. Detailed Implementation

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

[0039] Example 1:

[0040] This invention proposes an intelligent control system for heat treatment of drill pipes based on digital twins. The structure of the system is as follows: Figure 1 As shown, the data logic of the system is as follows: Figure 2 As shown; including:

[0041] The multi-source sensing module collects environmental parameters, medium parameters, and workpiece state parameters from the heat treatment production line. After spatiotemporal alignment and noise reduction preprocessing, it outputs a time-series input vector that reflects the physical working conditions.

[0042] The coupled simulation module receives the time-series input vector and reconstructs the heat exchange boundary conditions using the included medium parameters. Combining these boundary conditions, it drives the integrated temperature field submodule, microstructure evolution submodule, and stress field submodule to perform strongly coupled calculations based on the phase transition-induced strain mechanism, outputting a virtual state matrix containing the full-section temperature distribution, grain size distribution, phase transition volume fraction distribution, and stress distribution. The structure of the coupled simulation module is as follows: Figure 3 As shown;

[0043] The decision optimization module, based on the virtual state matrix, executes a multi-objective optimization algorithm in the multi-dimensional process parameter space to find the optimal balance point of strength and toughness index under the premise of satisfying microstructure and stress constraints, and generates a process instruction set containing the optimal heating power curve and dynamic spraying parameters.

[0044] The execution control module receives and parses the process instruction set, converts it into equipment control signals to drive the heating and spraying mechanisms, and optimizes the production of the next batch.

[0045] The environmental parameters include furnace temperature; the medium parameters include cooling medium flow rate and spray quenching pressure; the workpiece state parameters include workpiece surface temperature and timestamp.

[0046] The spatiotemporal alignment and noise reduction preprocessing includes: for electromagnetic interference present in industrial sites, a multi-scale wavelet transform algorithm is applied to decompose the collected furnace temperature, flow rate and pressure time series data, and high-frequency noise components are removed by setting an adaptive threshold while retaining low-frequency trend terms; for data asynchrony caused by inconsistent sampling frequencies of different types of sensors, a unified time reference axis is constructed, and a cubic spline interpolation algorithm is used to reconstruct and align the sparse sampled data, so that the time steps of all state variables are aligned with the calculation cycle of the coupled simulation module.

[0047] Environmental parameter acquisition: Thermocouples are installed in each temperature zone of the heating furnace to collect the temperature of the furnace chamber.

[0048] Medium parameter acquisition: Install electromagnetic flow meters and pressure transmitters in the pipeline of the quenching tank to collect the instantaneous flow rate and nozzle injection pressure of the cooling medium, respectively.

[0049] Workpiece status parameter acquisition: Install dual-color infrared thermometers at the furnace opening and quenching station to collect temperature sequences of key points on the workpiece surface; at the same time, record the time points when the workpiece enters or leaves each process, and synchronously obtain the raw material attribute parameters of the workpiece.

[0050] To address the issues of strong electromagnetic interference and inconsistent sampling frequencies among different sensors in industrial settings, the following steps are employed for spatiotemporal alignment and noise reduction:

[0051] Multi-scale wavelet transform noise reduction: For non-stationary industrial signals, a multi-scale wavelet transform algorithm is used to decompose the acquired data sequence. The system automatically calculates the noise level estimate at each scale and sets an adaptive threshold to remove high-frequency detail components (i.e., noise components) with amplitudes less than the threshold, and reconstructs low-frequency approximate components (i.e., trend terms). This effectively filters out environmental electromagnetic noise while preserving key signal abrupt changes such as the sudden drop in temperature upon entering water.

[0052] Spatiotemporal alignment and reconstruction: Construct a unified time reference axis; For sparse data with low sampling frequency, use cubic spline interpolation algorithm for resampling processing, and fill the time gap by solving spline functions with continuous second derivatives, so that the time steps of all state variables are strictly aligned with the calculation cycle of subsequent coupled simulation modules, ensuring data synchronization.

[0053] This invention solves the problem of poor data quality in industrial settings. First, compared with traditional filtering methods, multi-scale wavelet transform can effectively preserve key transient features during heat treatment and avoid signal distortion caused by excessive smoothing, providing a high-fidelity input source for digital twin models. Second, the spatiotemporal alignment mechanism based on cubic spline interpolation eliminates the data asynchrony between multi-source heterogeneous sensors, solves the confusion of causal relationships caused by time phase differences, and ensures the accuracy and physical consistency of boundary conditions in subsequent strongly coupled calculations.

[0054] The heat exchange boundary condition reconstruction process includes:

[0055] The coupled simulation module has a built-in boundary mapping unit based on fluid dynamics, which is used to transform the medium parameters in the time-series input vector into microscopic physical boundary conditions. The boundary mapping unit establishes a nonlinear mapping function, which takes the cooling medium flow rate and quenching pressure as independent variables and the comprehensive heat transfer coefficient of the workpiece surface as the dependent variable. The nonlinear mapping function not only includes the basic natural convection term, but also introduces a forced convection weight term dominated by flow rate and an impact kinetic energy weight term dominated by pressure.

[0056] At each simulation time step, the boundary mapping unit reads the current cooling medium flow rate and spray quenching pressure values, calculates the dynamically changing comprehensive heat transfer coefficient, and assigns the comprehensive heat transfer coefficient to the surface mesh nodes of the drill pipe geometry model to simulate the time-varying process of the cooling intensity of the pipe surface caused by hydraulic fluctuations in the production line and changes in nozzle status.

[0057] The core of the coupled simulation module lies in the determination of boundary conditions. In the spray quenching process, the boundary conditions are not static constants, but dynamic variables that fluctuate drastically with the fluid state. Therefore, the system has a built-in boundary mapping unit based on fluid dynamics.

[0058] As a core of logical computation, the boundary mapping unit acts as a translator between macroscopic production line data and microscopic physical simulation. Its main function is to receive the cooling medium flow rate and spray quenching pressure values ​​from the multi-source sensing module, and through the built-in nonlinear mapping function, transform the discrete process parameters into a comprehensive heat transfer coefficient distribution field covering the entire surface of the workpiece that can drive the solution of the heat conduction equation.

[0059] The physical construction logic of the nonlinear mapping function: In order to accurately describe the complex fluid-structure interaction heat transfer behavior during the spraying process, the boundary mapping element establishes a nonlinear mapping function that integrates multiple physical mechanisms, and its functional form is as follows:

[0060] ;

[0061] in, The overall heat transfer coefficient; Based on the baseline value of natural convection; , The weighting coefficients are obtained by inverting historical data using the least squares method. The instantaneous flow rate of the cooling medium is normalized data. The normalized data represents the injection pressure from the nozzle. , These two empirical indices are key parameters used to describe the nonlinear growth law of flow rate and heat exchange efficiency (diminishing marginal effect) and the influence of pressure on vapor film rupture. In actual implementation, the initial values ​​are first set with reference to fluid mechanics theory, and then corrected and established by combining standard sample calibration test or CFD simulation data.

[0062] Furthermore, the specific operational procedures for the standard sample calibration test include: measuring the natural cooling curve under no-spray conditions, and calculating the basic natural convection reference value using the lumped parameter method; fixing the injection pressure P, gradually increasing the flow rate Q in 10% increments, and recording the cooling curves under each operating condition; calculating the peak heat transfer coefficient for each operating condition using an inverse heat transfer algorithm, and determining the flow rate weight through logarithmic coordinate regression fitting. With index With a fixed flow rate Q, the pressure P is adjusted according to a gradient, and the pressure weights are determined similarly through fitting. With index Orthogonal experiments were conducted in the full parameter space to correct for the coupling effect between the two and finally determine the initial function parameters.

[0063] This function overcomes the drawbacks of traditional empirical formulas that only consider a single variable, and is constructed using a weighted superposition method. Specifically, it includes the following three core components:

[0064] Part 1: Basic Natural Convection Reference Values;

[0065] This section represents the basic heat exchange capacity between the workpiece surface and the surrounding medium under conditions of no external fluid impact or fluid stillness. It is mainly determined by the specific heat capacity and thermal conductivity of the medium, as well as the temperature difference between the workpiece surface and the medium. This ensures that the basic heat exchange based on natural cooling or radiative heat dissipation can still be calculated before or after the spray system is turned on or off, thus ensuring the continuity of the simulation process on the time axis.

[0066] Part Two: Flow-Driven Forced Convection Weighting Term;

[0067] This section aims to reflect the stripping effect of macroscopic fluid flow on the thermal boundary layer. Physically, the flow rate of the cooling medium directly determines the fluid velocity flowing over the workpiece surface, thus affecting the Reynolds number. This term describes the nonlinear growth law of flow rate and heat transfer coefficient through a power-law relationship: as the flow rate increases, the turbulence intensity of the fluid increases, the laminar sublayer close to the workpiece surface becomes thinner, and the thermal resistance is significantly reduced. This term uses an exponential parameter to simulate that the larger the flow rate, the faster the heat is removed.

[0068] Part Three: Pressure-dominated impact kinetic energy weighting term;

[0069] This section specifically addresses the vapor film barrier problem during the initial stage of high-temperature quenching, known as the Leidenfrost effect. At high temperatures, a heat-insulating vapor film easily forms on the workpiece surface, leading to a sharp drop in cooling efficiency. This weighting term maps the quenching pressure to the kinetic energy of the droplets or jet impacting the workpiece surface. When the impact force generated by the pressure exceeds the surface tension threshold of the vapor film, this term increases significantly, simulating the physical process of the high-pressure jet piercing and shattering the vapor film, allowing the liquid medium to directly contact the high-temperature metal surface. This mechanism ensures that the model can accurately capture abrupt changes in cooling modes caused by pressure fluctuations, namely the transition from film boiling to nucleation boiling.

[0070] Dynamic computation and spatiotemporal mapping process: Within each discrete time step of the simulation computation, the boundary mapping unit performs the following refinement operations:

[0071] Parameter reading and interpolation: First, read the current cooling medium flow rate and spray quenching pressure values ​​from the time-series input vector. If the current simulation time step is between two sensor sampling points, the unit will automatically call the interpolation algorithm to obtain the instantaneous estimated value at that moment.

[0072] Non-uniform calculation of local heat transfer coefficient; considering that the nozzles on the industrial spray ring are usually distributed in an array, the cooling intensity of the workpiece surface is not uniform in space, for example, the cooling is strong directly below the nozzle and weak between the nozzles; the boundary mapping unit introduces a spatial distribution factor, and combined with the above nonlinear mapping function, not only calculates the global heat transfer coefficient reference value that changes with time, but also, based on the geometric coordinates of the drill pipe surface, superimposes the local enhancement coefficient of the nozzle jet coverage area to generate a heat transfer coefficient distribution field that changes with position;

[0073] The process involves assigning and driving values ​​to the mesh nodes. Finally, the calculated comprehensive heat transfer coefficient, which includes spatiotemporal characteristics, is precisely assigned to each surface mesh node of the drill pipe geometry model using a mapping algorithm. This ensures that the boundary conditions used by the heat conduction solver to calculate the temperature distribution at the next moment are no longer uniform ideal values, but rather realistically reflect the complex physical boundaries of hydraulic fluctuations in the production line, nozzle injection angles, and jet impact forces. This achieves a high-fidelity simulation of the time-varying process of cooling intensity on the pipe surface.

[0074] Specifically, the boundary mapping unit performs the following spatial mapping steps when calculating the overall heat transfer coefficient:

[0075] First, a cylindrical coordinate system and a nozzle array coordinate model of the spray device are established on the drill pipe surface. Second, based on the current cooling medium flow rate and spray quenching pressure, a global heat transfer coefficient baseline value is calculated using a nonlinear mapping function. Subsequently, a spatial distribution factor function is constructed, which uses the straight-line distance between any grid node on the drill pipe surface and the nearest nozzle jet center as the independent variable to describe the physical law of cooling intensity attenuation with increasing distance. Finally, all grid nodes on the drill pipe surface are traversed, and the local enhancement coefficient of each node affected by the surrounding nozzle jet is calculated. This coefficient is then weighted and superimposed with the global heat transfer coefficient baseline value to generate a non-uniform heat transfer coefficient distribution field covering the entire surface and exhibiting peak and trough distribution characteristics with geometric position.

[0076] Step 1: Constructing the geometric mapping relationship; First, the drill pipe surface is discretized into a set of finite element mesh nodes, and the axial coordinates and circumferential angles of each node in the cylindrical coordinate system are determined. Simultaneously, the geometry of the physical spray ring is reconstructed in digital space, and the projected coordinate positions of all nozzle jet impact points on the workpiece surface are recorded. The system pre-sets the geometric parameters for uniform circumferential distribution and multi-row axial distribution of nozzles.

[0077] Step 2: Define the spatial distribution factor function; define a jet influence function to describe the attenuation of the cooling capacity of a single nozzle. Physical experiments show that cooling is strongest directly below the nozzle, and the cooling capacity gradually weakens with increasing distance from the center. This embodiment uses a Gaussian attenuation model to describe this phenomenon, that is, the cooling capacity decreases exponentially with increasing distance from the nozzle center. This function mainly includes two physical parameters: one is the effective coverage radius of the jet controlled by the nozzle orifice diameter and the jet diffusion angle, and the other is the center enhancement amplitude coefficient relative to the reference value;

[0078] Specifically, spatial distribution factor function for:

[0079] ;

[0080] in, The Euclidean distance from the grid node to the nearest nozzle jet center; The effective coverage radius of the jet; Enhance the amplitude coefficient around the center;

[0081] The central enhancement amplitude coefficient and the effective coverage radius of the jet were obtained through fixed-point injection calibration tests on instrumented standard drill pipe specimens. Specifically, thermocouples were embedded at the jet impact center and edge on the specimen surface to collect cooling curves, and the peak heat transfer coefficient at the impact center and the reference heat transfer coefficient in the non-impact zone were calculated using an inverse heat transfer algorithm. The central enhancement amplitude coefficient was calculated by subtracting the reference heat transfer coefficient from the peak heat transfer coefficient and then dividing by the reference heat transfer coefficient. Alternatively, computational fluid dynamics simulation can be used to numerically simulate the nozzle jet flow field and extract the wall Nusselt number distribution characteristics for fitting.

[0082] Step 3: Calculate and superimpose the local enhancement coefficient. At each simulation time step, for each grid node on the drill pipe surface, the following calculations are performed: First, search the nozzle set for the nearest effective nozzles to the grid node; typically, only nozzles within the effective coverage radius need to be considered. Second, calculate the local enhancement coefficient. The system linearly superimposes the cooling contributions of all effective nozzles in the neighborhood to the node. This means that the coefficient is larger in the overlapping region of multiple nozzle jets, while it approaches zero in the gap region between nozzles. Finally, generate the final heat transfer coefficient. The system uses a multiplicative superposition method, multiplying the global heat transfer coefficient baseline value by a correction term including the local enhancement coefficient. This calculation logic means that the higher the baseline flow velocity, the more the local cooling differences will be amplified.

[0083] Step 4: Generate a time-varying distribution field; finally, the system outputs a heat transfer coefficient matrix covering the entire surface. As the drill rod spirals forward on the production line, the coordinate positions of the mesh nodes relative to the nozzles change periodically. Boundary mapping elements update this relative positional relationship, thereby simulating spiral-shaped strong and weak cooling zones on the drill rod surface, realistically reproducing the thermal history of the workpiece as it passes through the spray ring.

[0084] This invention realistically reproduces the weak cooling effect at the nozzle gap by introducing a spatial distribution factor. The digital twin model can predict in advance which areas may experience insufficient martensitic transformation due to their location in the injection blind zone, thereby guiding process engineers to adjust the nozzle arrangement density or increase the injection pressure, fundamentally eliminating soft spot defects. Through non-uniform field calculations, the model can capture the extreme stresses in these tiny regions. This enables the system to accurately identify potential risk points for microcrack initiation, avoiding misjudgments caused by averaging calculations and significantly improving the safety boundary of process parameter optimization.

[0085] This invention breaks through the limitations of static boundary conditions in traditional simulations. By introducing a flow-dominated forced convection term and a pressure-dominated impact kinetic energy term, it can accurately describe the uneven local cooling phenomenon caused by nozzle pressure fluctuations or blockages when the drill pipe passes through the spray ring. Secondly, this dynamic mapping mechanism can capture the key physical event of vapor film rupture, namely the process of a sudden increase in heat transfer coefficient caused by high-pressure impact, which is crucial for controlling quenching cracks and hardness uniformity.

[0086] The temperature field submodule runs a transient heat conduction solver, using the comprehensive heat transfer coefficient as the boundary drive, to calculate the instantaneous temperature of each discrete node in the radial direction of the drill pipe from the surface to the core;

[0087] During the calculation process, a parameter observer based on the extended Kalman filter algorithm is run synchronously. The parameter observer sets the comprehensive heat transfer coefficient as the state variable to be estimated, reads the workpiece surface temperature in the time-series input vector, calculates the residual between the predicted surface temperature and the calculated surface temperature, and corrects the comprehensive heat transfer coefficient through the gain matrix.

[0088] The calibrated transient heat conduction solver outputs a temperature distribution vector containing the nodal temperature values ​​of the entire cross section.

[0089] The transient heat conduction solver operates as follows: The temperature field submodule is the core computational engine for coupled simulation. This module first discretizes the drill pipe radially from the outer surface to the inner surface into several concentric ring-shaped mesh nodes based on the finite difference method or finite element method. Within each computation cycle, the solver reads the comprehensive heat transfer coefficient provided by the boundary mapping elements as an external driving condition, and combines it with the thermophysical parameters of the drill pipe material as a function of temperature, including density, specific heat capacity, and thermal conductivity, to call the discretized numerical solution method for the transient heat conduction partial differential equation. This process calculates how heat is transferred between the surface and the core within the current time step, thereby updating the instantaneous temperature values ​​on each mesh node and forming the temperature distribution field across the entire cross-section.

[0090] The parameter observation process based on the Extended Kalman Filter (EKF) algorithm: To eliminate the accumulated model errors caused by material property deviations, boundary condition estimation errors, or numerical approximations, the system integrates a parameter observer based on the EKF algorithm. This observer treats the integrated heat transfer coefficient not merely as an input parameter, but expands it into a state variable to be estimated for the system. The specific online observation and correction process includes the following precise steps:

[0091] Step 1: Prior State Prediction; Before reading the measured data, the parameter observer first uses the thermal conduction physics model to deduce the prior state at the current moment based on the optimal estimated state at the previous moment. This step not only predicts the theoretical temperature value of the drill pipe surface, but also calculates the state error covariance matrix; this covariance matrix quantifies the uncertainty of the model prediction, and it will automatically accumulate the random disturbances inside the system over time according to the process noise model.

[0092] Step 2: Calculation of observation residuals; The system reads the measured surface temperature of the drill pipe collected by the infrared thermometer. The observer compares the predicted surface temperature calculated by the physical model with the measured value from the sensor, and calculates the difference between the two, which is the observation residual. This residual reflects the degree of deviation between the current physical model and the actual physical process.

[0093] Step 3: Kalman Gain Calculation and State Correction; The observer calculates the optimal Kalman gain matrix based on the prediction covariance matrix (representing the model's unreliability) and the measurement noise covariance matrix (representing the sensor's unreliability). The Kalman gain is essentially a weighting coefficient that determines whether the system should trust the model's calculated values ​​or the sensor's measured values ​​more. Subsequently, this gain matrix is ​​used to feed the observation residuals back into the state vector.

[0094] This feedback not only corrects the current temperature distribution, but more importantly, it directly corrects the augmented state variable, the overall heat transfer coefficient. This means that the system uses mathematical methods to calculate, in reverse, what the boundary heat transfer coefficient should be at this moment to produce the currently observed surface temperature; the resulting posterior estimate is the calibrated overall heat transfer coefficient. This improves the accuracy of the overall heat transfer coefficient and ensures the accuracy of subsequent calculations.

[0095] The reverse update strategy for the model between batches: The above process realizes the numerical correction in a single production. In order to make the model more and more accurate, the system also designed a reverse learning mechanism between batches.

[0096] After a complete production batch is completed, the system extracts the posterior integrated heat transfer coefficient sequence calibrated by Kalman filtering throughout the entire process, as well as the corresponding cooling medium flow rate sequence and quenching pressure sequence for that batch. The system uses the least squares method or gradient descent algorithm to construct a regression optimization task. The goal of this task is to refit the nonlinear mapping function in the boundary mapping unit, minimizing the total error between the calculated function value and the posterior calibration value, and inversely solving and updating the weight coefficients of natural convection, forced convection, and impact kinetic energy in the function. The updated weight coefficients are saved to the database as the basic parameters of the boundary mapping unit for the next production batch.

[0097] During long-term operation, even when faced with environmental drift issues that are difficult to quantify, such as the gradual scaling and blockage of quenching nozzles leading to changes in flow characteristics, or the aging of quenching media causing a decrease in cooling capacity, the system can automatically capture these changes and adjust the weight parameters within the model by backtracking data from each batch. This ensures that the digital twin model can maintain a high degree of consistency with the physical production line throughout its entire lifecycle without the need for frequent manual calibration.

[0098] This invention achieves deep fusion of physical models and sensor data by introducing an extended Kalman filter algorithm. It not only filters out high-frequency measurement noise from infrared thermometers caused by water mist interference, but also utilizes surface temperature observation information to deduce the internal structure of the workpiece through physical equations, providing a more accurate core temperature inference than simple measurement. This effectively solves the problem of unmeasurable internal conditions during heat treatment. Through online correction of the heat transfer coefficient and subsequent reverse updating of model parameters, an immune mechanism against environmental changes is established, ensuring the accuracy of temperature field calculations and the reliability of process control under long-term operation.

[0099] The microstructure evolution submodule performs tissue calculations to drive the heating and cooling stages;

[0100] Using the temperature distribution vector as input, the calculation logic of the high-temperature region and the phase transformation region is activated in parallel; in the high-temperature region, the average equivalent diameter of the austenite grains is accumulated and calculated using a grain growth kinetic model based on thermal activation energy, and a grain size distribution vector is generated.

[0101] In the phase transition region, the non-diffusion-type phase transition kinetic equation is used to quantify the exponential relationship between the martensite transformation amount and the supercooling, and generate a phase transition volume fraction distribution vector.

[0102] The microstructure evolution submodule is mainly responsible for simulating the phase transformation behavior and grain growth process of drill pipe materials during the entire heat treatment cycle. The specific calculation logic is divided into a heating and holding stage and a cooling and quenching stage.

[0103] Heating stage: cumulative calculation of austenite grain growth; during the heating and holding stage, the core task of the submodule is to calculate the degree of austenite grain coarsening, which has a decisive impact on the impact toughness of the final product.

[0104] Model construction logic: The system adopts a grain growth kinetics model based on thermal activation energy, which assumes that the migration rate of grain boundaries is controlled by both temperature and time; the model includes a grain growth exponent determined by material composition and a grain growth rate constant determined by temperature.

[0105] Calculation Execution Process: Within each simulation time step, the submodule reads the instantaneous temperature of the current node. If this temperature exceeds the austenitizing transformation temperature, the system begins to perform cumulative grain size calculations. The algorithm differentiates the entire heating process, calculates the grain size increment within each tiny time increment, and adds this increment to the grain size at the previous time step. Through this time integration method, the model can accurately reflect the influence of different heating rates and holding times on the final average equivalent grain diameter, outputting a feature vector containing the grain size distribution across the entire cross-section.

[0106] Cooling stage: Quantification of multiphase structure transformation; In the spray cooling stage, the submodule calculates the transformation ratio of supercooled austenite to low-temperature structures such as martensite and bainite based on the temperature drop process.

[0107] Non-diffusional phase transformation calculation: For the martensitic transformation, which is of most concern in the heat treatment of drill pipe, a non-diffusional phase transformation kinetic equation is used for calculation. This equation not only considers the influence of chemical composition on the martensitic transformation initiation temperature, but also focuses on quantifying the nonlinear exponential relationship between the martensitic transformation amount and the degree of undercooling.

[0108] Dynamic evolution logic: When the node temperature drops below the martensitic transformation initiation temperature, the submodule calculates the volume fraction of transformed martensite at that moment based on the current degree of supercooling, i.e., the difference between the current temperature and the initiation temperature. This calculation process is continuous and dynamic; as the temperature decreases, the martensite content gradually increases until the temperature falls below the transformation termination temperature or reaches room temperature. Finally, the submodule outputs a distribution vector containing the phase transformation volume fraction of each node across the entire cross-section, clearly distinguishing the microstructure differences between the surface and the core.

[0109] The stress field submodule constructs a mechanism for the transmission of microstructure to macroscopic mechanics, decomposing the total strain rate into elastic strain rate, plastic strain rate, thermal strain rate, phase transformation volume expansion strain rate, and phase transformation induced strain rate.

[0110] Among them, the phase transformation volume expansion strain rate is proportional to the rate of change of the phase transformation volume fraction and the volume expansion coefficient of the material; the phase transformation induced plastic strain rate is proportional to the rate of change of the phase transformation volume fraction and the current deviatoric stress tensor, and the analysis is based on normalized data.

[0111] When solving the thermo-elastic-plastic constitutive equation, the volume expansion caused by the change in crystal structure is superimposed on the geometric equation as an intrinsic strain source term to simulate the internal force confrontation and stress relaxation process caused by the asynchronous phase transformation between the surface and the core. The stress distribution vector containing the axial, tangential and radial stress values ​​of the entire cross section is calculated. A virtual state matrix is ​​constructed based on the temperature distribution vector, grain size distribution vector, phase transformation volume fraction distribution vector and stress distribution vector.

[0112] The stress field submodule builds a bridge from microstructural changes to macroscopic mechanical responses. Its core lies in dealing with the complex mechanical behavior of phase transformation-induced strain, and predicting residual stress and deformation by solving the thermo-elastic-plastic constitutive equations.

[0113] Refined decomposition of total strain rate: To accurately describe material deformation under complex conditions, the simple thermoelastic model was abandoned, and the total strain rate of the material was decomposed into five independent physical components for superposition calculation:

[0114] Elastic strain rate: follows Hooke's law, describes the reversible deformation of a material within its yield limit, and is related to the current stress state and the elastic modulus as a function of temperature.

[0115] Plastic strain rate: following the Mises yield criterion and flow law, it describes the irreversible permanent deformation of a material after it exceeds its yield limit.

[0116] Thermal strain rate: caused solely by temperature changes, it describes the volume change of a material due to thermal expansion and contraction, and is controlled by the coefficient of linear expansion.

[0117] Phase transformation volumetric expansion strain rate: This is one of the key terms in strongly coupled calculations. Since the specific volume of martensite is significantly greater than that of austenite, the microstructure transformation is accompanied by volumetric expansion; this strain rate is proportional to the rate of change of the phase transformation volume fraction and the volumetric expansion coefficient of the material, reflecting the macroscopic volumetric effect caused by the change in microstructure.

[0118] Phase transformation-induced plastic strain rate: This mechanism describes how, during a phase transformation, even when the macroscopic stress is below the material's yield strength, significant plastic flow occurs along the deviatoric stress direction due to internal weakening caused by microscopic lattice rearrangement. This strain rate is proportional to the rate of change of the phase transformation volume fraction and the current deviatoric stress tensor.

[0119] Simulation process of internal force resistance and stress relaxation: In the solution process, the stress field submodule uses the calculated phase transformation volume expansion strain and phase transformation induced plastic strain as intrinsic strain source terms, and superimposes them into the geometric equations for finite element solution.

[0120] Asynchronous phase transformation simulation: Due to the rapid cooling of the drill pipe surface and the slow cooling of the core, the surface first undergoes a martensitic phase transformation, resulting in volume expansion. This expansion is restrained by the untransformed core, generating compressive stress. Subsequently, the core undergoes a phase transformation and expansion, but is constrained by the hardened surface, generating tensile stress. The submodule dynamically reproduces this internal force conflict process caused by the asynchronous phase transformation between the surface and the core through step-by-step iterative calculations.

[0121] Stress redistribution calculation: Utilizing the phase transformation-induced plasticity mechanism, the abnormal plastic yielding of the material under high stress was simulated. This yielding effect significantly releases some of the peak stress, i.e., stress relaxation. Finally, the submodule calculates the residual stress values ​​distributed along the drill pipe axis, tangentially, and radially after cooling, forming a complete virtual stress state matrix.

[0122] This invention, by explicitly introducing phase transformation-induced plasticity and volume expansion mechanisms, can accurately simulate the compressive effect of surface martensite formation on the core austenite during quenching, as well as the subsequent stress reversal process. This allows the system to accurately identify the timing and location of tensile stress peaks, thereby precisely predicting potential risk areas for quenching cracks and ensuring process safety. Secondly, by strongly coupling the microstructure evolution with the macroscopic stress field, the virtual state matrix not only includes macroscopic residual stress data but also microscopic grain size and phase composition characteristics. This provides a direct physical basis for the subsequent decision optimization module to evaluate the impact toughness of the drill pipe, enabling direct quantification and optimization of product strength and toughness indicators.

[0123] The decision optimization module performs the following operations during the generation of the process instruction set:

[0124] A macro-micro performance mapping model is constructed, using the phase transformation volume fraction distribution, grain size distribution and stress distribution in the virtual state matrix as input features to predict the tensile strength and impact energy of the drill pipe after heat treatment.

[0125] Subsequently, a multi-objective evolutionary algorithm was applied to iteratively search within the feasible region of process parameters, with optimization objectives covering maximizing strength and toughness and minimizing residual stress.

[0126] After the algorithm generates the Pareto optimal solution set, it selects the best process scheme according to the preset engineering weights and converts it into heating temperature curve setting data and spray flow rate and pressure time series data to form the process instruction.

[0127] Construction of the macro-micro performance mapping model: In order to transform the microscopic physical quantities obtained from simulation into macroscopic mechanical performance indicators of engineering concern, the decision optimization module first constructs a high-precision macro-micro performance mapping model, which acts as a translator from virtual state to real performance.

[0128] Input feature extraction: The system extracts key feature data from the virtual state matrix output by the coupled simulation module, including the volume fraction distribution of martensite and bainite across the entire cross section, the average grain size of austenite, the amount of carbide precipitation, and the statistical characteristics of residual stress peak and distribution gradient.

[0129] Proxy model training: A nonlinear mapping relationship is constructed using support vector machine regression or ensemble learning algorithms. The training dataset for the model comes from two parts: historically accumulated actual production and physicochemical testing data, and calibrated high-fidelity simulation calculation data. Through training, the model can predict the tensile strength, yield strength, elongation, and low-temperature impact energy of drill pipe under a given microstructure at an extremely fast speed, replacing destructive physical experiments. In this embodiment, support vector regression is used to establish the nonlinear mapping model, and radial basis functions with strong nonlinear approximation capabilities are selected as the kernel function. During model training, a grid search method combined with five-fold cross-validation is used to iteratively optimize the penalty coefficient and kernel function parameters within a preset range to find the parameter combination with the minimum generalization error, thereby preventing overfitting or underfitting of the model.

[0130] The multi-objective optimization process based on evolutionary algorithms: After obtaining performance prediction capabilities, the module initiates a multi-objective evolutionary algorithm to find the optimal solution within the feasible region of the process parameters. The specific execution steps are as follows:

[0131] Decision variable coding: The process parameters to be optimized are coded, including the heating power and holding time of each temperature zone of the heating furnace, as well as the time-series control points of the cooling medium flow rate and pressure during the spray quenching stage.

[0132] Population Initialization and Evaluation: An initial population containing several sets of process schemes is randomly generated. For each individual in the population (i.e., each set of process schemes), a coupled simulation module is invoked for rapid deduction, and its fitness is evaluated using a performance mapping model. The fitness function contains three conflicting objectives: maximizing tensile strength, maximizing impact toughness, and minimizing final residual stress.

[0133] Non-dominated sorting and iteration: The algorithm employs a non-dominated sorting mechanism with an elitist strategy, dividing the population into different tiers. Individuals at the Pareto front represent the optimal set of solutions where, in the current state, it is impossible to improve one performance metric without decreasing another.

[0134] Evolutionary operations: The next generation of the population is generated through selection, crossover, and mutation operators. During the iteration process, the algorithm automatically eliminates inferior solutions that lead to excessive grain coarsening or insufficient surface hardness, while retaining and propagating superior solutions that balance strength and toughness with low stress. After hundreds of generations of iterative evolution, the algorithm ultimately outputs a Pareto optimal solution set containing multiple trade-offs.

[0135] Optimal process solution selection and instruction generation: After the algorithm converges, the decision module selects the unique optimal process solution from the Pareto solution set based on the preset engineering weights; for example, for deep well drill pipe, the impact toughness is assigned a higher weight; for ultra-deep well drill pipe, the tensile strength is assigned a higher weight; and the multi-attribute decision method is used to select the unique optimal process solution from the Pareto solution set.

[0136] Subsequently, the system transforms the discrete parameters in the scheme into continuous instructions that the equipment can execute. For example, the optimized temperature control point is fitted to a smooth heating curve setpoint, and the time-series changes in pressure are converted into ladder diagram instructions or recipe data packages for programmable logic controllers, forming the final process instruction set.

[0137] The execution control module is responsible for implementing decisions from the digital world into the physical world, and its core task is to provide feedforward control for the next batch of production.

[0138] Command parsing and safety verification: After receiving the process command set, the module first performs logical parsing and safety boundary verification. If the optimized parameters exceed the physical limits of the equipment, such as heating power exceeding the rated value or touching the safety red line, the system will automatically trigger the fuse mechanism, revert to safe and conservative parameters, and issue an alarm.

[0139] Inter-batch parameter distribution: After verification, during the interval between the end of the current batch production and the start of the next batch production, the module downloads the new temperature setting curve, spray pressure and flow formula to the underlying heating controller and hydraulic actuator through the industrial communication protocol.

[0140] Closed-loop optimization execution: When the next batch of drill pipe enters the production line, the equipment will operate strictly according to this new set of parameters that have been simulated and optimized throughout the entire process. This process constitutes a large closed-loop control at the batch level, ensuring that the process parameters of each batch of products are the optimal solution tailored to the feedback from the previous batch and the current equipment status.

[0141] In this embodiment, the feasible range of process parameters is defined as follows: the temperature setting range for each temperature zone of the heating furnace is 850 degrees Celsius to 980 degrees Celsius to ensure complete austenitization of the material and prevent overheating; the holding time setting range is 30 minutes to 90 minutes; the flow rate setting range for the cooling medium in the spray quenching stage is 50 cubic meters per hour to 200 cubic meters per hour, and the spray quenching pressure setting range is 0.2 MPa to 0.8 MPa. These boundary values ​​are set based on the rated power limit of the production equipment and the heat treatment process window of the material.

[0142] This invention utilizes a multi-objective evolutionary algorithm to automatically find equilibrium points that are difficult for human experience to detect within a vast parameter space covering the entire process of heating, heat preservation, and quenching. While ensuring ultra-high strength, it maximizes the preservation of material toughness and suppresses residual stress by finely adjusting the cooling path, generating an optimal process that balances multiple performance indicators and avoids the blindness of manual trial and error. Secondly, it employs an inter-batch optimization strategy to avoid the conflict between the time-consuming calculations of complex coupled simulations and the fast pace of production lines. The system utilizes the window period between batches for in-depth calculations, significantly improving batch stability and yield, and achieving continuous self-evolution of the heat treatment process.

[0143] Example 2:

[0144] This invention proposes an intelligent control system for heat treatment of drill pipes based on digital twins, comprising:

[0145] The multi-source sensing module collects environmental parameters, medium parameters, and workpiece state parameters from the heat treatment production line. After spatiotemporal alignment and noise reduction preprocessing, it outputs a time-series input vector that reflects the physical working conditions.

[0146] The coupled simulation module receives the time-series input vector and reconstructs the heat exchange boundary conditions using the included medium parameters. Combining the heat exchange boundary conditions, it drives the integrated temperature field submodule, microstructure evolution submodule, and stress field submodule to perform strongly coupled calculations based on the phase transformation induced strain mechanism, and outputs a virtual state matrix containing the full-section temperature distribution, grain size distribution, phase transformation volume fraction distribution, and stress distribution.

[0147] The decision optimization module, based on the virtual state matrix, executes a multi-objective optimization algorithm in the multi-dimensional process parameter space to find the optimal balance point of strength and toughness index under the premise of satisfying microstructure and stress constraints, and generates a process instruction set containing the optimal heating power curve and dynamic spraying parameters.

[0148] The execution control module receives and parses the process instruction set, converts it into equipment control signals to drive the heating and spraying mechanisms, and optimizes the production of the next batch.

[0149] The spatiotemporal alignment and noise reduction preprocessing includes: for electromagnetic interference present in industrial sites, a multi-scale wavelet transform algorithm is applied to decompose the collected furnace temperature, flow rate and pressure time series data, and high-frequency noise components are removed by setting an adaptive threshold while retaining low-frequency trend terms; for data asynchrony caused by inconsistent sampling frequencies of different types of sensors, a unified time reference axis is constructed, and a cubic spline interpolation algorithm is used to reconstruct and align the sparse sampled data, so that the time steps of all state variables are aligned with the calculation cycle of the coupled simulation module.

[0150] Specifically, in this embodiment, the Dobessi 4 wavelet is selected as the basis function for the furnace temperature and flow signals. This basis function has good tight support and regularity, which can effectively balance the smoothness of the signal and the ability to capture abrupt changes, avoiding the loss of key feature information at the moment of water entry during the denoising process. The number of wavelet transform decomposition layers is set to four. This number of layers is based on the Shannon sampling theorem, which aims to subdivide the original signal frequency band into a high-frequency noise region and a low-frequency feature region, ensuring that the high-frequency electromagnetic interference of the sensor can be completely separated into the first and second layer detail components.

[0151] The adaptive threshold is calculated using a general thresholding rule. First, the median of the absolute values ​​of the detail components at each scale is calculated, and this median is divided by a constant to obtain an estimate of the noise standard deviation. Then, the final threshold is calculated by combining the natural logarithm of the signal length. Wavelet coefficients with absolute values ​​below the threshold are considered noise and set to zero, while coefficients with absolute values ​​above the threshold are retained or subjected to soft thresholding.

[0152] When constructing a unified time reference axis, a cubic spline interpolation algorithm with natural boundary conditions is adopted; that is, the second derivative of the interpolation curve is required to be zero at both ends. This setting conforms to the physical fact that the system is in a steady state before and after the heat treatment process, which can prevent non-physical artificial oscillations at the beginning and end of the data sequence and ensure that the time-aligned data is smooth and accurate.

[0153] The coupled simulation module has a built-in boundary mapping unit based on fluid dynamics, which is used to transform the medium parameters in the time-series input vector into microscopic physical boundary conditions. The boundary mapping unit establishes a nonlinear mapping function, which takes the cooling medium flow rate and quenching pressure as independent variables and the comprehensive heat transfer coefficient of the workpiece surface as the dependent variable. The nonlinear mapping function not only includes the basic natural convection term, but also introduces a forced convection weight term dominated by flow rate and an impact kinetic energy weight term dominated by pressure.

[0154] At each simulation time step, the boundary mapping unit reads the current cooling medium flow rate and spray quenching pressure values, calculates the dynamically changing comprehensive heat transfer coefficient, and assigns the comprehensive heat transfer coefficient to the surface mesh nodes of the drill pipe geometry model to simulate the time-varying process of the cooling intensity of the pipe surface caused by hydraulic fluctuations in the production line and changes in nozzle status.

[0155] Furthermore, the boundary mapping unit incorporates a variable cross-section geometry adaptive correction module. For the unique "pipe body-transition zone-thickened end" variable wall thickness structure of the drill pipe, a geometry factor is introduced when calculating the comprehensive heat transfer coefficient. This factor is a function of the axial position z, constructed based on the local wall thickness variation rate and surface curvature of the drill pipe. When the simulation scans to the thickened transition zone at the pipe end, the change in wall thickness gradient is automatically identified, and based on the boundary layer separation theory in fluid mechanics, the forced convection weight term in this region is nonlinearly attenuated and corrected. Simultaneously, a turbulence enhancement coefficient is superimposed on the large planar region at the thickened end. The resulting heat transfer coefficient field exhibits a segmented, differentiated distribution characteristic strongly coupled with the drill pipe's geometric profile in the axial direction.

[0156] Geometric discretization: The drill pipe model is divided along the axial direction into the pipe body region, the transition cone surface region, and the joint thickening region.

[0157] Flow field correction: In the transition cone region, the cooling medium velocity separation is caused by the abrupt change in diameter. When calculating the local Reynolds number, a correction factor is introduced to simulate the local decrease in cooling capacity.

[0158] Heat capacity compensation: At the thickened end, due to the large heat capacity, in order to match the cooling rate, the system will automatically increase the nozzle pressure corresponding to this area during the optimization phase, and introduce a pressure enhancement factor simultaneously during boundary condition reconstruction.

[0159] This solution addresses common quality issues in drill pipe heat treatment, such as insufficient end hardness and excessive body hardness. Traditional models often ignore the disturbance of the flow field caused by geometric abrupt changes, leading to soft spots or cracks in the transition zone. This solution, through geometric adaptive correction, accurately reproduces the complex heat transfer behavior at the variable cross-section, ensuring the uniformity of microstructure transformation throughout the entire length of the drill pipe.

[0160] Furthermore, the boundary mapping unit may also include a helical motion pulse heat transfer model.

[0161] For the drilling rod undergoing high-speed rotation while advancing during quenching, this model abandons the static time-averaged approach and instead defines the heat transfer boundary of the workpiece surface as a pulse function that fluctuates frequently with time and rotation angle. The model calculates the duty cycle of the effective impact zone of the nozzle jet at any point on the drilling rod surface. Within the impact zone, the heat transfer coefficient is taken as the high-pressure impact value; in the non-impact zone (back surface), the heat transfer coefficient instantaneously switches to the flowing liquid film heat transfer value. This generates a dynamic boundary condition that oscillates in a rectangular wave on the time axis and is distributed in a spiral band shape in space.

[0162] This solution significantly improves the ability to predict spiral hardness bands on the surface. Traditional models, which use the average heat transfer coefficient, cannot explain why spiral hardness fluctuations appear on the drill pipe surface. This model, by reproducing the "rotational pulse" effect, can accurately predict "zebra stripe" defects caused by improper matching of rotational speed and feed rate, guiding process engineers to adjust the rotational speed to obtain a seamless, uniformly hardened layer.

[0163] Although embodiments of the invention have been shown and described, it will be understood by those skilled in the art that various changes, modifications, substitutions and alterations can be made to these embodiments without departing from the principles and spirit of the invention, the scope of which is defined by the appended claims and their equivalents.

Claims

1. A smart control system for heat treatment of drill pipes based on digital twins, characterized in that, include: The multi-source sensing module collects environmental parameters, medium parameters, and workpiece state parameters from the heat treatment production line. After spatiotemporal alignment and noise reduction preprocessing, it outputs a time-series input vector that reflects the physical working conditions. The coupled simulation module receives the time-series input vector and reconstructs the heat exchange boundary conditions using the included medium parameters. Combining the heat exchange boundary conditions, it drives the integrated temperature field submodule, microstructure evolution submodule, and stress field submodule to perform strongly coupled calculations based on the phase transformation induced strain mechanism, and outputs a virtual state matrix containing the full-section temperature distribution, grain size distribution, phase transformation volume fraction distribution, and stress distribution. The decision optimization module, based on the virtual state matrix, executes a multi-objective optimization algorithm in the multi-dimensional process parameter space to find the optimal balance point of strength and toughness index under the premise of satisfying microstructure and stress constraints, and generates a process instruction set containing the optimal heating power curve and dynamic spraying parameters. The execution control module receives and parses the process instruction set, converts it into equipment control signals to drive the heating and spraying mechanisms, and optimizes the production of the next batch.

2. The intelligent control system for heat treatment of drill pipe based on digital twin as described in claim 1, characterized in that: The environmental parameters include furnace temperature; the medium parameters include cooling medium flow rate and spray quenching pressure; the workpiece state parameters include workpiece surface temperature and timestamp.

3. The intelligent control system for heat treatment of drill pipes based on digital twins according to claim 1, characterized in that: The spatiotemporal alignment and noise reduction preprocessing includes: for electromagnetic interference present in industrial sites, a multi-scale wavelet transform algorithm is applied to decompose the collected furnace temperature, flow rate and pressure time series data, and high-frequency noise components are removed by setting an adaptive threshold while retaining low-frequency trend terms; for data asynchrony caused by inconsistent sampling frequencies of different types of sensors, a unified time reference axis is constructed, and a cubic spline interpolation algorithm is used to reconstruct and align the sparse sampled data, so that the time steps of all state variables are aligned with the calculation cycle of the coupled simulation module.

4. The intelligent control system for heat treatment of drill pipes based on digital twins according to claim 1, characterized in that: The heat exchange boundary condition reconstruction process includes: The coupled simulation module has a built-in boundary mapping unit based on fluid dynamics, which is used to transform the medium parameters in the time-series input vector into microscopic physical boundary conditions. The boundary mapping unit establishes a nonlinear mapping function, which takes the cooling medium flow rate and quenching pressure as independent variables and the comprehensive heat transfer coefficient of the workpiece surface as dependent variable. The nonlinear mapping function not only includes the basic natural convection term, but also introduces a forced convection weight term dominated by flow rate and an impact kinetic energy weight term dominated by pressure. At each simulation time step, the boundary mapping unit reads the current cooling medium flow rate and spray quenching pressure values, calculates the dynamically changing comprehensive heat transfer coefficient, and assigns the comprehensive heat transfer coefficient to the surface mesh nodes of the drill pipe geometry model to simulate the time-varying process of the cooling intensity of the pipe surface caused by hydraulic fluctuations in the production line and changes in nozzle status.

5. The intelligent control system for heat treatment of drill pipes based on digital twins according to claim 4, characterized in that: The temperature field submodule runs a transient heat conduction solver, using the comprehensive heat transfer coefficient as the boundary drive, to calculate the instantaneous temperature of each discrete node in the radial direction of the drill pipe from the surface to the core; During the calculation process, a parameter observer based on the extended Kalman filter algorithm is run synchronously. The parameter observer sets the comprehensive heat transfer coefficient as the state variable to be estimated, reads the workpiece surface temperature in the time-series input vector, calculates the residual between the predicted surface temperature and the calculated surface temperature, and corrects the comprehensive heat transfer coefficient through the gain matrix. The calibrated transient heat conduction solver outputs a temperature distribution vector containing the nodal temperature values ​​of the entire cross section.

6. The intelligent control system for heat treatment of drill pipes based on digital twins according to claim 5, characterized in that: The microstructure evolution submodule performs tissue calculations to drive the heating and cooling stages; Using the temperature distribution vector as input, the calculation logic of the high-temperature region and the phase transformation region is activated in parallel; in the high-temperature region, the average equivalent diameter of the austenite grains is accumulated and calculated using a grain growth kinetic model based on thermal activation energy, and a grain size distribution vector is generated. In the phase transition region, the non-diffusion-type phase transition kinetic equation is used to quantify the exponential relationship between the martensite transformation amount and the supercooling, and generate a phase transition volume fraction distribution vector.

7. The intelligent control system for heat treatment of drill pipes based on digital twins according to claim 6, characterized in that: The stress field submodule constructs a mechanism for the transmission of microstructure to macroscopic mechanics, decomposing the total strain rate into elastic strain rate, plastic strain rate, thermal strain rate, phase transformation volume expansion strain rate, and phase transformation induced strain rate. Among them, the phase transformation volume expansion strain rate is proportional to the rate of change of the phase transformation volume fraction and the volume expansion coefficient of the material; the phase transformation induced plastic strain rate is proportional to the rate of change of the phase transformation volume fraction and the current deviatoric stress tensor. When solving the thermo-elastic-plastic constitutive equation, the volume expansion caused by the change in crystal structure is superimposed on the geometric equation as an intrinsic strain source term to simulate the internal force confrontation and stress relaxation process caused by the asynchronous phase transformation between the surface and the core. The stress distribution vector containing the axial, tangential and radial stress values ​​of the entire cross section is calculated. A virtual state matrix is ​​constructed based on the temperature distribution vector, grain size distribution vector, phase transformation volume fraction distribution vector and stress distribution vector.

8. The intelligent control system for heat treatment of drill pipes based on digital twins according to claim 1, characterized in that: The decision optimization module performs the following operations during the generation of the process instruction set: A macro-micro performance mapping model is constructed, using the phase transformation volume fraction distribution, grain size distribution and stress distribution in the virtual state matrix as input features to predict the tensile strength and impact energy of the drill pipe after heat treatment. Subsequently, a multi-objective evolutionary algorithm was applied to iteratively search within the feasible region of process parameters, with optimization objectives covering maximizing strength and toughness and minimizing residual stress. After the algorithm generates the Pareto optimal solution set, it selects the best process scheme according to the preset engineering weights and converts it into heating temperature curve setting data and spray flow rate and pressure time series data to form the process instruction set.

Images