Dynamic Coupling Analysis on Thermo–Chemo–Mechanical Field and Fluid–Structure Interaction for Aero-Engine Turbine Blade with Functional Gradient Thermal Barrier Coatings

: In this study, an extended layerwise/solid-element (XLW/SE) method is developed for the thermo–chemo–mechanical (TCM) coupling problem of an aero-engine turbine blade with thermal barrier coatings (TBCs). The method consists of two parts, the extended layerwise (XLW) method and the three-dimensional (3D) solid-element (SE) method, which are adopted to formulate the governing equations of TBCs and substrate, respectively. Then, according to the compatibility conditions of displacement, temperature, concentration and internal force equilibrium at the TBCs/substrate interface, the governing equation of the ﬁnal blade structure is assembled. Through a time integration, the dynamic responses of displacement, temperature and concentration can be calculated. In addition, the ﬂuid–structure coupling analysis is conducted by using COMSOL. The nonuniform thermal load is imported into the XLW/SE method to calculate the mechanical response of blade structure. Finally, the corresponding computing program is compiled with C++. In numerical examples, the TCM coupling analysis is conducted on the blade structure with and without interfacial debonding and delamination damages. To validate the effectiveness of the proposed method, the dynamic TCM responses of the XLW/SE model is compared with those of a 3D elastic model generated by COMSOL, which shows that the two models are in good agreement.


Introduction
With the increasing output power of gas turbine and the increasing weight ratio of aircraft engine, the turbine inlet temperature is greatly improved.How to ensure the safe operation of aero-engine turbine blades under a high-temperature, high-pressure, oxidation and corrosion environment has long been the major focus in the field of aircraft engineering.In response to this problem, a variety of technologies have been adopted, for example, the blade casting technology [1] for improving the blade temperature resistance, the film cooling technology [2-7] with real-time cooling and thermal protection system and the TBCs technology with low thermal conductivity and high stability [8,9].
As a multilayer structure, TBCs can effectively reduce the substrate temperature, and maintain the machinability and reliability of the alloy substrate.The cross-sectional view of TBCs structure is shown in Figure 1 [10].The substrate is composed of nickel-based or cobalt-based superalloys with sufficient mechanical strength.The TBCs consist of two parts: a ceramic coat and a bond coat.The ceramic coat is composed of Y 2 O 3 and ZrO 2 (YSZ).Its outstanding thermal insulation performance can protect the substrate from oxidation, corrosion and erosion.The bond coat is an antioxidation layer composed of NiCrAlY or NiCoCrAlY alloy.It is used to firmly adhere the ceramic coat to the substrate, so as to reduce the mismatch stress caused by the different thermal expansion coefficients of the ceramic coat and substrate.The aero-engine turbine blade with TBCs is commonly used in a multiphysics environment.The variations of service environment may affect the life cycle of the turbine blade.In recent years, several studies have been devoted to the multi-physics coupling analysis of TBCs structure.For instance, Liu et al. [11] simulated the thermal stress of a turbine blade with TBCs under the conditions of air film cooling and thermally grown oxide (TGO).Wang et al. [12] simulated the thermal insulation of a turbine blade with TBCs in a cascade flow field.Zhou and Hashida [13] analyzed the thermal stress fields in TBCs, considering the coupled effects of temperature gradient and oxidation.Zhang et al. [14] studied the effect of Cr and Co oxides on TGO growth at 1000 • C. The thermodynamics and kinetics of oxidation were analyzed.
Some studies have focused on the effects of various damage forms on TBCs structure [15][16][17][18][19][20].For instance, Xiao et al. [21] developed a phase-field damage model to quantify the cumulative damage, crack nucleation, propagation and coalescence, and they also analyzed the impact of thermal cycling and growth stresses on the TBCs structure.Li and Shan [22] performed a TCM coupling analysis of TBCs structure with multiple interface and transverse cracks; however, they ignored the impacts of the metal substrate on the TBCs structure.
The multiphysics coupling problem of TBCs fracture is still the bottleneck of turbine blade design, optimization and production.The XLW method can be used as an efficient tool to address this problem.For instance, Li et al. [23][24][25][26][27][28][29] analyzed the layered damage and transverse cracks of composite laminates by using the XLW method, which achieved a higher accuracy and less computation compared with the 3D elastic model.In the finite-element governing equation, the degrees of freedom (DOFs) for the upper and lower surfaces of the shell-like structures can be adjusted.This facilitates the combination of the XLW method with other methods according to the compatibility conditions of TBCs/substrate interface [30][31][32][33][34][35][36][37].Li [38,39] developed a steady-state and dynamic thermomechanical XLW method for the laminated composite plates with multiple delaminations and transverse cracks.A similar study by Li and Shan [40] developed a chemo-mechanical XLW method for the multilayered porous media with different damage forms.
Based on the above review, most studies focused on the fracture analysis of coatings only.Some studies considered both the TBC and the substrate but simplified the structure as a 2D problem and did not approach the complex curved surface structure of the actual blade.Furthermore, the geometric dimensions of the coatings and blade differ greatly, which is a cross-scale problem.Hence, using traditional solid element to construct the grids increases the grid number and computation, which also increases the complexity of the finite element analysis.In addition, the real temperature field, which has to be determined by a thermal fluid-structure coupling analysis, was not subject to the existing fracture analysis model.Therefore, according to the existing research, this study proposes a computational framework for an aero-engine turbine with TBCs and develops the corresponding computing program.This method can simultaneously apply two methods for cross-scale modeling, and without any additional geometric assumptions.Moreover, the method can be used for a thermo-chemo-mechanical and fluid-structure coupling analysis of complex structures with damage, while ensuring the accuracy of the analysis.This is difficult for other methods to achieve, and the proposed damage analysis model can be used for a complex structure, rather than a simple plate and circular structures.
The rest of the paper is organized as follows.Section 2 introduces the TCM variational theorem.The TCM formulations for the substrate are defined by using SE method.In Section 3, the functionally gradient material (FGM) structure of TBCs is described.The coatings with delamination and transverse cracks is modeled by using the XLW method.In Section 4, the XLW/SE method for the TCM coupling analysis is presented.The governing equation of a turbine blade with TBCs is assembled.Then, a mixed semidiscretization method is applied to calculate thermal, chemical and mechanical responses.In Section 5, a fluid-structure coupling analysis is conducted to solve the flow field of the turbine blade.In Section 6, the flowchart of the computing program for the TCM and fluid-structure coupling analyses is presented.In Section 7, the numerical examples are given to verify the convergence and accuracy of the proposed XLW/SE method.Finally, a conclusion and future work are given in Section 8.

Thermo-Chemo-Mechanical Variational Theorem
According to the Hamilton principle, the energy forms of mechanical, chemical and thermal fields are considered in this study.The TCM mixed variational principle for the present solution domain is given by The boundary conditions are expressed as In Equation ( 1), the first row represents the energy generated by the mechanical field.ε αβ and σ αβ denote the mechanical strain and mechanical stress, respectively.u α (α = 1, 2, 3) denotes the displacement component in the x, y, z directions, respectively.ūα denotes the imposed displacement on boundary Γ u .f α and F α denote the mechanical body force and mechanical surface force on boundary Γ F , respectively.ρ denotes the mass density.
The second row represents the energy generated by the thermal field.θ, e α and ϑ 0 denote temperature change, temperature gradient and positive reference temperature, respectively.h α denotes the heat flux vector, and h denotes the heat flux on boundary Γ H . s T , η and θ denote the heat source, entropy density change rate and imposed temperature change on boundary Γ θ , respectively.
The third row represents the energy generated by the chemical field.c denotes the chemical ion concentration, and c denotes the imposed concentration on boundary Γ c .m α denotes the mass flux vector, and m denotes the mass flux on boundary Γ c .μ, τ and g α denote the chemical potential change rate, penalty factor and chemical ion concentration gradient, respectively.n α and n β denote the outward unit normal vectors.

Solid-Element Method for Thermo-Chemo-Mechanical Coupling Problem of Substrate Structure
The SE method is adopted to solve the 3D dynamic TCM coupling problems of the substrate structure.The displacement u α , temperature change θ and chemical ion concentration c can be expressed as where n e denotes the number of nodes per element.As the isoparametric elements were employed, n e = 8 in this study.N i denotes the shape function.The gradient relations of strain ε αβ , temperature gradient e α and chemical gradient components g α corresponding to u α , θ and c are, respectively, defined as Based on the linear constitutive law for the TCM problem, the constitutive equations of mechanical stress σ αβ , entropy density η, chemical potential µ, heat flux h α and mass flux m α are given by where λ αβ and denote the thermal expansion coefficient and thermal stress-temperature coefficient, respectively.= ρc v /ϑ 0 , where c v denotes the specific heat capacity.C αβkl , ν TC and ν CT denote the stiffness, thermo-chemical coefficient and chemo-thermal coefficient, respectively.k, κ αβ and R αβ denote the chemical potential constant, heat conduction coefficient and mechanical-chemical coefficient, respectively.ϕ αβ denotes the diffusion coefficient of anions and cations.Based on Equations ( 1) and (3), the finite formulations of the TCM coupling problem can be derived as where M, C and K denote the mass, damping and stiffness matrices, respectively.u and ü denote the velocity and acceleration vectors, respectively.θ and ċ denote the temperature and chemical concentration change rates, respectively.F u , F θ and F c denote the mechanical, thermal and chemical load vectors, respectively.The submatrix of each aforementioned element is expressed as where B denotes the strain-displacement matrix per finite element, which is expressed as ) 3. Thermo-Chemo-Mechanical Extended Layerwise Method for TBCs Structure

Functionally Gradient Material Structure
In TBCs, the bond coat is used as a transition layer to reduce the difference of thermal expansion coefficient between the metal substrate and the ceramic coat.Nevertheless, such structure yields a large thermal stress and residual stress during service, thus greatly shortening its life cycle.A TBCs structure with the gradient composition ratio, as shown in Figure 2, can eliminate the macroscopic interface and reduce the residual stress and thermal stress.For the FGM structure, the volume fractions of ceramic coat V f 1 and bond coat V f 2 are expressed as where n denotes the gradient index.h denotes the thickness of the FGM structure.z denotes the coordinate in the layer thickness direction, ranging from − h 2 to h 2 .In Figure 2, the top surface is pure ceramic coating (i.e., z = h 2 , V f 1 = 1, V f 2 = 0), while the bottom surface is pure bond coating (i.e., z = − h 2 , V f 1 = 0, V f 2 = 1).The middle part is the functionally graded structure composed of ceramic and bond materials.The material properties of the FGM structure change in gradient along the layer thickness direction.Based on Equation (12), the material property is expressed as where P(z) denotes the material property parameter and is the function of z.P 1 and P 2 are the material property parameters of the top and bottom surfaces, respectively.

Thermo-Chemo-Mechanical Extended Layerwise Method of TBCs Structure
The XLW method transforms the 3D fracture problem into a one-dimensional fracture problem in the z direction and a two-dimensional fracture problem in the x and y directions.The linear Lagrange interpolation function is used to model the dynamic displacement, temperature change and chemical ion concentration of TBCs with multiple delaminations.As shown in Figure 3 [29], in the z direction, the interpolation nodes are located at the bottom and top surfaces of TBCs and the middle surface of each mathematical layer.The displacement u α , temperature change θ and concentration c at point (x, y, z) can be discretized as where φ k is the linear Lagrange interpolation function.Ξ k is the shape function for modeling displacement discontinuities caused by delaminations.Θ k is the shape function for modeling strain discontinuities caused by interlaminar interfaces.k is the layer number, and N is the total number of mathematical layers.N D is the number of nodes.u αik , u αlk and u αrk are the standard nodal DOF, the additional nodal DOF for displacement discontinuities and the additional nodal DOF for strain discontinuities, respectively.i, l and r are the corresponding subscripts, and α = 1, 2, 3 are the x, y and z directions, respectively.(c ik , c lk , c rk ) and (θ ik , θ lk , θ rk ) are the nodal DOFs of chemical and temperature fields, respectively.
For the TBCs with transverse cracks in the x and y directions, the in-plane nodal displacement u αζk , temperature change θ ζk and concentration c ζk can be discretized as where ψ m denotes the two-dimensional Lagrange interpolation function.Λ s denotes the shape function for modeling the in-plane transverse cracks.Π hb denotes the shape function for modeling the in-plane transverse crack tips.Ũαζkm , Θζkm and Cζkm are the standard nodal DOFs.Ūαζks , Θζks and Cζks are the additional nodal DOFs yielded from the in-plane transverse cracks.Ûαζkhb , Θζkhb and Ĉζkhb are the additional nodal DOFs yielded from the in-plane transverse crack tips.ζ = i, l, r.Based on Equations ( 1), ( 14) and ( 15), the finite-element governing equation for the dynamic TCM-XLW method is expressed as The stiffness matrix K, mass matrix M and damping matrix C are given by Note that p = m, s, f and q = n, g, h, where are the numbers of in-plane nodes enriched by the transverse crack discontinuities and the transverse crack tips, respectively.b = 1, 2, • • • N f , and N f denotes the number of the crack-tip enrichment functions.The values of index pairs (m, n), (s, g) and ( f , h) are related to the shape functions ψ, Λ and Π, respectively.In Equation ( 16), K uu αβζηkepq can be obtained by replacing the shape function ψ corresponding to p and q.Other submatrices of the element stiffness, mass and damping matrices can be obtained in a similar way.The detailed expressions of the above submatrices are listed in Appendix A.

Transformation of Coordinate System
The element stiffness matrix in the local coordinate system was derived in the previous subsection.For the complexity of the turbine blade surface, the coordinate systems needs to be modified.Figure 4 displays the global coordinate system and local coordinate system of finite elements.In this study, x , y and z are the three directions of the global coordinate system.x, y and z denote the three directions of the local coordinate system.Note that x and x are parallel to each other.In two coordinate systems, the nodal displacement vectors are expressed as The transformation equations of two coordinate systems and the transformation matrix L are given by where λ is the cosine of the local nodal coordinates in x, y and z directions of the global coordinate system.For instance, in Figure 4, the local nodal coordinates of Node 1 and Node 4 are (x 1 , y 1 , z 1 ) and (x 4 , y 4 , z 4 ).The global nodal coordinates are (x 1 , y 1 , z 1 ) and (x 4 , y 4 , z 4 ).The direction cosines λ x , λ y and λ z are, respectively, expressed as and where s 14 is the distance between two nodes.By following the above method, the element stiffness matrix K e and load vector Q e in the global coordinate system can be transformed from the element stiffness matrix K e and load vector Q e in the local coordinate system.The transformation equations are given by where the transformation matrix T is expressed as Based on the above, the governing equation of the overall structure can be obtained.According to the nodal displacement vectors a i and a i in the two coordinate systems, the element stress σ ij can be calculated, which satisfies the following functional relationship and the projection of the local coordinate system on the global coordinate system is where e i , e i denotes the angle between two basis vectors e i and e i .

Modeling of Aero-Engine Turbine Blade with TBCs
The aero-engine turbine blade is a complex structure, as shown in Figure 5.The curvature and thickness of the blade are not fixed values at different positions.Moreover, the thickness difference between the TBCs and the substrate may increase the complexity of the finite-element modeling.To satisfy the actual condition, the proposed XLW/SE method was adopted to solve the dynamic TCM coupling problem of the turbine blade with TBCs.The thickness of the FGM structure was 150 µm.The cooling air-film holes of the blade were also considered.

Governing Equations of Aero-Engine Turbine Blade with TBCs
In this section, the governing equations of the TBCs and the substrate for the turbine blade are formulated.In order to ensure the compatibility of nodes and elements on the TBCs/substrate interface, the compatibility conditions of displacement, temperature and concentration at the interface are defined, as shown in Figure 6.The meshing of the TBCs is based on the outward free surface of the substrate.Finally, the governing equation of the overall blade structure is assembled.
Compatibility conditions of TCM coupling problem for the aero-engine turbine blade with TBCs.
With respect to the TBCs structure, the displacement, temperature and concentration DOFs can be divided into two parts: (i) the external DOFs of the TBCs/substrate interface and (ii) the internal DOFS of the TBCs structure excluding the interface.Referring to the governing equation for the TCM-XLW method in Equation ( 16), the governing equation of TBCs is given by where P represents TBCs.The subscripts 1 and 2 are the external and internal DOFs, respectively.The subscripts 12 and 21 denote the coupling of external and internal DOFs.θu and cu denote the mechanical-thermal coupling and mechanical-chemical coupling.K, C and M denote the stiffness, damping and mass matrices, respectively.u With respect to the substrate structure, the displacements, temperature and concentration DOFs can also be divided into two parts: (i) the external DOFs of the TBCs/substrate interface and (ii) the internal DOFs of the substrate excluding the interface.The governing equation of the substrate is given by where s represents the substrate.The compatibility conditions of displacement, temperature, concentration and interface force equilibrium at the TBCs/substrate interface are given by In Equations ( 28) and ( 29), the summation of the first rows is calculated as Based on the above, Equations ( 31) to (33), and the second, fourth and sixth rows of Equations ( 28) and ( 29), the final governing equation of the overall structure (i.e., substrate and TBCs) is assembled as where the submatrices of stiffness, damping, mass and the vector matrices of displacement, acceleration and load are as follows

Time Integration of Thermo-Chemo-Mechanical Responses
Referring to the Newmark and Crank-Nicolson integration method, a mixed semidiscretization method was developed to calculate the dynamic TCM responses of blade structure at each time step ∆t.The responses including temperature change Θ t+∆t , concentration C t+∆t , velocity Ut+∆t and displacement U t+∆t at time t + ∆t are given by where γ, β 1 and β 2 are the coefficients that determine the numerical stability and accuracy.In this study, the value of β 1 was set to 0.5 to ensure the unconditional stability without numerical dissipation.The values of γ and β 2 were set to 0.5 and 0.25 to ensure the numerical stability.The derivatives of the aforementioned responses are expressed as where . Based on the above, Equation (34) at time t + ∆t can be rewritten as Substituting Equation ( 37) into (38), we have where

Fluid-Structure Analysis for Aero-Engine Turbine Blade with Cooling Structure
Gas turbine blades are continuously impacted by high-temperature and high-pressure gases in the actual working condition.An efficient cooling method is needed to reduce the blade temperature and ensure its strength and service life.The film cooling technology is commonly used in turbine blade structure, which has multiple cooling holes and an internal cooling chamber, as shown in Figure 7 [41].In working state, the cooling gas flows through the internal passage to cool the blades from inside.Then, the cooling gas flows through the cooling holes and mixes with the external high temperature gases to reduce the wall temperature of turbine blades.The temperature change may result in a large thermal stress in the solid interior of turbine blades.This is the main reason for blade creep, fatigue and other forms of damage, which directly affect the safety and stability of aero-engine operation.Therefore, it is necessary to obtain an accurate thermal load distribution of turbine blades with the cooling structure and TBCs.
In this study, the Reynolds-averaged Naiver-Stokes (RANS) equations were used to solve the flow field.Since the major form of gas flow in an aero-engine turbine is turbulence, the shear stress transfer (SST) model was used as the turbulence model to complete the RANS equations and simulate the motion state of turbulence.The SST model was developed by combining the k − ω and k − ε models, which had fewer grids and better adaptability.The differential form of the governing equations for an unsteady fluid was expressed as where S M , S E and h * are the momentum source term, energy source term and total specific enthalpy, respectively.v α , p f and ρ f are the fluid velocity, pressure and density, respectively.µ f is the dynamic viscosity.

Computational Framework of Extended-Layerwise/Solid-Element Method
To implement the proposed computational framework, the corresponding computing program was compiled by using C++. Figure 8 presents the flowchart of the computational process, which can be divided into four steps.Details of each step are described as follows.(1) The models of TBCs and substrate are input into the XLW/SE method to start the program.For each model, the boundary conditions of temperature, concentration and displacement are defined.For the flow-thermal coupling problem, the SST model based on the RANS method is used as the turbulence model.By using the fluid dynamics solver provided by COMSOL, the nonuniform thermal load is obtained and then input into the XLW/SE method to calculate the mechanical response of blade (2) By following the SE method for the substrate described in Section 2 and the XLW method for TBCs from Section 3, the stiffness, mass and damping matrices of the substrate and TBCs are formulated, and then used to formulate their governing equations, respectively.
(3) By following the XLW/SE method from Section 4, the displacement, temperature and concentration DOFs are determined.The governing equation of the overall structure is assembled.
(4) The iteration of the time integration is performed.Referring to Equation ( 40), the equivalent stiffness matrix K and load vector Ft+∆t are calculated.Then, referring to Equation (39), the responses of temperature, concentration and displacement at time t + ∆t are calculated.The iteration is stopped when the stopping criterion is satisfied.

Validation of Extended-Layerwise/Solid-Element Method
To verify the convergence and accuracy of the proposed XLW/SE method, an aeroengine turbine blade with TBCs was used in the numerical examples, which was subjected to thermal and chemical loads.Each layer of the model was isotropic, and the boundary condition was bottom clamped.The outermost wall temperature was 800 K, and the chemical icon concentration was 80 mol/m 3 .The initial temperature and concentration were 0 K and 0 mol/m 3 , respectively.The TBCs and substrate material parameters that characterized the chemical, thermal and mechanical properties are listed in Table 1.The functional gradient index n was 0.5.To begin, the convergence of the meshing scheme and time step should be analyzed.In this numerical example, the time step ∆t = 1 s and three different meshing schemes were adopted, as shown in Figure 9.The numbers of elements and nodes for each meshing scheme are listed in Table 2.A random point A (16.0 cm, −12.2 cm, 0.28 cm) located inside the blade structure was selected.Its dynamic responses of displacement u 1 , u 2 and u 3 in the x, y and z directions are, respectively, shown in Figure 10.It can be observed that the displacement curves gradually converge with time, and the convergence trends of the three meshing schemes are consistent.Since the curves of the second and third meshing schemes coincide, it can be confirmed that the second meshing scheme is convergent.For the second meshing scheme, the displacement curves of three different time steps ∆t are presented Figure 11 shows.It can be observed that when ∆t = 1 s, the displacement curve shows a convergent trend and can provide the convergence result.Therefore, the second meshing scheme with the time step ∆t = 1 s is employed in the following cases; see the numerical examples below.To validate the effectiveness of the proposed XLW/SE method, a 3D elastic model developed in COMSOL 5.6 was compared with the XLW/SE model.Firstly, by following the SE method from Section 2 and the XLW method from Section 3, the 3D models of the substrate and TBCs were simulated, respectively.Then, by following the XLW/SE method from Section 4, the responses of displacement, temperature and concentration of point A were obtained.Comparing the response curves and maximum displacement of the XLW/ME model with those of the 3D elastic model, as shown in Figure 12 and Table 3.The corresponding nephograms at 250 s and 1000 s are presented in Figure 13.According to the peak value and steady-state value in Figure 12 and the percent errors in Table 3, it is found that the TCM responses of the TCM coupling model is in good agreement with the 3D elastic model, which verifies the effectiveness of the proposed XLW/SE method.

TBCs with Debonding at Interface
The spalling of TBC during service results in the generation of a debonding area and even the failure of the TBC system.This section simulates the effect of TBC debonding on TCM responses (i.e., displacement, temperature and concentration distribution) in the layer thickness direction.The working conditions of the blade structure were the same as those in the previous section.Since the material of TBCs was functionally graded, the top layer and bottom layer were pure ZrO 2 and NiCoCrAlY, respectively.Two cases were considered: (i) Case 1: aero-engine turbine blade without debonding at the interface; (ii) Case 2: aeroengine turbine blade with debonding at the interface.The interface between the coating and substrate can exhibit debonding damage, as shown in Figure 14.
A random point B (12.0 cm, 7.3 cm and 3.8 cm) on the substrate behind the damaged area was selected, and its dynamic responses of displacement, temperature and concentration are shown in Figure 15.The results show that the presence of a debonding crack affects the normal conduction of temperature and concentration and prolongs the time for the temperature and concentration to reach stability.Meanwhile, the existence of the debonding has a significant impact on the displacement changes in all the directions of point B. In Figure 15, the displacement curves of two cases are basically coincident at the beginning.However, due to the existence of interfacial debonding, a difference of displacement between the two cases gradually appears, and the steady-state value of displacement for the blade with debonding is lower than that without damage.However, the interface debonding has a significant effect on the transient value of temperature change and concentration, but no influence on the steady-state value.

TBCs with Delamination and Interfacial Debonding
In the complex service environment, delamination and interfacial debonding cracks often occur simultaneously in TBCs with functionally graded structure.This complex damage problem was considered and the aforementioned thermal barrier coatings structure was employed.Three cases were considered: (i) Case 1: the turbine blade without damage; (ii) Case 2: the turbine blade with debonding; (iii) Case 3: the turbine blade with delamination and interfacial debonding.The delamination at the interface between the 2nd and 3rd layers of TBCs and the size of the delamination area were the same as in the previous example, as shown in Figure 14.
At point B, the dynamic responses of displacement, temperature and concentration are shown in Figure 16.The corresponding nephograms of the TCM responses at 250 s and 1000 s are presented in Figure 17.By comparing Case 1 and Case 3, it is found that the delamination prevents the conduction of temperature and concentration, leading to a delay for the response curves to stabilize, and the displacement amplitude is smaller than that of Case 1.By comparing Case 2 and Case 3, it can be found that there are differences in the time-varying curves of u 1 , u 2 and u 3 at Point B, and this difference is gradually increasing with time.At the same time, the temperature and concentration response curves are also further affected, and the conduction process becomes slower, which is caused by the combined effect of delamination and debonding.The delamination has a significant effect on the transient value of the temperature change and concentration, but no influence on the steady-state value.

Coupling Analysis of Thermo-Chemo-Mechanical Field and Fluid-Structure Interaction
To ensure the structural stability of the turbine blade in a high-temperature environment, the gas-film holes and the internal cooling cavities can effectively reduce heat transfer.Moreover, in the flow field environment, when the gas-film cooling jet encounters the hightemperature and high-speed gas, a complex nonuniform temperature field is formed.This section conducts the coupling analysis of the TCM field and fluid-structure interaction of the turbine blade with the cooling structure.For this, we used another local turbine blade model with TBCs, whose geometric model and mesh model are shown in Figure 18.Firstly, commercial finite element software was used to calculate the flow-thermal coupling.The temperature of the outer mainstream gas was 1550 K, and the airflow velocity was 128 m/s.The temperature of the inner cold airflow was 850 K, and the airflow velocity was 106 m/s.The surface of the model structure, the inner flow channel wall and the inner wall interface of the film hole were nonslip walls, and the external flow field was set as periodic boundary conditions.An SST turbulence model was adopted for the turbulence model, and a turbulence intensity of 0.05 was selected.Only the internal cold airflow was considered to be mixed with the external mainstream gas through the cooling hole of the turbine blade.The nonuniform wall temperature load calculated by commercial software was introduced into the XLW/SE method to further calculate and analyze turbine blades.The uneven temperature field of the multifilm cooling hole wall is shown in Figure 19.As can be seen from Figure 19, because the blade adopts the cooling method of a porous film, the temperature of most wall areas in the jet direction of the cooling hole is generally reduced, and a better thermal protection effect is obtained.However, there is a continuous local high-temperature region in the leading edge and upper-and lower-edge regions of the wall facing the incoming flow, which is due to the fact that there are no more cooling holes in this region, so that it is difficult to organize an effective cooling.Two cases were considered: (i) Case 1: turbine blades with uniform temperature; (ii) Case 2: turbine blade with uneven temperature.For point C (1 cm, 5.01 cm and 3.24 cm) at the base interface near the cooling hole, the distributions of displacement, temperature vs. time are shown in Figure 20.It can be found that compared with a uniform temperature, the turbine blade at a nonuniform temperature is affected by air-film cooling, resulting in inconsistent temperature values in each part of the blade wall surface.The convergence rate of displacement and temperature with time is affected, and the amplitude of each physical quantity decreases to a certain extent.The corresponding nephograms at 500 s and 2000 s are presented in Figure 21.equation of the overall blade structure was assembled by defining the compatibility conditions of displacement, temperature, concentration and internal force equilibrium at the TBCs/substrate interface.The time integration was performed to calculate the TCM responses.The fluid-thermal coupling analysis was conducted to solve the flow field of the turbine blade with a cooling structure, and the obtained nonuniform thermal load was also input into the XLW/SE method.In the numerical examples, a computing program was adopted to calculate the dynamic responses of a turbine blade with and without debonding and delamination damages.By comparing the results of the XLW/SE model and a 3D elastic model developed by COMSOL, the effectiveness of the proposed XLW/SE method was validated.

Point C 2 cm
From a practical point of view, the XLW/SE method enabled us to solve the coupling effects of the complex curved surface structure with TBCs and analyze the impacts of damage and a nonuniform temperature field on the deformation displacement, concentration diffusion and temperature conduction.In addition, this study applied the proposed method to a real-world engine turbine blade, so as to obtain the dynamic coupling solution of the temperature field, concentration field, flow field and deformation field.Furthermore, for the subsequent damage prediction problem, the damage location and size could be further judged by the change of temperature and concentration, which had a certain practical application value.
The research work in this article is a foundational work, which needs further improvements to ensure the normal operation of turbine engines under harsher environmental conditions, for instance, the introduction of Eu3+ into yttria-stabilized-zirconia (YSZ) TBCs to improve the thermal insulation performance and interfacial fracture toughness of the coatings.In future work, the mechanical properties of this new material will be analyzed by a numerical calculation method (i.e., XLW/SE method), and verified through experiments.
In addition, this study mainly focused on temperature conduction, ion diffusion and deformation displacement.On the basis of existing research, the interface oxidation problem and calcium-magnesium-aluminum-silicate (CMAS) corrosion problem of thermal barrier coatings with complex curved surface structure can be studied, which can provide a guarantee for further revealing the failure mechanism of TBCs.Furthermore, the analysis results can provide guidance for improving the structural design of turbine blades.
The submatrices in the stiffness matrix of the chemical field are expressed as abζηke are the temperature-stress coefficients and concentration-stress coefficients, respectively.D1 abζηke and D1 abζηke are the heat conduction coefficients and diffusion coefficients, respectively.I ζηke is the density component coefficient.

Figure 4 .
Figure 4. Global coordinate system and local coordinate system of finite elements.
Displacement,temperature and concentration condition: Interal force balance condition: Substrate TBCs 1 1 chemical external load and chemical interaction force.

Figure 7 .
Figure 7. Schematic diagram of turbine blade with cooling structure and thermal barrier coating [41].
model construction by XLW method Substrate model construction by SE method Flow-thermal coupling analysis by Comsol Final governing equation assembly and solving

Figure 11 .
Figure 11.Displacements for three different time steps.

Figure 13 .
Figure 13.Nephograms of TCM responses of nondestructive TBCs structure by using XLW/SE method and COMSOL finite element software.

Figure 15 .
Figure 15.Displacement, temperature and concentration for TBCs structure with debonding at interface.
Displacement in the x direction (b) Displacement in the y direction (c) Displacement in the z direction (d) Temperature (e) Concentration

Figure 16 .
Figure 16.Displacement, temperature and concentration for TBCs structure with debonding and delamination.

Figure 17 .
Figure 17.Nephograms of TCM responses in Case 2 and Case 3.

Figure 18 .
Figure 18.(a) Geometrical model of turbine blade with multifilm cooling holes; (b) 3D mesh model of turbine blade with multifilm cooling holes.

Figure 19 .
Figure 19.Uneven temperature field of turbine blade with multifilm cooling holes.

Table 1 .
TCM Parameters of TBCs and substrate.

Table 2 .
Number of elements and nodes for three different meshing schemes.

Table 3 .
Maximum displacements of the 3D elastic model and XLW/SE model.