Finite strip-Riccati transfer matrix method for buckling analysis of tree-branched cross-section thin-walled members

Thin-walled components are gaining wide application in the field of modern engineering structures. The buckling analysis of thin-walled structures has thus become an important research topic. Here, we developed a new method named as finite strip-Riccati transfer matrix method (Riccati FSTMM) for buckling analysis of tree-branched cross-section thin-walled members. The method integrates Riccati transfer matrix method (Riccati TMM) for tree multi-body system with semi-analytical finite strip method (SA-FSM). Compared to SA-FSM, Riccati FSTMM features a smaller matrix and higher calculation efficiency, with no need for global stiffness matrix. In addition, by arranging uniformly distributed middle nodal lines inside strip elements, we developed the high order finite strip-Riccati transfer matrix method (Riccati HFSTMM) for buckling analysis of tree-branched cross-section thin-walled members. This method further improves the efficiency and accuracy of Riccati FSTMM. We tested the two proposed methods with two numerical examples, and demonstrated their superior reliability and efficiency over the finite element method (FEM).


Introduction
In recent decades, plates and thin-walled members have been widely used in engineering fields such as military industry, aerospace domain, and building construction. Thin-walled structures are light weight and have high stiffness to weight ratio. Given these engineering properties as well as energy efficiency, thin-walled structures are expected to gain wide range of applications.
Structural buckling is an important aspect in designing and studying of thin-walled structures. The investigation of buckling instability can be carried on by numerical, analytical, and experimental technologies. 1 Nowadays, there are various methods for buckling analysis of thinwalled members, such as finite element method (FEM), 2 finite strip method (FSM), 3 boundary element method, 4 finite difference method, 5 Direct Strength Method (DSM), 6 Generalized Beam Theory (GBT), 7

etc.
In modern engineering fields, FEM is a powerful method that can solve most engineering problems. However, FEM requires a large number of unknowns to approximate accurate results for discretizing structures in all directions, and is thus time consuming and resource intensive. For regular geometry shape thin-walled structures, finite strip method (FSM) is an effective method to replace FEM in structural analysis, as it can reduce the element number considerably and improve computational efficiency greatly. Nowadays, a variety of finite strip methods have been proposed for the buckling problems of thin-walled members. Based on semi-energy approach, the semi-analytical finite strip method (SA-FSM) was developed for the buckling, 8 shear buckling, 9 dynamic buckling of plate structures. 10 And SA-FSM has been widely used in computer software (CUFSM, 11 THIN-WALL 12 ), with the purpose of developing the signature curves of the buckling stress versus buckling halfwavelength for thin-walled members. Based on SA-FSM, the modal decomposition of thin-walled members with rounded corners can be discussed by the constrained finite strip method (cFSM), 13 which is to transform the nodal basis system of FSM into a modal basis system. More recently, a technique 13 was proposed to eliminate the major problems caused by the rounded corners in the original cFSM. Spline finite strip method (SFSM) 14 utilizes piece-wise polynomial shape functions to represent the longitudinal field, which can adapt to more complex boundary conditions and loading compared to SA-FSM. Taking advantages of SFSM, constrained spline finite strip method (cSFSM) 15 allowing the modal decomposition have been proposed by integrating GBT principles in SFSM. Besides, the compound strip method (CSM) can overcome the difficulty in traditional FSM capturing the influence of stiffeners. 16 Furthermore, the elastic buckling analyses of channel sections subject to pure shear stresses can be conducted by SAFSM and SFSM. 17 And the buckling analyses of both square plates and lipped channel sections with central circular and square holes subject to shear have been studied by SFSM. 18 To our knowledge, all these kinds of finite strip methods require global stiffness matrix, which could be a lot of computations.
The transfer matrix method for multibody systems (MSTMM), which features no global dynamics equations, high programming, low order of system matrix, and high computational speed, is a powerful numerical analysis approach for multibody systems with various topological structures. 19 To further enhance the efficiency of SA-FSM in analyzing the buckling problem of thin-walled members, semi-analytical finite strip transfer matrix method (FSTMM) 20 was advanced by introducing MSTMM in SA-FSM. The method requires no global stiffness, thus greatly reduces the computational load while deriving all the benefit of SAFSM. Previously, FSTMM has been applied in the buckling analysis of thin plates 21 and open-branched section thin-walled members 22 with simply-simply (SS) supported boundary condition of loaded edges. For thin-walled members with complicated cross section under other support conditions, the buckling analysis by FSTMM is under study. In addition, owing to successive multiplication of transfer matrices for largescale systems, some elements in overall transfer matrix tend to be large value, which may result in numerical instabilities. [23][24][25] This problem can be dealt with Riccati transfer matrix method (Riccati TMM). 26 However, the characteristic equations established by Riccati TMM may have asymmetric poles. When zero search technique based on the sign change (such as bisection method) is used to solve such characteristic equations, these asymmetric poles will be mistaken as the solution. In order to improve the search efficiency of eigenvalues, the techniques to eliminate the poles of characteristic equation were proposed. 25,27 Besides, a recursive eigenvalue search algorithm for TMM of linear flexible multibody systems can handle this problem. 25 In order to extend the application range of FSTMM and lift its calculation ability, we developed a new method named as finite strip-Riccati transfer matrix method (Riccati FSTMM) for buckling analysis of tree-section thin-walled members. This method can efficiently solve the buckling problem of complicated tree-section thinwalled members under various support conditions. In this paper, we introduced the Riccati TMM for tree multibody systems 26 in FSTMM, and the transfer matrix of branch in tree-section thin-walled member was defined. On this basis, we introduced the high order finite strip method 28 is introduced, and presented the high order finite strip-Riccati transfer matrix method (Riccati HFSTMM) for buckling analysis of tree-section thinwalled members. Compare with Riccati FSTMM, Riccati HFSTMM can significantly improve computational accuracy by adding little computation. The rest of the paper is organized as follows: the general theorem of the Riccati FSTMM for buckling analysis of tree-branched cross-section thin-walled members in Section ''Finite strip-Riccati transfer matrix method for buckling analysis of tree-branched crosssection thin-walled members,'' the Riccati HFSTMM in Section ''High order finite strip-Riccati transfer matrix method for buckling analysis of treebranched cross-section thin-walled members,'' two numerical examples and the comparison between our proposed methods and FEM in Section ''Examples and analysis,'' and the conclusions in Section ''Conclusions.'' Finite strip-Riccati transfer matrix method for buckling analysis of tree-branched cross-section thin-walled members Degree of freedom and shape function Take a T-section member as shown in Figure 1(a) for example, discretizing the member into s strips along the longitudinal direction, therefore, the total number of nodal lines is s + 1. Both the global coordinate denoted as X-Y-Z and local coordinate denoted as x-yz are introduced. Figure 1(b) shows two membrane degrees of freedoms (DOFs) u and v, and two bending DOFs w and u of each nodal line. And the longitudinal displacements can be represented by shape functions according to the corresponding boundary conditions. Two support conditions are mentioned in this paper, that is, one end simply supported and the other end clamped (SC) and the other one end clamped and the other end free (CF). The displacement function of SC can be expressed as equation (1), and the one of CF can be written as equation (2).
where a is the length of the member, p is the axial halfwave number, y is the longitudinal coordinate in local coordinate.

Fundamental stiffness matrix
For the detailed derivation of the elastic stiffness matrix and the geometric stiffness matrix of the strip, refer to Hou et al. 16 The (8 3 8) block elements k pq e of the elastic stiffness matrix k e can be given as the combination of membrane and bending terms, where k pq eM and k pq eB are the (4 3 4) membrane and bending elastic stiffness matrices correspondingly.
When the boundary condition is not simply supported at both load edges, the elastic stiffness matrix k e is not orthogonal. And the submatrices of k e vary as the halfwave number changes. In that case, it is difficult to represent k e with a specific matrix. Therefore, a generalized elastic stiffness matrix can be established as follows, where, the (8 3 8) block elements k pq e of the (8m 3 8m) generalized elastic matrix k e vary with p and q.
Similarly, the (8 3 8) block elements k pq g of the geometric stiffness matrix k g can be expressed as the combination of membrane and bending terms, where k pq gM and k pq gB are the (4 3 4) membrane and bending geometric stiffness matrices, respectively.
Accordingly, the generalized geometric stiffness matrix can be written as, where, the (8 3 8) block elements k pq g of the (8m 3 8m) generalized geometric matrix k g vary with p and q.

Control equations of strip element
According to virtual work principal, the control equations of strip element can be obtained, where, k e and k g is generalized elastic stiffness matrix and geometric stiffness matrix respectively, D m denotes the generalized displacements vector, and R m is the generalized internal forces. Assuming that the initial axial force is T a , the real axial force can be expressed as T = lT a (l is the buckling coefficient). As a result, the geometric stiffness matrix can be represented as k g = lk g (T a ) , substitute it into the equation (7) to get where the coefficient matrix of generalized displacements vector D m is k = k e À lk g(T a ) , and Þ . For convenience of calculations, rearrange the generalized displacements vector D m and generalized internal forces R m , where the subscript r represent re-arrangement, and the rest quantities of R m r are consistent with P i . Rewrite equation (8) by above reorganization, and the generalized control equations of the strip element can be expressed as, where the coefficient matrix K (8m 3 8m) can be deduced by the corresponding determinant transformation of coefficient matrix k in equation (8).

State vector and transfer matrix
The state vector of nodal line l covers generalized displacements and generalized internal forces, Where i, j denote the input and output nodal lines of each strip, D M, l = ½u l , v l T and D B, l = ½w l , u l T are the displacement vectors in nodal line l, R M, l = ½P l , Q l T and R B, l = ½F l , M l T define the generalized internal force vectors in the line l, the subscript n is the serial number of strips. By equation (11) and the definition of state vectors of nodal line in equation (12), the transfer equation of the strip n is given, To facilitate the expression of the transfer matrix method, the (8m 3 8m) coefficient matrix K in equation (11) can be represented as the following block matrix form, where K ij (i, j = 1, 2, 3, 4) are the (2m 3 2m) submatrices of the coefficient matrix K.
As mentioned above, the expression of the (8m 3 8m) transfer matrix Z n of the strip n can be obtained, Tree delivery strategy and Riccati transformation As illustrated in Figure 2, the transfer path of treesection thin-walled member is leaf-to-root, and each node is graded in sequence along the pathway. Thereout, the root node is on the top level while the leaf nodes are on the bottom level. In addition, the number of leaf nodes differs at different levels and leaf nodes are nonoverlapping in the tree. Furthermore, intermediate node can have multiple subordinate knots but only one ancestor node in path. Since the root node of a multi-branches tree is unique, the transfer path can be determined once the root nodal line is located.
In dynamics system, transfer matrix method for multi-body system (MSTMM) classifies the structural elements into two types: SISO element and MISO element. The former is the element with single input and single output end, whereas the latter is the element with multiple input ends and single output end. And the SISO element can be treated as a particular case of MISO element. As shown in Figure 3, two kinds of equations are involved for the MISO element n with L inputs, Transfer equations of strip element. For tree-section thinwalled members, the transfer equations of MISO element n with L inputs are described below, where, i k denotes the input nodal line of the kth input strip connected to the branch point, I k denotes the kth input end of the MISO element n, and I denotes the input nodal line of the MISO element n.
Geometric equation of strip element. Different from singlebranch thin-walled members, tree-section thin-walled structures possess the geometric equation defining the geometric relationship between input ends of the MISO element, where, H I k , n is a constant matrix. The dimension of state vectors of the element with L inputs and single output is 4m 3 1. We introduce Riccati transformation here and break the state vectors d into 2m 3 1 known quantities and 2m 3 1 unknown variables.
where A 1 , A 2 represent the known quantities and unknown variables respectively. Therefore, equations (16) and (17) can be rewritten as, where the (4m 3 4m) matrix Z ab (a, b = 1, 2) can be obtained by reordering the elements of Z k , H 1, I k , n , and H 2, I k , n are the submatrices of H I k , n .
Given the angle between strip n and strip n + 1 shown in Figure 4, the relationship of the vectors on the common nodal line of two adjacent strips can be indicated by the following transformation,    where a is the angle between conjunctional strips, T n represents the transformation matrix of strip n and strip n + 1 at the joint nodal line. The combination of equations (19) and (21) can be given as, where, I denotes the input nodal line of the MISO element n, the (4m 3 4m) matrix Z 00 ab (a, b = 1, 2) is the combination of the transfer matrix Z ab and the transformation matrix T n . Equation (22) can be decomposed as, For the element n with L inputs and single output, continuous conditions of displacements and equilibrium conditions of forces should be considered. Given that the displacement of each end connected to the branch point is the same, only the displacement of the first input end will be taken. And according to equilibrium conditions, the sum of the inner forces of L inputs is equal to the internal forces of the output end, therefore, all the forces of input ends are reserved. Based on the above analysis, equation (24) can be simplified as, The Riccati transformation can be introduced here, where S i k is the Riccati transfer matrix of the input nodal line of the kth input strip. Hence, equations (23) and (25) can be rewritten as, Moreover, it can be obtained from equation (20), where H I k , n = (H 1 S + H 2 ) À1 I k , n (H 1 S + H 2 ) I 1 , n , As regards the geometric equation of the tree-section thin-walled member in this paper, Combine equations (22) and (29), Substitute the equation (31) into equation (27), Combine equations (28) and (32), Hence, the Riccati recursive equation of MISO element with L inputs can be given, And we can remove 2k input ends to get the Riccati recursive equation of SISO elements, For the buckling analysis of tree-section thin-walled members only the case that unloaded edges are free (ff) is considered, and the constraint condition is 0 0 0 0 u v w u ½ T . Thus, the initial value of the matrix S is defined as, where the subscript i denotes the input nodal lines of tree-section thin-walled member, the subscript f is the leaf strip elements in transfer path of tree-section thinwalled member. Along the path, repeatedly use the equations (34) and (35) and end up with S s of the terminal nodal line of the last strip.
where A 1, j, s is a known zero vector, A 2, j, s is the unknown quantity, the subscript s represents the last strip in the path. To make the non-zero solutions of the equation (37) possible, the following condition must be satisfied.
The above equation is the characteristic equation for buckling of tree-section thin-walled member under generalized boundary condition.
High order finite strip-Riccati transfer matrix method for buckling analysis of tree-branched cross-section thin-walled members Degree of freedom of high order strip As shown in Figure 5, we introduce the numbering system and right-handed coordinate system similar to the previous section. In Figure 5(b), the dotted line represents the middle nodal line of strip, and i, in, j denote three nodal lines in a same strip. There are two membrane degrees of freedoms (DOFs) u i and v i , and two bending DOFs w i and u i in each nodal line. The displacement function of strip with one middle nodal line is identical to that presented in Section ''Degree of freedom and shape function,'' and it will not be repeated again.

Fundamental stiffness matrix
For the detailed derivation of the elastic stiffness matrix and the geometric stiffness matrix of high-order finite strip, refer to Rui et al. 19 The ( where, the (12 3 12) block elements k pq e of the (12m 3 12m) generalized elastic matrix k e of high-order strip vary with p and q.
Similarly, the (12 3 12) block elements k pq g of the geometric stiffness matrix k g of high-order strip can be expressed as the combination of membrane and bending terms, k where k pq gM and k pq gB are the (6 3 6) membrane and bending geometric stiffness matrices, respectively.
And the generalized geometric stiffness matrix can be written as, k g = k where, the (12 3 12) block elements k pq g of the (12m 3 12m) generalized geometric matrix k g of highorder strip vary with p and q.

Control equations of high order strip
According to virtual work principal, the control equations of high-order strip can be obtained, where, k e and k g is generalized elastic stiffness matrix and geometric stiffness matrix correspondingly, D m and R m denote the generalized displacements vector and the generalized internal forces of high-order strip, respectively. Assuming that the initial axial force is T a , the real axial force can be expressed as T = lT a (l is the buckling coefficient). As a result, the geometric stiffness matrix can be represented as k g = l k g (T a ) , substitute it into the equation (43) where the coefficient matrix of generalized displacements vector Þ . For ease of calculation, rearrange the generalized displacements vector D m and the generalized internal where the subscript r denotes re-arrangement, where, the (8m 3 8m)Z e n is the transfer matrix of highorder strip element, the subscript I and O represent the input and output nodal line respectively, the subscript n denotes the serial number of strips.
The stiffness matrix K in equation (47) can be rewritten as the following block matrix form,  From equation (52), it can be obtained as follows, The combination of equations (52) and (53) can be given as, where Combine the equations (53) and (54) and equation (55) is obtained, The expression of the bending term can be given by equation (51) in the same way, where The combination of equations (55) and (56) can be written as, where, where, the (12m 3 8m) matrix z n can be expressed as following form, Riccati recursive equations of high order strip elements The tree delivery strategy and Riccati transformation have been provided in Section ''Finite strip-Riccati transfer matrix method for buckling analysis of tree-branched crosssection thin-walled members'' and the detailed process will not be repeated here. The derivation of the Riccati recursive equations of high-order strip elements is entirely consistent with the part in Section ''Tree delivery strategy and Riccati transformation.'' Substitute the transfer matrix Z n in Section ''Finite strip-Riccati transfer matrix method for buckling analysis of tree-branched crosssection thinwalled members'' with the transfer matrix Z e n of high order strip and repeat the procedure of derivation, and the Riccati recursive equation of high-order strip element n with L inputs and single output can be obtained, where the subscript I denotes the input nodal lines, the subscript f is the leaf strip elements. Along the path, repeat the equations (61) and (62) and end up with S s of the terminal nodal line of the last strip.
where A 1, O, s is a known zero vector, A 2, O, s is the unknown quantity, the subscript s represents the last high-order strip in the path. To make the non-zero solutions of equation (64) possible, the following condition must be satisfied Equation (65) is the characteristic equation for buckling analysis of tree-section thin-walled members by Riccati HFSTMM.

Examples and analysis
In this section, two examples are given to validate and illustrate the reliability of Riccati FSTMM and Riccati HFSTMM, including a T-section member and a H-section member. The material parameters throughout this paper are as follow, Young modulus E = 2 3 10 5 N=mm 2 , Poisson's ratio v = 0:3, and shear modulus G = E=2(1 + v). The results by these two methods will be compared with the ones by FEM, that proofs the validity of the proposed methods. And the computational efficiency of two proposed methods will be discussed in Comparative analysis.

Illustrations of T cross-section members
The dimensions of T-shaped member are presented in Figure 6(a). The section height is 100 mm, the flange width is 80 mm, the flange lip length is 20 mm, the plate thickness is 2 mm, the initial force T a = 50N=mm 2 . As shown in Figure 6(b) and (d), the member analyzed by Figure 6. Section diagram of T-section thin-walled member.
Riccati FSTMM is divided into 7 and 13 strips along the loaded edge, the input and output nodal lines are denoted by solid points, and the transfer path can be determined by tree delivery strategy. Figure 6(c) shows the division and path of the member solved by Riccati HFSTMM with seven strips, the hollow points denote the middle lines. Use the proposed methods and FEM to calculate T-section members under the boundary conditions of SC and CF at the loaded edges, and obtain their buckling curves respectively. Figure 7 plots the relationship diagram of the buckling coefficient l and member length a of the member of SC at the loaded edges. There are two distinct buckling regions existing in this example, local and global buckling modes. The buckling shape transfers from local to global modes at a = 700 mm. And Figure 7 gives the local and global buckling shape of this example when the length is 120 and 800 mm individually.
The characteristic curves of the member of CF at the loaded edges is presented in Figure 8. The buckling shape transfers from local to global modes at a = 300 mm. And this figure shows the local and global buckling shape for the example at a = 120 mm and a = 800 mm.
It can be found from Figures 7 and 8 that, when Tsection member is divided into seven strips, the fit between the buckling curve obtained by Riccati FSTMM and the one by FEM is not perfect. When the number of strips increases to 13, the curve is in better agreement with FEM. While Riccati HFSTMM is used, the buckling curve can coincide with the one by FEM well when the member is divided into seven strips.

Illustrations of H cross-section members
As shown in Figure 9(a), the geometrical properties of H-shaped are as follows, the section height is 100 mm, the web is 50 mm, the flange lip length is 20 mm, the plate thickness is 2 mm, and the initial force T a = 50N=mm 2 . The section diagram and transfer path of this member solved by Riccati FSTMM with 9 strips and 19 strips are presented in Figure 9(b) and (d), correspondingly. The member analyzed by Riccati HFSTMM is subdivided into 19 strips, and its path is illustrated in Figure 9(c). Calculated the member with Riccati FSTMM, Riccati HFSTMM and FEM in turn, and obtain the buckling curves of H-section member under the boundary conditions of SC and CF at the loaded edges, as shown in Figures 10 and 11 respectively.
The signature curves of the buckling coefficient l versus length a for H-shaped member of SC are plotted in Figure 10. Note that the distortional region exists in this loading condition. The buckling shape transfers from local to distortional mode at a = 400 mm, and turn into global mode at a = 1430 mm. The buckling shapes corresponding to local, distortional, and global modes with lengths 140, 800, and 3000 mm are presented in this figure.
As shown in Figure 11, three kinds of buckling modes can be found in H-section member of CF. The buckling shape transfers from local to distortional mode at a = 130 mm, and turn into global buckling mode at a = 300 mm. And Figure 11 gives the local, distortional, and global buckling shapes for the  member of CF when the length is 70, 300, and 1080 mm respectively. It should be noticed that, there are a little difference between these four curves when the member length is relatively small, and the discussion about this phenomenon will be given in the next section.

Comparative analysis
Comparative analysis for computational precision. In order to show directly the computational accuracy of different methods, the results by Riccati FSTMM with two kinds of division and the ones by Riccati HFSTMM are compared with the results by FEM. The calculation formula for relative error is given below, where, l e is the buckling coefficient calculated by proposed methods, l FEM is the buckling coefficient obtained by FEM. The relative errors of T-section and H-section members of SC and CF by different methods are presented in Tables 1 to 4 respectively. As shown in Tables 1 and 2, there is large error between the results of T-section member by Riccati FSTMM with seven strips and the ones by FEM. And the results by Riccati FSTMM with 13 strips are approximate to the value by FEM. For the Riccati HFSTMM with 7 strips, however, the errors are smaller than the ones by Riccati FSTMM with 13 strips.   As for H-section member, we can notice from Tables 3 and 4 that, the relative error between Riccati FSTMM with nine strips and FEM is comparing large. And the error is narrowing when the strips in Riccati FSTMM increases to 19. Nevertheless, the relative error of Riccati HFSTMM with nine strips is minimal.
Comparative analysis for computational precision. This part records the computational time of the T-section member of SC in length a = 500 mm in different methods. From Table 5, it can be observed firstly that in Riccati FSTMM and Riccati HFSTMM of the modals of 13 strips need more computational time than the models of 7 strips, from 11.83 to 8.79 s in Riccati FSTMM and from 11.98 to 9.02 s in Riccati HFSTMM, respectively. Second, when the same number of strips can be used, the computing time of Riccati HFSTMM is a little longer than that of Riccati FSTMM, and in the case of 7(13) strips, Riccati HFSTMM requires 9.02 s (11.98 s) and Riccati FSTMM needs 8.79 s (11.83 s). Finally, we can notice that Riccati HFSTMM with 7 strips demands less computing time than Riccati FSTMM with 13 strips, and it owns high accuracy as shown in Table 1.
To sum up, it can be concluded that the proposed methods are both valid, and Riccati HFSTMM is more accurate than Riccati FSTMM for its efficiency independent of the strip width. Moreover, Riccati HFSTMM needs less computational time with necessary accuracy. Last but not least, the small difference exists between Riccati HFSTMM or Riccati FSTMM with enough strips and FEM for the reason that shear stain is included in FEM but neglect in proposed methods.

Conclusions
In this article, we presented the finite strip-Riccati transfer matrix method (Riccati FSTMM) for buckling analysis of tree-section thin-walled members under generalized boundary conditions. On this basis, with a middle nodal line arranged in strip element we developed the high order finite strip-Riccati transfer matrix method (Riccati HFSTMM) for buckling analysis of tree-section thin-walled member. In order to test these two methods, we analyzed the buckling problem of Tsection and H-section members under two boundary conditions. We compared the results with that of FEM and showed superiority of our proposed methods. Specifically, the computational accuracy of Riccati FSTMM increases with the number of strips, and the accuracy of Riccati HFSTMM can be assured even when the member is divided into less strips.