A method for predicting erosion characteristics and surface temperature of rough coatings based on adjacency grid

Abstract

This invention discloses a method for predicting the erosion characteristics and temperature distribution of a rough thermal barrier coating surface based on adjacent grids. The method includes establishing a grid model, using a particle enrichment model suitable for thermal barrier coatings, using a thin-wall thermal resistance model that supports the dynamic thickness of the coating, establishing a particle trajectory correction method based on the slope of the adjacent grids, and setting a time magnification factor. The method considers the influence of the coating roughness characteristics in the direction of particle flow on the particle's future collision trajectory. The presence of existing particle trajectories in the opposite direction is determined based on the coating roughness characteristics, thereby accurately predicting the actual collision trajectory of the particles on the coating surface. The particle erosion model and the dynamic thin-wall thermal resistance model are combined to predict the erosion characteristics and temperature distribution of the coating surface at any time in the future. The method of the invention can predict the surface life of thermal barrier coatings to a limited extent, providing data support and guidance for thermal protection system design and equipment maintenance.

Description

Technical Field

[0001] The present invention belongs to the technical field, and in particular relates to a method for predicting the erosion characteristics and surface temperature of a rough coating based on an adjacent grid. Background Art

[0002] As core components of aircraft engines, turbine blades are subject to a range of erosion effects on their surfaces under extreme conditions, including high temperatures, high pressures, and complex thermal environments. These conditions include not only the direct impact of high-temperature combustion gases, but also multiple factors such as particle erosion, alternating thermal loads, and centrifugal loads. With increasing service life, the microscopic roughness features of the coating surface formed by high-temperature oxidation, mechanical wear, and particle impact become increasingly prominent. The evolution of this surface morphology further exacerbates the non-uniform erosion and local temperature distortion of the coating by changing the particle collision trajectory and energy transfer path. Therefore, accurately predicting the coating erosion behavior and the resulting performance degradation becomes particularly important. However, existing numerical simulation methods are mostly based on the assumption of ideal smooth surfaces and fail to couple the geometric constraints of dynamic roughness morphology in the calculation of particle motion trajectories, resulting in systematic deviations in the prediction of erosion rate and thermal load distribution.

[0003] Existing research has shown that surface roughness significantly influences particle erosion behavior and the failure process. When coatings are severely worn or delaminated, surface temperatures can rise dramatically, even exceeding the blade's maximum temperature tolerance, leading to ablation. This directly reduces aerodynamic efficiency and poses a potential threat to engine safety. Therefore, in-depth analysis of the erosion degradation mechanism of rough thermal barrier coatings has important theoretical and practical implications for predicting coating life and revealing blade temperature distribution. Summary of the Invention

[0004] The purpose of the present invention is to provide a method for predicting the erosion characteristics and surface temperature of a rough coating based on adjacent grids, and to provide a method for predicting the erosion characteristics and surface temperature distribution of a rough coating by searching adjacent grids around the impact grid and recalculating the unit normal vector of the coating based on the local coating thickness, so as to solve the technical problem of the lack of a method for predicting the erosion characteristics and surface temperature of a rough thermal barrier coating.

[0005] To achieve the above objectives, the present invention adopts the following technical solution: a method for predicting the erosion characteristics and surface temperature distribution of rough coatings based on adjacent grids, comprising:

[0006] Step 1: Scan the laboratory coating sample using an electron microscope to obtain coating thickness point cloud data at a preset sampling rate and simultaneously obtain coating roughness;

[0007] Step 2: Establish a mesh model: Determine the mesh size based on the shape and size of the erosion surface and the electron microscope sampling rate. Use ICEM software to establish a physical model of the erosion surface and divide it into hexahedral meshes.

[0008] Step 3: Use the erosion surface grid coordinates in the erosion surface physical model to match the position information contained in the coating thickness point cloud data, and use UDM (user-defined memory) to store the coating thickness information for later use;

[0009] Step 4: Establish a particle trajectory correction method based on the slope of the adjacent grid: When a particle collides with a grid on the erosion surface, find the four adjacent grids of the hexahedral grid unit. At the same time, define the two adjacent grids in the same direction as the particle flow as slope correction units, and the two adjacent grids in the opposite direction as "shadow area" judgment units. Use the coating thickness information stored in the slope correction unit and the coating thickness at the impact point to recalculate the coating normal at the impact point. Use the coating thickness information of the "shadow area" judgment unit and the coating thickness at the impact point to determine whether the particle trajectory exists. The judgment method is to detect whether the coating stored in the grid unit opposite to the particle flow direction obscures the particle incident trajectory. If obscured, it means that the area is a shadow area, and the effect of the particle on the coating is canceled;

[0010] Step 5: Compile the particle trajectory correction method UDF (user-defined function) established in step 4, compile the particle erosion model UDF suitable for thermal barrier coating erosion prediction, and compile the dynamic thin-wall thermal resistance model UDF suitable for coating thickness variation;

[0011] Step 6: Import the mesh model into Ansys-Fluent software, set the numerical simulation boundary conditions, load the particle trajectory correction method UDF, particle erosion model UDF, and dynamic thin-wall thermal resistance model UDF, set the iterative hook function, and input the amplification factor into the hook function;

[0012] Step 7: Perform steady-state calculations on the computational domain. Once the flow field is stable and filled with particles, activate the particle trajectory correction method, particle erosion model, and dynamic thin-wall thermal resistance model. The particle erosion model is used to calculate the coating erosion amount. Combined with the mesh area at the particle erosion site, the coating residual thickness and erosion rate are calculated.

[0013] Step 8: Input the residual thickness of the coating in step 7 into the dynamic thin-wall thermal resistance model to obtain the surface temperature distribution of the coating;

[0014] Step 9: Determine the current surface temperature of the coating by obtaining the coating surface temperature distribution in step 8, and then calculate the coating insulation temperature or the cooling efficiency of the thermal barrier coating coating area based on the current surface temperature of the coating.

[0015] Furthermore, in step 1, the scanning accuracy of the electron microscope on the coating can be adjusted autonomously according to actual needs, and the coating point cloud information obtained by scanning needs to be converted into actual coating thickness information.

[0016] Furthermore, in step 2, the mesh model needs to meet a minimum mesh quality higher than 0.25, and the mesh division size needs to be coordinated with the point cloud data interval (sampling rate).

[0017] Furthermore, in step 3, the coating thickness point cloud data needs to be matched with the grid model position, and the matching determination method is:

[0018] ;

[0019] in, is the position coordinate of the model, is the location information contained in the point cloud data, is the device sampling rate in μm; when the discriminant is true, the coating thickness information is stored in the UDM of the local grid. At the same time, if the model meshing used in the calculation does not require high precision, its grid size may be larger than the sampling rate. In this case, the above discriminant still applies. When the scanned coating thickness information is missing, this prediction method is still feasible, and the coating thickness information can be unified and stored in the UDM. On the other hand, in actual use, the model area to be simulated and predicted may be larger than the sample. In this case, it is necessary to use the coating roughness characteristic information combined with the autocorrelation function to derive coating thickness distribution information suitable for the simulation model.

[0020] Furthermore, in step 4, the core method for finding adjacent meshes is to call the mesh surface sweeping macro and mesh cell surface search macro built into Ansys-Fluent to identify the adjacent mesh surfaces that collide with the mesh surface. The current collision mesh surface is on the boundary surface and has a corresponding mesh unit, namely c0. The sweep macro is called to find six mesh surfaces on unit c0, one of which is the collision mesh surface. The corresponding c0_new_i and c1_new_i units are found based on the remaining five mesh surfaces, where i∈(0-4). When the mesh surfaces on c0_new_i and c1_new_i are on the same model surface as the collision mesh surface, the adjacent mesh surfaces of the collision mesh surface are considered to have been found. In principle, there are four adjacent mesh surfaces, but when the collision mesh surface is on the model boundary, there may be fewer than four adjacent mesh surfaces. Once the adjacent mesh surfaces that meet the requirements are found, they are divided into slope correction units in the same direction as the particle and "shadow area" judgment units in the opposite direction of the particle according to the angle with the particle velocity. The coating thickness information at the collision location and the coating thickness at the adjacent mesh are used to refit the rough coating plane, and the plane normal vector is calculated. The particle collision trajectory is corrected based on this normal vector. When the coating in the opposite direction of the particle flow obscures the incident trajectory of the particle, the influence of the particle is eliminated and the path is deleted.

[0021] Furthermore, in step 5, the particle erosion model suitable for thermal barrier coating erosion prediction is:

[0022] ;

[0023] ;

[0024] in: ;

[0025] in, is the total erosion amount, is the particle mass, is the collision angle, is the collision velocity, is the critical angle, is the critical speed, is the cutting erosion coefficient, is the deformation erosion coefficient, where BB is the mass loss of the coating under high-angle erosion, AA is the cutting amount when the collision angle is less than the critical angle, and CC is the cutting amount when the collision angle is greater than the critical angle;

[0026] At particle collision angle Less than critical angle When the collision angle is Greater than or equal to the critical angle When the particle collision speed is Less than critical speed When the particles do not cause deformation and quality loss, the cutting erosion coefficient , deformation erosion coefficient Affected by the physical properties of particles and coatings;

[0027] The critical angle is related to the coating temperature:

[0028] ;

[0029] in is the coating temperature;

[0030] The critical speed is related to the physical properties of particles and coatings:

[0031] ;

[0032] in, is the yield stress of the coating; is the particle density; q p Poisson's ratio of particles; is the Young's modulus of the granular material; q m is the Poisson’s ratio of the target; The Young's modulus is the one that takes into account the porosity and Poisson's ratio of the target material; and the cutting erosion coefficient is Deformation erosion coefficient Expressed as:

[0033] ;

[0034] ;

[0035] in is the shape factor of the particle, is the out-of-plane Vickers hardness, is the in-plane Vickers hardness; and is the temperature dependence;

[0036] The temperature dependence and Described as:

[0037] ;

[0038] ;

[0039] H v is the Vickers hardness of the coating, V p is the particle velocity, Is the true density of the coating, which is different from the relative density of the coating The relationship is ,in The density of dense YSZ (yttrium zirconium oxide) thermal barrier coating is 5650 kg / m 3 , is the shape factor, 、 、 as well as Fit parameters for the temperature dependence.

[0040] Furthermore, in step 5, the porosity of the coating is:

[0041] ;

[0042] ;

[0043] in, is the internal porosity, is the inter-column porosity, corresponding to, is the internal relative density, is the relative density between columns, is the relative density of the coating;

[0044] The UDF of the coating's Young's modulus is:

[0045] ;

[0046] in, is the density of the thermal insulation coating, is the coating Poisson’s ratio, E m is the Young's modulus of the omnidirectional dense coating and E m =210Gpa;

[0047] When the operating temperature is greater than 800°C, the UDF of the change in Vickers hardness and service time is:

[0048] ;

[0049] ;

[0050] When the operating temperature is from room temperature to 800℃, the UDF of the change of Vickers hardness and service time is:

[0051] ;

[0052] ;

[0053] in, 、 、 refers to the in-plane fitting parameters, and 、 、 refers to the out-of-plane fitting parameters, and It is the hardness value estimated by the Vickers hardness calculation formula when working at a higher temperature.

[0054] Furthermore, in step 5, the calculation core of the dynamic thin-wall thermal resistance model is to equate the thermal barrier coating with different thicknesses and consistent thermal conductivity to an equivalent coating with a thickness equal to the initial coating thickness and a varying thermal conductivity, and to specify the equivalent coating thickness as the initial coating thickness H.

[0055] The calculation formulas for the normal and tangential thermal conductivity of the coating are derived based on the definition formula of thermal resistance;

[0056] Among them, thermal resistance The expression is:

[0057] ;

[0058] in is the length of the heat transfer path, is the thermal conductivity, is the cross-sectional area perpendicular to the direction of heat transfer;

[0059] When calculating the equivalent value of the thermal conductivity of the coating in the normal direction, the coating thickness changes while the bottom area remains unchanged, and its thermal resistance remains unchanged. Therefore, the expression of the thermal conductivity of the coating in the normal direction is derived as follows:

[0060] ;

[0061] When calculating the equivalent value of the thermal conductivity of the coating in the tangential direction, the coating thickness remains unchanged while the side area decreases, and its thermal resistance remains unchanged. Therefore, the expression of the thermal conductivity of the coating in the tangential direction is derived:

[0062] ;

[0063] in, is the normal equivalent thermal conductivity, is the tangential equivalent thermal conductivity, is the residual thickness of the coating.

[0064] Furthermore, in step 6, the iteration hook function specifies the service time represented by each iteration. , and add the time magnification factor , as the Ansys-Fluent software is gradually iterated, the service time gradually accumulates, where the service time is described as:

[0065] ;

[0066] in, For the length of service, is the length of service, is the time magnification factor, is the number of iterations, is the service time represented by a single iteration.

[0067] Furthermore, in step 7, the expression of the erosion rate is:

[0068] ;

[0069] in is the erosion rate, is the amount of coating erosion caused by the current particles, is the cross-sectional area of ​​the current grid perpendicular to the heat transfer direction, is the mass of the collision particle, Length of service;

[0070] The expression of the residual thickness of the coating is: ;

[0071] in is the residual thickness of the coating, is the initial coating thickness, is the amount of coating erosion caused by the current particles, is the true density of the coating.

[0072] Furthermore, in step 8, the dynamic thin-wall thermal resistance model is loaded using the biaxial anisotropic material thermal conductivity model.

[0073] Furthermore, in step 9, the cooling efficiency of the thermal barrier coating coating area is The calculation formula is:

[0074] ;

[0075] in, and are the gas inlet temperature and the cold air inlet temperature, is the outer surface temperature of the coating;

[0076] Coating insulation temperature The calculation formula is:

[0077] ;

[0078] Where, is the outer surface temperature of the coating, is the inner surface temperature of the coating.

[0079] The present invention provides a method for predicting the erosion characteristics and surface temperature of a rough coating based on an adjacent grid, which has the following beneficial effects:

[0080] (1) The prediction method of the present invention achieves accurate and rapid prediction of the erosion characteristics of rough coatings by utilizing the UDF sweep macro built into Ansys-Fluent, and develops a calculation method for accurately identifying the actual unit normal vector of the coating at the impact position, while determining whether the position is in the "shadow area". This method has a high margin. When the impact grid is located at the model boundary, the number of adjacent grids is unpredictable, but it can still be effectively judged. At the same time, this method has a high accuracy rate, and the effective recognition rate has reached 95.6% after testing. This method successfully achieves accurate and rapid prediction of the erosion morphology of rough coatings and simulates the particle erosion characteristics of the rough coating surface. The method of the present invention is only applied in steady-state flow fields, and the simulation value before stabilization has no effect on the final numerical calculation results. At the same time, the application of this method avoids the tedious small-scale coating physical modeling and meshing, has low requirements for computing resources, saves labor costs, and better meets actual engineering needs.

[0081] (2) The prediction method of the present invention uses a particle erosion model suitable for thermal barrier coatings, taking into account the properties of the thermal barrier coating porosity, Young's modulus, and Vickers hardness that change with service life. Characteristic parameters such as coating erosion amount, erosion rate, and coating residual thickness are calculated. At the same time, the calculated coating residual thickness data is coupled with a thin-wall thermal resistance model that supports dynamic thickness. The thin-wall thermal resistance model is successfully applied to the setting of the material anisotropic thermal conductivity by means of equivalent coating thermal conductivity, achieving real-time prediction of coating and blade surface temperature and cooling efficiency.

[0082] (3) The present invention simplifies the calculation process and saves computing resources by avoiding the establishment of a real coating physical model. At the same time, the concept of amplification factor is introduced, which greatly reduces the time required for numerical simulation and achieves fast and accurate coating temperature prediction.

[0083] (4) The method of the present invention can simulate the coating erosion characteristics under any working conditions and has the advantages of fast simulation speed, economy, time and labor saving, and simplification of complex problems. BRIEF DESCRIPTION OF THE DRAWINGS

[0084] Figure 1 Flowchart of the method for predicting coating erosion characteristics according to the present invention;

[0085] Figure 2 This is a flow field division result diagram in an embodiment of the present invention;

[0086] Figure 3 The particle impact grid of the present invention and the distribution diagram of the adjacent grids in the direction of its velocity component;

[0087] Figure 4The coating thickness distribution diagram of the particle impact grid and the adjacent grid in the direction of its velocity component in the present invention;

[0088] Figure 5 Schematic diagram of the corrected plane Plane' and the corrected normal vector of the present invention;

[0089] Figure 6 Schematic diagram of two-dimensional wall normal vector correction according to the present invention;

[0090] Figure 7 This is a diagram showing the correction angle calculation method of the present invention;

[0091] Figure 8 This is a schematic diagram of a smooth plane collision of the present invention;

[0092] Figure 9 Schematic diagram of the judgment of the shadow area by the adjacent grid slope correction method of the present invention;

[0093] Figure 10 This is a diagram showing the verification results of the adjacent grid slope correction method of the present invention;

[0094] Figure 11 This is an erosion rate diagram before blade trajectory correction in Example 1 of the present invention;

[0095] Figure 12 This is an erosion rate diagram after blade trajectory correction in Example 1 of the present invention;

[0096] Figure 13 Figure 1 of the residual thickness of the coating on the blade surface according to the first embodiment of the present invention;

[0097] Figure 14 This is a temperature distribution diagram of the coating surface in Example 1 of the present invention;

[0098] Figure 15 This is the distribution diagram of the autocorrelation value of the coating thickness information in Example 2;

[0099] Figure 16 This is a simulated coating thickness distribution diagram derived from Example 2;

[0100] Figure 17 This is the residual thickness distribution diagram of the coating on the flat plate surface in Example 2;

[0101] Figure 18 This is the temperature distribution diagram of the coating on the flat plate surface in Example 2. DETAILED DESCRIPTION

[0102] The present invention will be further explained below with reference to the accompanying drawings.

[0103] like Figure 1 As shown, the present invention provides a method for predicting the erosion characteristics and surface temperature of a rough coating based on an adjacent grid, comprising:

[0104] Step 1: Use an electron microscope to scan the laboratory coating sample to obtain coating thickness point cloud data at a specific sampling rate and simultaneously obtain coating roughness.

[0105] Step 2: Establish a mesh model: Determine the mesh size based on the shape and size of the erosion surface and the electron microscope sampling rate. Use ICEM software to establish a physical model of the erosion surface and divide it into hexahedral meshes.

[0106] Step 3: Use the erosion surface grid coordinates to match the position information contained in the coating thickness point cloud data, and use UDM to store the coating thickness information for later use.

[0107] Step 4: Establish a particle trajectory correction method based on the slope of adjacent grids: When a particle collides with a grid on the erosion surface, search for the four adjacent grids of the hexahedral grid unit. Simultaneously, define the two adjacent grids in the same direction as the particle flow as slope correction units, and the two adjacent grids in the opposite direction as "shadow area" judgment units. Use the coating thickness information stored in the slope correction unit and the coating thickness at the impact point to recalculate the coating normal at the impact point. Use the coating thickness information of the "shadow area" judgment unit and the coating thickness at the impact point to determine whether the particle trajectory exists.

[0108] Step 5: Add a particle erosion model suitable for thermal barrier coating erosion prediction and a dynamic thin-wall thermal resistance model suitable for coating thickness changes.

[0109] Step 6: Import the above mesh model into Ansys-Fluent software, set the numerical simulation boundary conditions, load the particle trajectory correction method, particle erosion model and dynamic thin-wall thermal resistance model, set the iterative hook function, and input the amplification factor into the hook function.

[0110] Step 7: Perform steady-state calculations on the flow field. After the flow field stabilizes and is filled with particles, activate the above numerical model. The coating erosion amount is calculated using the particle erosion model. The coating residual thickness and erosion rate are calculated based on the grid area at the particle erosion site.

[0111] Step 8: Input the residual thickness of the coating in step 7 into the thin-wall thermal resistance model that supports the dynamic physical properties of the coating to obtain the surface temperature distribution of the coating;

[0112] Step 9: Determine the current surface temperature of the coating by obtaining the coating surface temperature distribution in step 8, and then calculate the coating insulation temperature or the cooling efficiency of the thermal barrier coating coating area based on the current surface temperature of the coating.

[0113] Among them, coating erosion characteristics refer to the changes in coating morphology and temperature distribution caused by coating quality loss due to erosion wear, including coating erosion rate, coating erosion mass, coating residual thickness, and blade surface and coating surface temperature.

[0114] After completing the above steps, the particle trajectory correction method based on the slope of the adjacent grid will be numerically verified and compared with the erosion rate to verify the accuracy of its model. The verified numerical prediction method will be put into practical application to realize the erosion situation and surface temperature distribution of the turbine blade surface coating at any time.

[0115] In step 1, the roughness characteristics of the target coating need to be collected and the electron microscope sampling rate needs to be determined.

[0116] In step 2, the erosion surface mesh size must match the sampling rate of the coating's roughness features. If high reproduction accuracy is not required, the point cloud data can be appropriately thinned. The mesh must also meet the following requirements: a minimum mesh quality greater than 0.25 is required, meaning that the Jacobian value of each hexahedron in the hexahedral mesh is greater than 0.25.

[0117] In step 3, the method for determining the position matching of the coating thickness point cloud data is:

[0118] ;

[0119] in, is the position coordinate of the model, is the location information contained in the point cloud data, is the device sampling rate, expressed in μm. When the discriminant is true, the coating thickness information is stored in the local grid's UDM. Furthermore, if the model mesh used for the calculation does not require high accuracy, the mesh size may be larger than the sampling rate; in this case, the discriminant still applies. It should be noted that this prediction method is still feasible when scanned coating thickness information is missing; the coating thickness information can be unified and stored in the UDM.

[0120] On the other hand, in actual use, the model area requiring simulation prediction may be larger than the sample. In this case, it is necessary to use the coating's roughness characteristics and the autocorrelation function to derive coating thickness distribution information suitable for the simulation model. The sample coating roughness data can be approximated as a random signal. Mathematically, random signals can be expressed as random processes. Specifically, changes in coating surface thickness generally conform to certain statistical laws, manifesting as a spatially random process that can be described by an autocorrelation function. The existing coating thickness signal is converted through a filter (Fourier transform) and then output as an autocorrelated random surface of a specified form. This process is relatively efficient and consistent with the actual coating distribution.

[0121] The autocorrelation function (Autocorrelation Fun) describes the statistics of the signal or data changing with the spatial (or temporal) position, that is, the correlation of the coating signal between any two points can be expressed by the autocorrelation function. For a one-dimensional signal, it can be defined as:

[0122] ;

[0123] in represents the average of the entire signal, is the translation amount. hour, , is the root mean square of the roughness thickness. Defined as:

[0124] ;

[0125] in, When sampling the coating sample surface, the interval is The discrete thickness signal needs to be further processed and converted into coating thickness information.

[0126] Taking the one-dimensional coating sampling line as an example, the number of discrete points is , the discrete interval is The RMS thickness of the coating sampling line is:

[0127] ;

[0128] in is the mean value of thickness fluctuation, defined as:

[0129] ;

[0130] After normalizing the autocorrelation function:

[0131] ;

[0132] When two points on the rough surface coincide, The maximum value is 1. When increases, the normalized autocorrelation function value decreases, that is, when the distance between two points on the rough surface is farther, the correlation is smaller.

[0133] For a two-dimensional coating surface, the autocorrelation function changes to:

[0134] ;

[0135] It represents the spatial correlation after the coating is offset. For a discrete system, it can be transformed into:

[0136] ;

[0137] in, is the value of the signal at that position.

[0138] The power spectrum, also known as the spectral density function, describes the intensity distribution of different frequency components in a signal. The correlation function and the power spectrum density function are Fourier transforms of each other. For a one-dimensional signal, its spectral density function can be expressed through the Fourier transform as follows:

[0139] ;

[0140] For 2D coatings: ;

[0141] The spectral density function describes the intensity of different spatial frequency components on the coating surface, explaining the frequency-domain characteristics of the coating thickness. The obtained frequency-domain spectral density is multiplied by white noise to simulate the undulations of the coating surface. Different spectral density gain logics predict different roughness states.

[0142] After reading the surface thickness information of the existing coating sample, the stored information is reconstructed into a two-dimensional surface for subsequent visualization and frequency domain analysis. The reconstructed data provides a suitable data structure for further operations such as Fourier transform, autocorrelation function calculation, and spectral density analysis. The degree of autocorrelation of the coating thickness information is then calculated using this two-dimensional data.

[0143] In step 4, the core of the particle trajectory correction method based on the slope of adjacent meshes lies in calling the built-in mesh surface sweeping macro and mesh cell face search macro in Ansys-Fluent to identify adjacent mesh surfaces of the impacting mesh surface. Since the impacting mesh surface is located at the boundary, a corresponding mesh cell, c0, exists. The sweeping macro locates six mesh surfaces on cell c0 (one of which is the impacting mesh surface). Based on the remaining five mesh surfaces, the corresponding c0_new_i and c1_new_i cells (i∈(0-4)) are searched. When the mesh surfaces on c0_new_i and c1_new_i lie on the same model surface as the impacting mesh surface, the adjacent mesh surfaces of the impacting mesh surface are considered to have been found. In principle, there are four adjacent mesh surfaces, but if the impacting mesh surface is at the model boundary, there may be fewer than four. Once the adjacent mesh surfaces that meet the requirements are found, they are divided into slope correction cells in the same direction as the particle and "shadow area" cells in the opposite direction of the particle based on their angle with the particle velocity. The coating thickness information at the impact position and the coating thickness at the adjacent grid are used to comprehensively calculate the coating unit normal vector at that position. At the same time, it is determined whether the particle path is obscured. If so, the impact caused by the particle is canceled and the path is deleted.

[0144] In step 5, the coating porosity is calculated as:

[0145] ;

[0146] ;

[0147] in, is the internal porosity, is the inter-column porosity, corresponding to, is the internal relative density, is the relative density between columns, is the relative density of the coating;

[0148] The UDF of the coating's Young's modulus is:

[0149] ;

[0150] in, is the density of the thermal insulation coating, is the coating Poisson’s ratio, E m is the Young's modulus of the omnidirectional dense coating and E m =210Gpa;

[0151] The fitting formula of the variation of Young's modulus of thermal barrier coating with temperature is as follows:

[0152] ;

[0153] in =132.27, =813.8K, =1574.48K.

[0154] When the operating temperature is greater than 800°C, the UDF of the change in Vickers hardness and service time is:

[0155] ;

[0156] ;

[0157] When the operating temperature is from room temperature to 800℃, the UDF of the change of Vickers hardness and service time is:

[0158] ;

[0159] ;

[0160] in, 、 、 refers to the in-plane fitting parameters, and 、 、 refers to the out-of-plane fitting parameters, and It is the hardness value estimated by the Vickers hardness calculation formula when working at a higher temperature.

[0161] The above-mentioned in-plane and out-of-plane Vickers hardness fitting parameters are shown in Table 1.

[0162] Table 1 Vickers hardness fitting parameters

[0163]

[0164] The specialized particle erosion model for thermal barrier coatings is:

[0165] ;

[0166] ;

[0167] in: ;

[0168] in, is the total erosion amount, is the particle mass, is the collision angle, is the collision velocity, is the critical angle, is the critical speed, is the cutting erosion coefficient, is the deformation erosion coefficient, where BB is the mass loss of the coating under high-angle erosion, AA is the cutting amount when the collision angle is less than the critical angle, and CC is the cutting amount when the collision angle is greater than the critical angle;

[0169] At particle collision angle Less than critical angle When the collision angle is Greater than or equal to the critical angle When the particle collision speed is Less than critical speed When the particles do not cause deformation and quality loss, the cutting erosion coefficient , deformation erosion coefficient Affected by the physical properties of particles and coatings;

[0170] The critical angle is related to the coating temperature:

[0171] ;

[0172] in is the coating temperature;

[0173] The critical speed is related to the physical properties of particles and coatings:

[0174] ;

[0175] in, is the yield stress of the coating; is the particle density; q p Poisson's ratio of particles; is the Young's modulus of the granular material; q m is the Poisson’s ratio of the target; The Young's modulus is the one that takes into account the porosity and Poisson's ratio of the target material; and the cutting erosion coefficient is Deformation erosion coefficient Expressed as:

[0176] ;

[0177] ;

[0178] in is the shape factor of the particle, is the out-of-plane Vickers hardness, is the in-plane Vickers hardness; and is the temperature dependence;

[0179] The temperature dependence and Described as:

[0180] ;

[0181] ;

[0182] H v is the Vickers hardness of the coating, V p is the particle velocity, Is the true density of the coating, which is different from the relative density of the coating The relationship is ,in is the density of dense YSZ thermal barrier coating, where Take 5650kg / m 3 , is the shape factor, 、 、 as well as are the temperature-dependent fitting parameters, which are shown in Table 2.

[0183] The cutting erosion coefficient Deformation erosion coefficient The Vickers hardness used in the deformation erosion coefficient is slightly different from the Vickers hardness used in the deformation erosion coefficient. , and the Vickers hardness used in the cutting erosion coefficient is the in-plane hardness However, the thermal barrier coating used in the present invention is an isotropic material, so the difference in hardness between the two is not large.

[0184] Table 2 Temperature dependence fitting parameters

[0185]

[0186] In step 5, referencing the dynamic thin-wall thermal resistance model includes the following steps: in the material settings, set the coating thermal conductivity property to change with the current coating thickness, and then derive the calculation formulas for the normal and tangential thermal conductivity of the coating based on the thermal resistance definition formula. The derived calculation formulas are written as UDFs and loaded into the biaxial anisotropic thermal conductivity model in the material thermal conductivity to obtain a thin-wall thermal resistance model that supports the dynamic physical properties of the coating.

[0187] Among them, thermal resistance The expression is:

[0188] ;

[0189] in is the length of the heat transfer path, is the thermal conductivity, is the cross-sectional area perpendicular to the direction of heat transfer;

[0190] When calculating the equivalent value of the thermal conductivity of the coating in the normal direction, the coating thickness changes while the bottom area remains unchanged, and its thermal resistance remains unchanged. Therefore, the expression of the thermal conductivity of the coating in the normal direction is derived as follows:

[0191] ;

[0192] When calculating the equivalent value of the thermal conductivity of the coating in the tangential direction, the coating thickness remains unchanged while the side area decreases, and its thermal resistance remains unchanged. Therefore, the expression of the thermal conductivity of the coating in the tangential direction is derived:

[0193] ;

[0194] in, is the normal equivalent thermal conductivity, is the tangential equivalent thermal conductivity, is the residual thickness of the coating.

[0195] In step 6, the iteration hook function specifies the service time represented by each iteration , and add the time magnification factor , as the Ansys-Fluent software is gradually iterated, the service time gradually accumulates, where the service time is described as:

[0196] ;

[0197] in, For the length of service, is the length of service, is the time magnification factor, is the number of iterations, is the service time represented by a single iteration.

[0198] In step 7, the expression for the erosion rate is:

[0199] ;

[0200] in is the erosion rate, is the amount of coating erosion caused by the current particles, is the cross-sectional area of ​​the current grid perpendicular to the heat transfer direction, is the mass of the collision particle, Length of service;

[0201] The expression of the residual thickness of the coating is: ;

[0202] in is the residual thickness of the coating, is the initial coating thickness, is the amount of coating erosion caused by the current particles, is the true density of the coating.

[0203] In step 9, the cooling efficiency of the thermal barrier coating coating area The calculation formula is:

[0204] ;

[0205] in, and are the gas inlet temperature and the cold air inlet temperature, is the outer surface temperature of the coating;

[0206] Coating insulation temperature The calculation formula is:

[0207] ;

[0208] Where, is the outer surface temperature of the coating, is the inner surface temperature of the coating.

[0209] First embodiment:

[0210] (1) Scan the sample and obtain its point cloud data.

[0211] (2) Establish a mesh model. Based on the shape and size of the coating-coated blades, establish a microscopic three-dimensional geometric model of the particle erosion flow field. Set the erosion surface mesh size according to the actual situation in the calculation flow field and the sampling rate of the coating thickness point cloud data. At the same time, set the boundary layer mesh at the fluid-solid coupling interface in the calculation domain. The specific size is calculated by Y+. Then set periodic boundaries on the upper and lower sides of the flow field channel to simulate the deposition of the blade channel. Finally, check and adjust the mesh parameters to make the overall mesh quality of the model reach above 0.25. The flow field division results are as follows: Figure 2 shown.

[0212] (3) Matching the coating data point cloud with the mesh model.

[0213] Due to the lack of point cloud data of typical blade surface coatings, the coating thickness is uniformly set to 300 μm.

[0214] (4) A particle trajectory correction method based on the slope of the adjacent grid is established. Among the four adjacent grids (non-model boundaries) found, two adjacent grids with the same direction as the particle flow direction are used as the calculation factors of the local coating unit normal vector. The distribution of the adjacent grids in the direction of the particle impact grid and its velocity component is as follows: Figure 3 As shown, the coating thickness distribution of the particle impact grid and the adjacent grid in the direction of its velocity component is as follows: Figure 4 As shown. The impact point grid and the centroid of the two adjacent grids can form a plane based on the principle of three points forming a plane. This plane completely coincides with the model plane, and the plane's normal vector is the normal vector of the model at that point. Based on this, the saved coating thickness information can be introduced, so that each of the three points grows a certain length along the plane's normal vector direction to obtain a new point. The growth length depends on the coating thickness at that position. The three newly obtained points then take their plane's normal vectors, which is also the corrected normal vector. Figure 5 Schematic diagram of the corrected plane and the corrected normal vector.

[0215] The specific calculation of the corrected collision angle is as follows: take a two-dimensional plane as an example. Figure 6 , Figure 7 As shown in the figure. Without changing the incident velocity, the direction of the wall normal vector is changed to achieve the purpose of correcting the particle exit velocity. The particle exit velocity can be obtained by calculating the wall normal vector correction and the restitution coefficient.

[0216] The consideration of the “shadow area” is similar to the above method. The principle is as follows: Figure 8 , Figure 9 As shown, taking a two-dimensional plane as an example, the flow collision behavior of particles in the flow field is as follows Figure 8As shown in the figure, the particles collide with the flow on a specific grid surface. However, after considering the rough morphology of the virtual coating, there is a situation where the inclination angle between the coating thickness on the opposite side of the original impact surface and the coating thickness of the grid surface is greater than the angle between the particle velocity and the plane. In this case, the particles should not continue to collide with the original grid surface, as shown in the figure. Figure 9 shown.

[0217] (5) Compile the UDF function and load it into Ansys-Fluent software.

[0218] (6) Set the calculation boundary conditions, open the energy and discrete term models, set the iteration hook function to link the number of iteration steps with the service time, and perform particle collision trajectory analysis in the flow field.

[0219] (7) When the computational domain reaches a steady state and the particles fill the flow field, the particle trajectory correction method, particle erosion model, and dynamic thin-wall thermal resistance model are activated. The angle between the calculated corrected normal vector and the vector Z[0,0,1] of the standard sinusoidal wave surface is used as a verification quantity. The effectiveness and accuracy of the method are verified by comparing the changing trend of the angle between the normal vector of the known coating fluctuation sinusoidal function and the vector z with the normal vector angle of the real roughness model according to the sinusoidal curve fluctuation. Figure 10 To verify the results, it is shown that the particle trajectory correction method of the present invention has a higher correction accuracy.

[0220] (8) The residual thickness of the coating is obtained and connected with the dynamic thin-wall thermal resistance model to obtain the surface temperature distribution. The erosion rate (g / kg, the coating mass loss caused by particles per unit mass) within a certain service time, the number of particle collisions, and the residual thickness of the coating are stored in the user-defined storage UDM and exported. Then, compared with the erosion results of the smooth surface, it is found that the erosion prediction method considering the surface roughness is more consistent with the actual coating erosion results, which proves the accuracy of the method of the present invention. The erosion rate distribution without particle trajectory correction is shown in Figure 11, and the results after using it are shown in Figure 12. Figure 12 As shown, it can be seen that the correction result is obvious.

[0221] (9) Example Result Analysis: This example studies the particle erosion characteristics of a coating applied on the surface of a typical turbine blade. The predicted residual thickness of the coating on the blade surface after 1000 hours of service is as follows: Figure 13 The temperature distribution is shown as Figure 14As shown, this demonstrates that the proposed prediction method is capable of simulating blade surface coating erosion over any service life. This method saves significant computing resources and time while maintaining high simulation accuracy. It can predict characteristic parameters such as coating thickness, temperature changes, blade cooling efficiency, and coating insulation temperature over the life of turbine blades. This provides effective data support for preventing blade ablation caused by coating wear, significantly impacting turbine blade protection.

[0222] Second embodiment:

[0223] (1) Obtain coating sample point cloud data.

[0224] (2) Establish a grid model.

[0225] In this example, a CMC plate coated with a thermal barrier coating is used as the prediction object, and the grid division method is the same as in Example 1. Autocorrelation analysis is performed on the coating thickness data obtained by scanning. Figure 15 The autocorrelation degree of the coating thickness information is shown. Among them, the values ​​of the x-axis and y-axis represent the offset in the x- and y-directions when calculating the autocorrelation function, that is, the distance between the two data points in this direction. The ACF value represents the similarity of the coating thickness at this distance. The larger the autocorrelation function value, the higher the similarity of the data points at this offset, and vice versa. The ACF value is generally between -1 and 1. In the coating thickness calculation, it varies according to the nature of the data. At the same time, the power spectral density of the coating thickness information is calculated as the noise source gain. The values ​​of the power spectral density in the x- and y-directions are the frequency components obtained by Fourier transforming the original signal.

[0226] Figure 16 The simulated coating thickness distribution caused by the spectral density function of single gain is given, and it is found that the arithmetic average roughness (Ra) of the sample is 0.010692mm, which is close to the arithmetic average roughness (Ra) of the sample coating of 0.011462mm.

[0227] (3) Matching the coating data point cloud with the mesh model.

[0228] (4) Compile and load the particle collision trajectory correction method.

[0229] (5) The particle erosion model and the dynamic thin-wall thermal resistance model are used in the same manner as in Example 1.

[0230] (6) Set the hook function and input the amplification factor into the hook function.

[0231] (7) Perform steady-state calculations on the computational domain. After the flow field stabilizes and is filled with particles, activate the particle trajectory correction method, particle erosion model, and dynamic thin-wall thermal resistance model. The coating erosion amount is obtained after calculation using the particle erosion model. Combined with the grid area at the particle erosion location, the coating residual thickness and erosion rate are calculated.

[0232] (8) Obtain the coating temperature distribution.

[0233] (9) Result analysis: This example studies the particle erosion characteristics of a flat plate coated with a thermal barrier coating, and predicts the residual thickness of the surface coating after 1000 hours of continuous erosion. Figure 17 The temperature distribution is shown as Figure 18 shown.

[0234] The above is only a preferred embodiment of the present invention. It should be pointed out that for ordinary technicians in this technical field, several improvements and modifications can be made without departing from the principles of the present invention. These improvements and modifications should also be regarded as within the scope of protection of the present invention.

Claims

1. A method for predicting the erosion characteristics and surface temperature of rough coatings based on adjacency grids, characterized in that: include: Step 1: Scan the laboratory coating sample using an electron microscope to obtain coating thickness point cloud data at a preset sampling rate and simultaneously obtain coating roughness; Step 2: Establish a mesh model: Determine the mesh size based on the shape and size of the erosion surface and the electron microscope sampling rate. Use ICEM software to establish a physical model of the erosion surface and divide it into hexahedral meshes. Step 3: Use the erosion surface grid coordinates in the erosion surface physical model to match the position information contained in the coating thickness point cloud data, and use UDM to store the coating thickness information for later use; Step 4: Establish a particle trajectory correction method based on the slope of adjacent grids: When a particle collides with a grid on the erosion surface, find the four adjacent grids of the collided hexahedral grid cell; The two adjacent grids in the same direction as the particle flow are defined as slope correction units, and the two adjacent grids in the opposite direction are defined as shadow area judgment units. The coating thickness information stored in the slope correction unit and the coating thickness at the impact point are used to recalculate the coating normal at the impact point. Use the shadow area to determine the coating thickness of the unit and the coating thickness at the impact point to determine whether the particle trajectory exists. The judgment method is to detect whether the coating stored in the grid unit opposite to the particle flow direction obscures the particle incident trajectory. If the obscuration indicates that the area obscured by the particle incident trajectory is a shadow area, the influence of the particle on the coating is cancelled. Step 5: Compile the particle trajectory correction method UDF established in step 4, compile the particle erosion model UDF suitable for thermal barrier coating erosion prediction, and compile the dynamic thin-wall thermal resistance model UDF suitable for coating thickness variation; Step 6: Import the mesh model into Ansys-Fluent software, set the numerical simulation boundary conditions, load the particle trajectory correction method UDF, particle erosion model UDF, and dynamic thin-wall thermal resistance model UDF, set the iterative hook function, and input the amplification factor into the hook function; Step 7: Perform steady-state calculations on the computational domain. Once the flow field is stable and filled with particles, activate the particle trajectory correction method, particle erosion model, and dynamic thin-wall thermal resistance model. The particle erosion model is used to calculate the coating erosion amount. Combined with the mesh area at the particle erosion site, the coating residual thickness and erosion rate are calculated. Step 8: Input the residual thickness of the coating in step 7 into the dynamic thin-wall thermal resistance model to obtain the surface temperature distribution of the coating; Step 9: Determine the current surface temperature of the coating by obtaining the coating surface temperature distribution through step 8, and then calculate the coating insulation temperature or the cooling efficiency of the thermal barrier coating coating area based on the current surface temperature of the coating; In step 5, the particle erosion model suitable for thermal barrier coating erosion prediction is: ； ； in: ; in, is the total erosion amount, is the particle mass, is the collision angle, is the collision velocity, is the critical angle, is the critical speed, is the cutting erosion coefficient, is the deformation erosion coefficient, where BB is the mass loss of the coating under high-angle erosion, AA is the cutting amount when the collision angle is less than the critical angle, and CC is the cutting amount when the collision angle is greater than the critical angle; At particle collision angle Less than critical angle When the collision angle is Greater than or equal to the critical angle When the particle collision speed is Less than critical speed When the particles do not cause deformation and quality loss, the cutting erosion coefficient , deformation erosion coefficient Affected by the physical properties of particles and coatings; The critical angle is related to the coating temperature: ； in is the coating temperature; The critical speed is related to the physical properties of particles and coatings: ； in, is the yield stress of the coating; is the particle density; q p Poisson's ratio of particles; is the Young's modulus of the granular material; q m is the Poisson’s ratio of the target; The Young's modulus is the one that takes into account the porosity and Poisson's ratio of the target material; and the cutting erosion coefficient is Deformation erosion coefficient Expressed as: ； ； in is the shape factor of the particle, is the out-of-plane Vickers hardness, is the in-plane Vickers hardness; and is the temperature dependence; The temperature dependence and Described as: ； ； H v is the Vickers hardness of the coating, V p is the particle velocity, Is the true density of the coating, which is different from the relative density of the coating The relationship is ,in is the density of dense YSZ thermal barrier coating, is the shape factor, 、 、 as well as Fit parameters for the temperature dependence.

2. The method for predicting erosion characteristics and surface temperature of rough coatings based on adjacent grids according to claim 1, characterized in that: In step 2, the minimum mesh quality is higher than 0.

25.

3. The method for predicting erosion characteristics and surface temperature of rough coatings based on adjacent grids according to claim 1, characterized in that: In step 3, the method for determining the position matching of the coating thickness point cloud data is: ； in, is the position coordinate of the model, is the location information contained in the point cloud data, is the device sampling rate in μm; when the discriminant is true, the coating thickness information is stored in the UDM of the local grid.

4. The method for predicting erosion characteristics and surface temperature of rough coatings based on adjacency grids according to claim 1, characterized in that: In step 4, the particle trajectory correction method based on the slope of the adjacent grid uses the sweep macro built into Ansys-Fluent to find the four adjacent grid surfaces of the collision grid surface. The search is as follows: the grid surface currently colliding is on the boundary surface, and there is a corresponding grid unit, namely the c0 unit; call the sweep macro to find the six grid surfaces on the c0 unit, one of which is the grid surface where the collision occurs, and find the corresponding c0_new_i and c1_new_i units based on the remaining five grid surfaces, i∈(0-4); when the grid surfaces on c0_new_i and c1_new_i are on the same model surface as the grid surface where the collision occurs, it is considered The adjacent mesh surfaces of the impact mesh surface are found; in principle, there are four adjacent mesh surfaces, but when the impact mesh surface is at the model boundary, there will be less than four adjacent mesh surfaces; after finding the adjacent mesh surface that meets the requirements, it is divided into slope correction units in the same direction as the particle and shadow area judgment units in the opposite direction of the particle according to the angle with the particle velocity; the rough coating plane is refitted using the coating thickness information at the impact position and the coating thickness at the adjacent mesh, and the plane normal vector is calculated, and the particle collision trajectory is corrected based on the normal vector; when the coating in the opposite direction of the particle flow obscures the incident trajectory of the particle, the influence of the particle is cancelled and the incident trajectory of the particle is deleted.

5. The method for predicting erosion characteristics and surface temperature of rough coatings based on adjacency grids according to claim 1, characterized in that: In step 5, the dynamic thin-wall thermal resistance model reduces the thermal barrier coating with different thickness and consistent thermal conductivity to an equivalent coating with the thickness of the initial coating and a changing thermal conductivity, and stipulates that the equivalent coating thickness is the initial coating thickness. ; The calculation formulas for the normal and tangential thermal conductivity of the coating are derived based on the definition formula of thermal resistance; Among them, thermal resistance The expression is: ； in is the length of the heat transfer path, is the thermal conductivity, is the cross-sectional area perpendicular to the direction of heat transfer; When calculating the equivalent value of the thermal conductivity of the coating in the normal direction, the coating thickness changes while the bottom area remains unchanged, and its thermal resistance remains unchanged. Therefore, the expression of the thermal conductivity of the coating in the normal direction is derived as follows: ； When calculating the equivalent value of the thermal conductivity of the coating in the tangential direction, the coating thickness remains unchanged while the side area decreases, and its thermal resistance remains unchanged. Therefore, the expression of the thermal conductivity of the coating in the tangential direction is derived: ； in, is the normal equivalent thermal conductivity, is the tangential equivalent thermal conductivity, is the residual thickness of the coating.

6. The method for predicting erosion characteristics and surface temperature of rough coatings based on adjacency grids according to claim 1, characterized in that: In step 6, the iteration hook function specifies the service time represented by each iteration , and add the time magnification factor , as the Ansys-Fluent software is gradually iterated, the service time gradually accumulates, where the service time is described as: ； in, For the length of service, is the length of service, is the time magnification factor, is the number of iterations, is the service time represented by a single iteration.

7. The method for predicting erosion characteristics and surface temperature of rough coatings based on adjacency grids according to claim 1, characterized in that: In step 7, the expression for the erosion rate is: ； in is the erosion rate, is the amount of coating erosion caused by the current particles, is the cross-sectional area of ​​the current grid perpendicular to the heat transfer direction, is the mass of the collision particle, Length of service; The expression of the residual thickness of the coating is: ; in is the residual thickness of the coating, is the initial coating thickness, is the amount of coating erosion caused by the current particles, is the true density of the coating.

8. The method for predicting erosion characteristics and surface temperature of rough coatings based on adjacency grids according to claim 1, characterized in that: In step 9, the cooling efficiency of the thermal barrier coating coating area The calculation formula is: ； in, and are the gas inlet temperature and the cold air inlet temperature, is the outer surface temperature of the coating; Coating insulation temperature The calculation formula is: ； Where, is the outer surface temperature of the coating, is the inner surface temperature of the coating.

Images