Bimetallic Thin-Walled Box Beam Thermal Buckling Response

A beam model for thermal buckling analysis of a bimetallic box beam is presented. The Euler–Bernoulli–Vlasov beam theory is employed considering large rotations but small strains. The nonlinear stability analysis is performed using an updated Lagrangian formulation. In order to account for the thermal effects of temperature-dependent (TD) and temperature-independent (TID) materials, a uniform temperature rise through beam wall thickness is considered. The numerical results for thin-walled box beams are presented to investigate the effects of different boundary conditions, beam lengths and material thickness ratios on the critical buckling temperature and post-buckling responses. The effectiveness and accuracy of the proposed model are verified by means of comparison with a shell model. It is revealed that all of the abovementioned effects are invaluable for buckling analysis of thin-walled beams under thermal load. Moreover, it is shown that the TD solutions give lower values than the TID one, emphasizing the importance of TD materials in beams.


Introduction
Thin-walled beams and structures are increasingly used in engineering branches, in standalone forms and as a stiffeners for plate-and shell-like structures, due to their high strength and light weight. However, these structures show susceptibility to local buckling and buckling failure [1,2]. Buckling analysis and the post-buckling response of such weightoptimized structures have been the topic of many research papers, such as [3][4][5][6], especially in the field of composite materials [7][8][9][10][11].
If the thermal environment is considered, the stability of structures has received significant attention in recent years: Duan et al. [12] performed thermal analysis of a beam element, Saha and Ali [13] presented a post-buckling mathematical model of a slender road under uniform temperature rise, while Cui and Hu [14] analyzed the thermal buckling and vibration of a beam. Jeyaraj et al. [15] investigated experimental and theoretical non-uniform heating of an isotropic beam. Burgreen and Mannit [16] and Burgreen and Regal [17] analysed the thermal buckling of bimetallic beams. In the case of composite beams, Aydogdu [18] obtained critical buckling temperatures of composite beams, Luan et al. [19] presented an analytical solution for buckling and vibration of FG beams, Kiani and Eslami [20,21] investigated buckling analysis under different types of thermal loads, while Giunta [22] analyzed FG beams under thermal/mechanical load using the Carrera unified formulation. However, there are not many papers about thermal buckling analyses of thin-walled structures: Libresceu [23] studied stability problems in a hightemperature environment and Ziane et al. [24] studied analytical methods for buckling and vibration responses of porous beams under thermomechanical loads.
In the present work, thermal buckling analysis of a thin-walled bimetallic box beam and frame structures is presented. The material is assumed to be linear, elastic and isotropic. The model is based on Euler-Bernoulli-Vlasov theory and on assumptions of large rotations and small strains. It is also assumed that the cross-section is not deformed in its own plane and that there are no shear strains in the middle surface. The nonlinear displacement field, which includes nonlinear displacement terms due to large rotation effects, is implemented. Using the UL description, the element geometric stiffness is derived. As an incremental iterative solution scheme, the Newton-Raphson method is used. Furthermore, this paper is a continuation of the research in which thermal buckling analysis of temperatureindependent materials was conducted [25], which has now been further expanded with temperature-dependent materials' properties. As far as the authors are aware, there is no beam model solution for thermal buckling analysis of thin-walled beam-type structures with temperature-dependent materials' properties.
The numerical results for thin-walled box beams are presented to investigate the effects of different boundary conditions, namely clamped-clamped, simply supported and clamped-simply supported, beam lengths and material thickness ratios on the critical buckling temperature and post buckling responses. In order to demonstrate the accuracy of the numerical algorithm, benchmark examples using shell FEM code were developed. Numerical results show that the abovementioned effects have a huge impact on the buckling analysis.

Kinematics
Two sets of coordinate systems related to the angle of orientation β are considered. The first one is a Cartesian (z, x, y) coordinate system where the z-axis coincides with the longitudinal beam that passes through the centroid O of each cross section, while the xand y-axes are principal axes. The second one is a contour coordinate system where the s-axis is tangential to the middle surface directed along the contour line of the cross-section while the n-axis is perpendicular to the s-axis.
The field of incremental displacement measures of a cross section are defined as [8]: where w 0 , u s and v s are the rigid-body translations of the cross-section centroid in the z-, xand y-direction, respectively, while ϕ z , ϕ x and ϕ y are the rigid-body rotations about the aforementioned axis; θ is a warping parameter of the cross-section.
In the case of small rotations, the incremental displacement field consists of the firstorder displacement values: where u z , u x and u y are the linear displacement increments of an arbitrary point on the crosssection defined by the x and y coordinates and the warping function ω(x, y). When the large rotations are considered, nonlinear displacement increments are expressed as follows: and should be added to those from Equation (2). Considering the nonlinear displacement field, the Green-Lagrange strain tensor components can be written as: where e ij and η ij are the linear and nonlinear strain components corresponding to the linear displacement, while e ij is the linear strain component corresponding to the nonlinear displacement due to the large rotations. The contour mid-line displacement w, and v can be seen more detail in [26]. Due to the in-plane rigidity hypothesis of the cross-section, the non-zero strain components are [7]: The stress resultants of the beam can be defined as: where F z represents axial force, M x and M y are bending moments with respect to xand y-axis, respectively, M t is the torsion moment and M ω is the warping moment (bimoment). t is the thickness of the closed section contour and F s is the St. Venant circuit flow [26].

Constitutive Equations
Consider a bimetallic beam made of two different metals. A beam wall with a core thickness of Ti-6Al-4V λt on the outer surface and SUS304 in the inner part of crosssection beam wall is shown in Figure 1. It is assumed that the layers of materials are perfectly bonded.  If the thermo-elastic material properties are considered as a fun T, they can be calculated for each material, as described in [19,27]: where P represents Young's modulus and thermal expansion coe , 1 , and are temperature-dependent coefficients listed in metals [28]. For simplicity, Poisson's ratio is assumed to be consta Table 1. Temperature dependent coefficients [19]. If the thermo-elastic material properties are considered as a function of temperature T, they can be calculated for each material, as described in [19,27]: where P represents Young's modulus E and thermal expansion coefficient α, while P 0 , P −1 , P 1 , P 2 and P 3 are temperature-dependent coefficients listed in Table 1 for different metals [28]. For simplicity, Poisson's ratio ν is assumed to be constant, ν = 0.3. Table 1. Temperature dependent coefficients [19].
It is assumed that the temperature of the whole beam is uniform and increased from the current ambient temperature T 0 to the critical value in incremental steps of 1 • C. If the axial displacements are prevented, the temperature at a point T(n, z) may be raised to T + ∆T, in the way that the beam buckles. ∆T is the temperature rise. The temperature that is read as the critical buckling temperature is the temperature difference compared to the ambient temperature. The process can be described as quasi-adiabatic since the heat exchange between the environment and the beam is neglected.
The stress-strain relations of the bimetallic beam can be written as: where Using Equations (5), (10) and (11), the beam forces can be expressed in a matrix form as: where R ij represents the thin-walled beam stiffness, as shown in Appendix A. N T z , M T x , M T y and M T ω are thermal force and thermal moments, respectively:

Finite Element Formulation
A two-nodded beam element with 14 degrees of freedom is shown in Figure 2. The nodal displacements and nodal force vectors are as follows: where the superscript e denotes the eth finite element. It should be noted that nodal displacement w and nodal forces F z , M x and M y are defined in the centroid O, while other nodal components are defined in the shear center.
where the superscript e denotes the eth finite element. It should be noted that nodal displacement and nodal forces , and are defined in the centroid , while other nodal components are defined in the shear center. Applying the principle of the virtual work, the incremental equilibrium equations of a beam element in linearized form are: where the equations from the left side consist of incremental virtual elastic strain energy: and the incremental virtual geometric potential: Applying the principle of the virtual work, the incremental equilibrium equations of a beam element in linearized form are: where the equations from the left side consist of incremental virtual elastic strain energy: and the incremental virtual geometric potential: On the right side of the equations, the terms represent the virtual work completed by external forces at the end and at the beginning of the present increment: In these equations, S ij is the second Piola-Kirchoff stress tensor, t i denotes the surface tractions, C ijkl presents the stress-strain tensor and the symbol δ indicates virtual quantities. By applying the linear interpolation functions for w 0 displacement and cubic interpolations for w s , u s and v s , one can obtain: Materials 2022, 15, 7537 where f e is the nodal force vector, k e E is the elastic stiffness matrix and k e G is the geometric stiffness matrix of the beam element. Nonlinear equilibrium equations are solved using the Newton-Raphson method as incremental iterative approach [29,30], and the explicit form of the terms given in nonlinear components were described previously in [31].

Results and Discussion
In numerical examples, the thin-walled box beam with height h = 100 mm, width b = 150 mm and thickness t = 10 mm is considered (Figure 3). For verification purposes, the critical buckling temperatures were obtained by shell FEM commercial code [32]. In order to simulate the bimetallic material, the beam walls were divided into two layers of different metals with a variable thickness ratio λ. Note that for λ = 0, the beam wall is fully SUS304, while as the index thickness ratio λ increases, the beam wall becomes fully Ti-6Al-4V.
where f is the nodal force vector, k is the elastic stiffness ma metric stiffness matrix of the beam element. Nonlinear equilibriu using the Newton-Raphson method as incremental iterative appr plicit form of the terms given in nonlinear components were desc

Results and Discussion
In numerical examples, the thin-walled box beam with hei = 150 mm and thickness = 10 mm is considered (Figure 3 poses, the critical buckling temperatures were obtained by shel [32]. In order to simulate the bimetallic material, the beam walls layers of different metals with a variable thickness ratio λ. Note wall is fully SUS304, while as the index thickness ratio λ increases fully Ti-6Al-4V.

Box Beam
In the first example, the eigenvalue results of the box beam for different boundary conditions, which are clamped-clamped (C-C), clamped-simply supported (C-S) and simply supported (S-S), beam lengths of L 1 = 6 m, L 2 = 8 m and L 3 = 10 m and different material thickness ratios λ are given in Tables 2-4. Critical buckling temperatures are given for the first two flexural buckling modes. As the thickness ratio λ increases, the critical buckling temperatures increase as well due to the material properties of TI-6Al-4V. Furthermore, longer beams obtain lower critical buckling temperatures. As expected, clamped-clamped beams exhibit the highest eigenvalues for every beam length. Good agreement of the present results and solutions derived from the 2D model for both flexural buckling modes is achieved. In order to further perform the nonlinear stability analysis of the box beam, a perturbation force of intensity ∆F = 500 N is introduced in the y-axis direction at the midpoint of the clamped-clamped beam. Temperature-displacement curves of the shortest beam with the comparison of temperature-dependent (TD) and temperature-independent (TID) materials are shown in Figures 4-6. The results are given for different thickness ratios: λ = 0.2, λ = 0.5 and λ = 0.8. It can be seen that curves match very well with the critical buckling temperatures achieved in the eigenvalue manner. As expected, TD solutions obtained lower critical buckling temperatures. In the case of λ = 0.2, the difference in critical temperature is around 4%, for λ = 0.5 the difference is 3.3%, and for λ = 0.8 it is 1.5%. In order to further perform the nonlinear stability analysis of the box beam, a perturbation force of intensity ∆ = 500 N is introduced in the y-axis direction at the midpoint of the clamped-clamped beam. Temperature-displacement curves of the shortest beam with the comparison of temperature-dependent (TD) and temperature-independent (TID) materials are shown in Figures 4-6. The results are given for different thickness ratios: = 0.2, = 0.5 and = 0.8. It can be seen that curves match very well with the critical buckling temperatures achieved in the eigenvalue manner. As expected, TD solutions obtained lower critical buckling temperatures. In the case of = 0.2, the difference in critical temperature is around 4%, for = 0.5 the difference is 3.3%, and for = 0.8 it is 1.5%.

L-Frame
Furthermore, the model is tested for thermal buckling analysis of an L-frame with the length of both legs being = 5 m and cross-section described in the previous chapter ( Figure 7). The frame is fixed at points A and C, while at point B, in-plane translations are prevented.

L-Frame
Furthermore, the model is tested for thermal buckling analysis of an L-frame with the length of both legs being L = 5 m and cross-section described in the previous chapter ( Figure 7). The frame is fixed at points A and C, while at point B, in-plane translations are prevented.  To verify the results, the critical buckling temperatures for the full SUS304, ∆T cr(SUS304) = 16.09°C, and full Ti-6Al-4V ∆T cr(Ti−6Al−4V) = 35.52°C, sections are solved by a shell commercial code. To perform nonlinear analysis, a small perturbation force ∆F = 50 N acting in the z-axis direction at point B is applied. The results are shown in Figure 8 for pure metals and for λ = 0.2t, λ = 0.5t and λ = 0.8t. The good recognition of the critical values can be noted. It can be observed that with an increase in the proportion of Ti-6Al-4V material, higher critical temperatures are achieved. SUS304, ∆ ( ) = 16.09 ℃, and full Ti-6Al-4V ∆ ( ) = 35.52 ℃, sections are solved by a shell commercial code. To perform nonlinear analysis, a small perturbation force ∆ = 50 N acting in the z-axis direction at point B is applied. The results are shown in Figure 8 for pure metals and for = 0.2 , = 0.5 and = 0.8 . The good recognition of the critical values can be noted. It can be observed that with an increase in the proportion of Ti-6Al-4V material, higher critical temperatures are achieved.

Conclusions
A thin-walled beam model capable of thermal buckling analysis has been presented. By means of the updated Lagrangian formulation, the incremental equilibrium equations have been developed using the nonlinear displacement field of the cross-section, taking into account the effects of large rotations. The reliability of the present model was verified by studying the benchmark examples, and the values obtained with the proposed model are in good agreement with those of the shell model. The effects of boundary conditions, the length of the beam and material thickness ratio on the critical buckling temperature

Conclusions
A thin-walled beam model capable of thermal buckling analysis has been presented. By means of the updated Lagrangian formulation, the incremental equilibrium equations have been developed using the nonlinear displacement field of the cross-section, taking into account the effects of large rotations. The reliability of the present model was verified by studying the benchmark examples, and the values obtained with the proposed model are in good agreement with those of the shell model. The effects of boundary conditions, the length of the beam and material thickness ratio on the critical buckling temperature and post-buckling response are of great importance. Additionally, it is shown that the TD solutions provide lower values than those of the TID solutions. The model was found to be efficient in predicting eigenvalues and nonlinear buckling behavior.