## 2 Traditional Matrix Structural Analysis (MSA)  
The traditional matrix structural analysis approach mainly focuses on the analysis of Euler-Bernoulli beam-columns. For first-order and second-order analysis (without and with considering the geometric nonlinearity), the element stiffness matrix of an Euler-Bernoulli beam-column element is formulated (e.g., [16; 17]) as Equations (1) and (2), respectively, showing the relationship between the element-end translational deformations/rotation angles (\(\delta_{1}\), \(\theta_{1}\), \(\delta_{2}\), \(\theta_{2}\)) and shear forces/bending moments (\(V_{1}\), \(M_{1}\), \(V_{2}\), \(M_{2}\)) (see Figure 1). Based on the element stiffness matrix, structural analyses can be conveniently conducted, including the traditional first-order analysis for the element deformations and internal forces using Equation (1) (e.g., [16; 17]), as well as elastic buckling and second-order stability analyses using Equation (2) (e.g., [13; 14; 15]). \[[\mathbf{K_{EB,lst}}]=\left[ \begin{array}{cccc}T_{\text{EB}}\frac{EI}{L^{3}} &Q_{\text{EB}}\frac{EI}{L^{2}}&-T_{\text{EB}}\frac{EI}{L^{3}}&Q_{\text{EB}} \frac{EI}{L^{2}}\\ &S_{\text{EB}}\frac{EI}{L}&-Q_{\text{EB}}\frac{EI}{L^{2}}&C_{\text{EB}}\frac{EI} {L}\\ &&T_{\text{EB}}\frac{EI}{L^{3}}&-Q_{\text{EB}}\frac{EI}{L^{2}}\\ \text{sym.}&&4\frac{EI}{L}\end{array} \right] \tag{1}\]  
\[[\mathbf{K_{EB,2 nd}}]=\left[ \begin{array}{cccc}T_{\text{EB}}\frac{EI}{L^{3}} &Q_{\text{EB}}\frac{EI}{L^{2}}&-T_{\text{EB}}\frac{EI}{L^{3}}&Q_{\text{EB}} \frac{EI}{L^{2}}\\ &S_{\text{EB}}\frac{EI}{L}&-Q_{\text{EB}}\frac{EI}{L^{2}}&C_{\text{EB}}\frac{EI} {L}\\ &&T_{\text{EB}}\frac{EI}{L^{3}}&-Q_{\text{EB}}\frac{EI}{L^{2}}\\ \text{sym.}&&5_{\text{EB}}\frac{EI}{L}\end{array} \right]\qquad \begin{array}{ c}T_{\text{EB}}=\lambda_{\text{EB}}^{3}\sin\lambda_{\text{EB}}/\phi_{\text{EB}} \\ &Q_{\text{EB}}=\lambda_{\text{EB}}^{2}(1-\cos\lambda_{\text{EB}})/\phi_{\text{EB }}\\ &S_{\text{EB}}=\lambda_{\text{EB}}(\sin\lambda_{\text{EB}}-\lambda_{\text{EB}} \cos\lambda_{\text{EB}})/\phi_{\text{EB}}\\ &C_{\text{EB}}=\lambda_{\text{EB}}(\lambda_{\text{EB}}-\sin\lambda_{\text{EB}} \cos\lambda_{\text{EB}})/\phi_{\text{EB}}\\ &\phi_{\text{EB}}=2-2\cos\lambda_{\text{EB}}-\lambda_{\text{EB}}\sin\lambda_{ \text{EB}}\end{array} \right. \tag{2}\]  
where \(EI\) denotes the cross-sectional flexural rigidity; \(L\) denotes the element length; \(T_{\text{EB}}\), \(Q_{\text{EB}}\), \(S_{\text{EB}}\), and \(C_{\text{EB}}\) are element stability functions of factor \(\lambda_{\text{EB}}\); \(\lambda_{\text{EB}}\) is an axial force factor for Euler-Bernoulli beam, defined as Equation (3); \(P_{\text{E}}=\pi^{2}EI/L^{2}\) is the Euler buckling axial force. \[\lambda_{\text{EB}}=\sqrt{\frac{PL^{2}}{EI}}=\pi\sqrt{\frac{P}{P_{\text{E}}}} \tag{3}\]For the analysis of Timoshenko beam-columns, the cross-sectional shear deformations are considered, and the element stiffness matrix is a little more complex than that of Euler-Bernoulli beam-columns. The first-order element stiffness matrix (without considering geometric nonlinearity) has been widely known, as shown in Equation (4) [18; 19]. However, for the second-order analysis (including the elastic buckling analysis), the exact solution of the element stiffness matrix considering geometric nonlinearity has not been derived. [PERSON] et al. [18] provided a simplified form (Equation (5)) for this second-order element stiffness matrix (also denoted as stability stiffness matrix) of Timoshenko beam-columns, based on a modification of the second-order element stiffness matrix of Euler-Bernoulli beam-columns (Equation (2)). However, this simplified form uses modified element stability functions, which could have considerable errors for elements with large axial force. Thus, conducting (elastic buckling and second-order) analyses of axial-loaded Timoshenko beam-columns based on this simplified form could lead to incorrect results with large errors.