Advanced Analysis of Space Steel Frames

1. Introduction

This chapter presents advanced analysis methods for space steel frames which consider both geometric and material nonlinearities. The geometric nonlinearities come from second-order $P-\Delta$ and $P-\delta$ effects (see Fig. 1.) as well as geometric imperfections, while the material nonlinearities are due to gradual yielding associated with residual stresses and flexure. The $P-\Delta$ effect results from the axial force acting through the relative displacement of the ends of the member, so it is referred to as a member chord rotation effect. The $P-\Delta$ effect is accounted in the second-order analysis by updating the configuration of the structure during the analysis process. The $P-\delta$ effect is caused by the axial force acting through the lateral displacement of the member relative to its chord, so it is referred to as a member curvature effect. The $P-\delta$ effect can be captured by using stability functions. Since the stability functions are derived from the closed-form solution of a beam-column subjected to end forces, they can accurately capture the $P-\delta$ effect by using only one element per member. Another way to capture the $P-\delta$ effect without using stability functions is to divide the member into many elements, and consequently, the $P-\delta$ effect is transformed to the $P-\Delta$ effect.

Geometric imperfections result from unavoidable errors during the fabrication or erection. There are three methods to model the geometric imperfections: (1) the explicit imperfection modeling, (2) the equivalent notional load, and (3) the further reduced tangent modulus. The explicit imperfection modeling for braced and unbraced members is illustrated in Fig. 2(a). For braced members, out-of-straightness is used instead of out-of-plumbness. This is due to the fact that the $P-\Delta$ effect due to the out-of-plumbness is vanished by braces. The limitation of this method is that it requires the determination of the direction of geometric imperfections which is often difficult in a large structural system. In the equivalent notional load method, the geometric imperfections are replaced by equivalent notional lateral loads in proportion to the gravity loads acting on the story as described in Fig. 2(b). The drawback of this method is that the gravity loads must be known in advance to determine the notional loads before analysis. Another way to account for the geometric imperfections is to further reduce the tangent modulus. The advantage of this method over the explicit imperfection modeling and equivalent notional load methods is its convenience and simplicity because it eliminates the inconvenience of explicit imperfection modeling and equivalent notional load methods.

Figure 1.

The $P-\delta$ and $P-\Delta$ effects of a beam-column

Figure 2.

Geometric imperfection methods

Residual stresses are created in the hot-rolled sections due to uneven cooling of the cross-section. Typical residual stress pattern for a hot-rolled wide flange section is illustrated in Fig. 3. When a member is subjected to a compressive force, the fibers which have the highest values of compressive residual stress will yield first, and the fibers with the tensile stress will yield last. It means that the yielding over the cross-section is a gradual process. Hence, the stress-strain curve for a stub column is smooth instead of linear elastic-perfectly plastic in the case of coupon as shown in Fig. 4(a). The gradual yielding over the cross-section is caused not only by residual stress but also by flexure as shown in Fig. 4(b). Although the stress-strain relationship of steel is assumed to be linear elastic-perfectly plastic, the moment-curvature relationship has a smooth transition from elastic to fully plastic. This is because the section starts to yield gradually from extreme fibers which have the highest stresses. Material nonlinearities can be taken into account using various methods based on the degree of refinement used to represent yielding. The elastic plastic hinge method allows a drastic simplification, while the plastic zone method uses the greatest refinement.

Figure 3.

Typical residual stress pattern for a hot-rolled wide flange section

In the current design approach, the strength and stability of a structural system and its members are treated separately, and hence, the information about the failure modes of a structural system is not provided. This disadvantage is overcome by using a second-order inelastic analysis called “advanced analysis”. Advanced analysis indicates any methods that efficiently and accurately capture the behavior and the strength of a structural system and its component members. This chapter will present two advanced analysis methods: (1) the refined plastic hinge method and (2) the fiber method. In these methods, the geometric nonlinearities are captured using the stability functions, while the material nonlinearities are considered using the refined plastic hinge model and fiber model. The benefit of employing the stability functions is that it can accurately capture geometrical nonlinear effects by using only one element per member, and hence, this leads to a high computational efficiency as demonstrated by the works of Thai and Kim (2008; 2009; 2011b; 2011c; 2011d; 2012).

Figure 4.

Gradual yielding of steel member

2. Advanced analysis

2.1. Stability functions accounting for second-order effects

Considering a beam-column element subjected to end moments and axial force as shown in Fig. 5. Using the free-body diagram of a segment of a beam-column element of length x, the external moment acting on the cut section is

$$M_{ext}=M_A+Py-\frac{M_A+M_B}{L}x=-EI{y}^{″}$$E1

where$E$, $I$, and $$L$$ are the elastic modulus, moment of inertia, and length of an element, respectively.

Figure 5.

Beam-column with double-curvature bending

Using$$k^2=P/EI$$, Eq. (1) is rewritten as

$$y^{″}+k^{2}P=\frac{M_A+M_B}{EIL}x-\frac{M_A}{EI}$$E2

The general solution of Eq. (2) is

$$y=C_1\sin kx+C_2\cos kx+\frac{M_A+M_B}{EIL{k}^2}x-\frac{M_A}{EI{k}^2}$$E3

The constants $$C_1$$ and $$C_2$$ are determined using the boundary conditions $$y\left(0\right)=y\left(L\right)=0$$

$$C_1=-\frac{M_A\cos kL+M_B}{EI{k}^2\sin kL}$$and $$C_2=\frac{M_A}{EI{k}^2}$$

Substituting Eq. (4) into Eq. (3), the deflection $y$ can be written as

$$y=-\frac{1}{EI{k}^2}\left[\frac{\cos kL}{\sin kL}\sin kx-\cos kx-\frac{x}{L}+1\right]M_A-\frac{1}{EI{k}^2}\left[\frac{1}{\sin kL}\sin kx-\frac{x}{L}\right]M_B$$E5

and rotation $y^{\prime }$ is given as

$$y^{\prime }=-\frac{1}{EIk}\left[\frac{\cos kL}{\sin kL}\cos kx+\sin kx-\frac{1}{kL}\right]M_A-\frac{1}{EIk}\left[\frac{1}{\sin kL}\cos kx-\frac{1}{kL}\right]M_B$$E6

The end rotation $${\theta }_A$$ and $${\theta }_B$$ can be obtained as

$$\begin{array}{l}{\theta }_A=y^{\prime }\left(0\right)=-\frac{1}{EIk}\left[\frac{\cos kL}{\sin kL}-\frac{1}{kL}\right]M_A-\frac{1}{EIk}\left[\frac{1}{\sin kL}-\frac{1}{kL}\right]M_B\text{(a)}\\ {\theta }_B=y^{\prime }\left(L\right)=-\frac{1}{EIk}\left[\frac{1}{\sin kL}-\frac{1}{kL}\right]M_A-\frac{1}{EIk}\left[\frac{\cos kL}{\sin kL}-\frac{1}{kL}\right]M_B\text{(b)}\end{array}$$E7

Eq. (7) can be written in matrix from as

$$\left\{\begin{array}{c}{M}_A\\ M_B\end{array}\right\}=\frac{EI}{L}\left[\begin{array}{cc}{S}_1& S_2\\ S_2& S_1\end{array}\right]\left\{\begin{array}{c}{\theta }_A\\ {\theta }_B\end{array}\right\}$$E8

where $S_1$ and $S_2$ are the stability functions defined as

$$\begin{array}{l}{S}_1=\frac{kL\left(\sin kL-kL\cos kL\right)}{2-2\cos kL-kL\sin kL}\text{(a)}\\ S_2=\frac{kL\left(kL-\sin kL\right)}{2-2\cos kL-kL\sin kL}(\text{b)}\end{array}$$E9

$S_1$and $S_2$ account for the coupling effect between axial force and bending moments of the beam-column member. For members subjected to an axial force that is tensile rather than compressive, the stability functions are redefined as

$$\begin{array}{l}{S}_1=\frac{kL\left(kL\cosh kL-\sinh kL\right)}{2-2\cosh kL-kL\sinh kL}(\text{a)}\\ S_2=\frac{kL\left(\sinh kL-kL\right)}{2-2\cosh kL-kL\sinh kL}\text{(b)}\end{array}$$E10

Eqs. (9) and (10) are indeterminate when the axial force is zero (i.e.$$kL=0$$). To overcome this problem, the following simplified equations are used to approximate the stability functions when the axial force in the member falls within the range of -2.0 ≤ $\rho$ ≤ 2.0

$$\begin{array}{l}{S}_1=4+\frac{2{\pi }^2\rho }{15}-\frac{(0.01\rho +0.543){\rho }_{}^2}{4+\rho }-\frac{(0.004\rho +0.285){\rho }_{}^2}{8.183+\rho }\text{(a)}\\ S_2=2-\frac{{\pi }^2\rho }{30}+\frac{(0.01\rho +0.543){\rho }_{}^2}{4+\rho }-\frac{(0.004\rho +0.285){\rho }_{}^2}{8.183+\rho }\text{(b)}\end{array}$$E11

where$$\rho =P/P_e=P/({\pi }^{2}EI/L^2)={\left(kL/\pi \right)}^2$$. For most practical applications, it gives excellent correlation to the "exact" expressions given by Eqs. (9) and (10). However, for $\rho$ other than the range of -2.0 ≤ $\rho$ ≤ 2.0, the conventional stability functions in Eqs. (9) and (10) should be used. The incremental member force and deformation relationship of a three-dimensional beam-column element under axial force and end moments can be written as

$$\left\{\begin{array}{c}\Delta P\\ \Delta M_{yA}\\ \Delta M_{yB}\\ \Delta M_{zA}\\ \Delta M_{zB}\\ \Delta T\end{array}\right\}=\left[\begin{array}{cccccc}\frac{EA}{L}& 0& 0& 0& 0& 0\\ 0& S_{1y}\frac{E{I}_y}{L}& S_{2y}\frac{E{I}_y}{L}& 0& 0& 0\\ 0& S_{2y}\frac{E{I}_y}{L}& S_{1y}\frac{E{I}_y}{L}& 0& 0& 0\\ 0& 0& 0& S_{1z}\frac{E{I}_z}{L}& S_{2z}\frac{E{I}_z}{L}& 0\\ 0& 0& 0& S_{2z}\frac{E{I}_z}{L}& S_{1z}\frac{E{I}_z}{L}& 0\\ 0& 0& 0& 0& 0& \frac{GJ}{L}\end{array}\right]\left\{\begin{array}{c}\Delta \delta \\ \Delta {\theta }_{yA}\\ \Delta {\theta }_{yB}\\ \Delta {\theta }_{zA}\\ \Delta {\theta }_{zB}\\ \Delta \varphi \end{array}\right\}$$E12

where$\Delta P$, $$\Delta M_{yA}$$, $$\Delta M_{yB}$$, $$\Delta M_{zA}$$, $$\Delta M_{zB}$$, and $\Delta T$ are the incremental axial force, end moments with respect to $y$ and $z$ axes, and torsion, respectively;$$\Delta \delta$$, $$\Delta {\theta }_{yA}$$, $$\Delta {\theta }_{yB}$$, $$\Delta {\theta }_{zA}$$, $$\Delta {\theta }_{zB}$$, and $$\Delta \varphi$$ are the incremental axial displacement, the end rotations, and the angle of twist, respectively; $S_{1n}$and $S_{2n}$ are stability functions with respect to $n$ axis $$\left(n=y,z\right)$$ given in Eqs. (9) and (10); and$EA$, $E{I}_n$, and $GJ$ denote the axial, bending, and torsional stiffness, respectively.

2.2. Refined plastic hinge model accounting for inelastic effects

The refined plastic hinge model is an improvement of the elastic plastic hinge one. Two modifications are made to account for a smooth degradation of plastic hinge stiffness: (1) the tangent modulus concept is used to capture the residual stress effect along the length of the member, and (2) the parabolic function is adopted to represent the gradual yielding effect in forming plastic hinges. The inelastic behavior of the member is modeled in terms of member force instead of the detailed level of stresses and strains as used in the plastic zone method. As a result, the refined plastic hinge method retains the simplicity of the elastic plastic hinge method, but it is sufficiently accurate for predicting the strength and stability of a structural system and its component members.

2.2.1. Gradual yielding due to residual stresses

The Column Research Council (CRC) tangent modulus concept is employed to account for the gradual yielding along the member length due to residual stresses. The elastic modulus $E$ (instead of moment of inertia$I$) is reduced to account for the reduction of the elastic portion of the cross-section since the reduction of the elastic modulus is easier to implement than a new moment of inertia for every different section. The rate of reduction in stiffness is different in the weak and strong direction, but this is not considered since the dramatic degradation of weak-axis stiffness is compensated for by the substantial weak-axis plastic strength. This simplification makes the present method more practical. The CRC tangent modulus $E_t$ can be written as

$$\begin{array}{l}\begin{array}{c}{E}_t=1.0E\text{for}P\le 0.5P_y(\text{a})\\ E_t=4\frac{P}{P_y}\left(1-\frac{P}{P_y}\right)E\text{for}0.5P_y<P\le P_y\text{(b)}\\ E_t=0\text{for}P>P_y\text{(c)}\end{array}\\ \end{array}$$E13

Equation (13) is plotted in Fig. 6. The tangent modulus $E_t$ is reduced from the elastic value when$P>0.5P_y$.

2.2.2. Gradual yielding due to flexure

The tangent modulus concept is suitable for the member subjected to axial force, but not adequate for cases of both axial force and bending moment. A gradual stiffness degradation model for a plastic hinge is required to represent the partial plastification effects associated with flexure. The parabolic function is used to represent the smooth transition from elastic stiffness at the onset of yielding to the stiffness associated with a full plastic hinge. The parabolic function $\eta$ representing the gradual stiffness degradation is obtained based on a calibration with plastic zone solutions of simple portal frames and beam-columns. It should be noted that only a simple relationship for $\eta$ is required to describe the degradation in stiffness associated with flexure. Although more complicated expressions for $\eta$ can be proposed, simple expression for $\eta$ is needed for keeping the analysis model simple and straightforward.

Figure 6.

Stiffness reduction due to residual stress

Figure 7.

Stiffness degradation function

The value of parabolic function $\eta$ is equal to 1.0 when the element is elastic, and zero when a plastic hinge is formed. The parabolic function $\eta$ can be expressed as (see Fig. 7.)

$$\begin{array}{c}\eta =1.0\text{for}\alpha \le 0.5\\ \begin{array}{c}\eta =4\alpha \left(1-\alpha \right)\text{for}0.5<\alpha \le 1.0\\ \eta =0\text{for}\alpha >1\end{array}\end{array}\begin{array}{c}(\text{a})\\ \text{(b)}\\ \text{(c)}\end{array}$$E14

where $\alpha$ is the force-state parameter which can be expressed by AISC-LRFD or modified Orbison yield surfaces as (seeFig. 8.).

For AISC-LRFD yield surface (AISC, 2005)

$$\begin{array}{l}\alpha =p+\frac{8}{9}{m}_y+\frac{8}{9}{m}_z\text{for}p\ge \frac{2}{9}{m}_y+\frac{2}{9}{m}_z(\text{a)}\\ \alpha =\frac{p}{2}+m_y+m_z\text{for}p<\frac{2}{9}{m}_y+\frac{2}{9}{m}_z\text{(b)}\end{array}$$E15

For modified Orbison yield surface (McGuire et al., 2000)

$$\alpha =p^2+m_z^2+m_y^4+3.5p^2{m}_z^2+3.0p^6{m}_y^2+4.5m_z^4{m}_y^2$$E16

where$p=P/P_y$, $m_z=M_z/M_{pz}$(strong-axis), $m_y=M_y/M_{py}$(weak-axis); P_y, M_yp, M_zp are axial load, and plastic moment capacity of the cross-section about $y$ and $z$ axes.

Figure 8.

Plastification surface

When the force point moves inside or along the initial yield surface$$\left(\alpha \le 0.5\right)$$, the element remains fully elastic (i.e. no stiffness reduction,$\eta =1.0$). If the force point moves beyond the initial yield surface and inside the full yield surface$$\left(0.5<\alpha \le 1.0\right)$$, the element stiffness is reduced to account for the effect of plastification at the element end. The reduction of element stiffness is assumed to vary according to the parabolic function in the Eq. (15b). When member forces violate the plastic strength surface$$\left(\alpha >1.0\right)$$, the member forces will be scaled down to move the force point return the yield surface based on incremental-iterative scheme.

When the parabolic function for a gradual yielding is active at both ends of an element, the incremental member force and deformation relationship in Eq. (12) is modified as

$$\left\{\begin{array}{c}\Delta P\\ \Delta M_{yA}\\ \Delta M_{yB}\\ \Delta M_{zA}\\ \Delta M_{zB}\\ \Delta T\end{array}\right\}=\left[\begin{array}{cccccc}\frac{E_{t}A}{L}& 0& 0& 0& 0& 0\\ 0& k_{iiy}& k_{ijy}& 0& 0& 0\\ 0& k_{ijy}& k_{jjy}& 0& 0& 0\\ 0& 0& 0& k_{iiz}& k_{ijz}& 0\\ 0& 0& 0& k_{ijz}& k_{jjz}& 0\\ 0& 0& 0& 0& 0& \frac{GJ}{L}\end{array}\right]\left\{\begin{array}{c}\Delta \delta \\ \Delta {\theta }_{yA}\\ \Delta {\theta }_{yB}\\ \Delta {\theta }_{zA}\\ \Delta {\theta }_{zB}\\ \Delta \varphi \end{array}\right\}$$E17

where

$$\begin{array}{l}{k}_{iiy}={\eta }_A(S_1-\frac{S_2^2}{S_1}(1-{\eta }_B))\frac{E_t{I}_y}{L}\text{(a)}\\ k_{ijy}={\eta }_A{\eta }_B{S}_2\frac{E_t{I}_y}{L}(\text{b})\\ k_{jjy}={\eta }_B(S_1-\frac{S_2^2}{S_1}(1-{\eta }_A))\frac{E_t{I}_y}{L}(\text{c)}\\ k_{iiz}={\eta }_A(S_3-\frac{S_4^2}{S_3}(1-{\eta }_B))\frac{E_t{I}_z}{L}\text{(d)}\\ k_{ijz}={\eta }_A{\eta }_B{S}_4\frac{E_t{I}_z}{L}\text{(e)}\\ k_{jjz}={\eta }_B(S_3-\frac{S_4^2}{S_3}(1-{\eta }_A))\frac{E_t{I}_z}{L}\text{(f)}\end{array}$$E18

where ${\eta }_A$ and ${\eta }_B$ are the values of parabolic functions at the ends A and B, respectively.

2.3. Fiber model accounting for inelastic effects

The concept of fiber model is presented in Fig. 9. In this model, the element is divided into a number of monitored sections represented by the integration points. Each section is further divided into $m$ fibers and each fiber is represented by its area $A_i$ and coordinate location corresponding to its centroid$\left(y_i,z_i\right)$. The inelastic effects are captured by tracing the uniaxial stress-strain relationship of each fiber on the cross sections located at the selected integration points along the member length.

The incremental force and deformation relationship, Eq. (12), which accounts for the $P-\delta$ effect can be rewritten in symbolic form as

$$\left\{\Delta F\right\}=\left[K_e\right]\left\{\Delta d\right\}$$E19

where

$$\left\{\Delta F\right\}={\left[\begin{array}{cccccc}\Delta P& \Delta M_{yA}& \Delta M_{yB}& \Delta M_{zA}& \Delta M_{zB}& \Delta T\end{array}\right]}^T$$E20

$$\left\{\Delta d\right\}={\left[\begin{array}{cccccc}\Delta \delta & \Delta {\theta }_{yA}& \Delta {\theta }_{yB}& \Delta {\theta }_{zA}& \Delta {\theta }_{zB}& \Delta \varphi \end{array}\right]}^T$$E21

Figure 9.

Fiber hinge model

$$K_e=\left[\begin{array}{cccccc}\frac{EA}{L}& 0& 0& 0& 0& 0\\ 0& S_{1y}\frac{E{I}_y}{L}& S_{2y}\frac{E{I}_y}{L}& 0& 0& 0\\ 0& S_{2y}\frac{E{I}_y}{L}& S_{1y}\frac{E{I}_y}{L}& 0& 0& 0\\ 0& 0& 0& S_{1z}\frac{E{I}_z}{L}& S_{2z}\frac{E{I}_z}{L}& 0\\ 0& 0& 0& S_{2z}\frac{E{I}_z}{L}& S_{1z}\frac{E{I}_z}{L}& 0\\ 0& 0& 0& 0& 0& \frac{GJ}{L}\end{array}\right]$$E22

in which the axial stiffness$EA$, bending stiffness$E{I}_n$, and torsional stiffness $GJ$ of the fiber element can be obtained as

$$\begin{array}{l}EA={{\displaystyle \sum _{j=1}^h{w}_j\left({\displaystyle \sum _{i=1}^m{E}_i{A}_i}\right)}}_j\text{(a)}\\ E{I}_y={{\displaystyle \sum _{j=1}^h{w}_j\left({\displaystyle \sum _{i=1}^m{E}_i{A}_i{z}_i^2}\right)}}_j(\text{b)}\\ E{I}_z={{\displaystyle \sum _{j=1}^h{w}_j\left({\displaystyle \sum _{i=1}^m{E}_i{A}_i{y}_i^2}\right)}}_j(\text{c})\\ GJ={{\displaystyle \sum _{j=1}^h{w}_j\left[{\displaystyle \sum _{i=1}^m\left(y_i^2+z_i^2\right)G_i{A}_i}\right]}}_j\text{(d)}\end{array}$$E23

in which $h$ is the total number of monitored sections along an element; $m$is the total number of fiber divided on the monitored cross-section; $w_j$is the weighting factor of the $j^{th}$ section; $E_i$and $G_i$ are the tangent and shear modulus of $i^{th}$ fiber, respectively; $$y_i$$and $$z_i$$ are the coordinates of $i^{th}$ fiber in the cross-section. The element stiffness matrix is evaluated numerically by the Gauss-Lobatto integration scheme since this method allows for two integration points to coincide with the end sections of the elements. Since inelastic behavior in beam elements often concentrates at the end of member, the monitoring of the end sections of the element is advantageous from the standpoint of accuracy and numerical stability. By contrast, the outermost integration points of the classical Gauss integration method only approach the end sections with increasing order of integration, but never coincide with the end sections and, hence, result in overestimation of the member strength (Spacone et al., 1996).

Section deformations are represented by three strain resultants: the axial strain $\epsilon$ along the longitudinal axis and two curvatures ${\chi }_z$ and ${\chi }_y$ with respect to $z$ and $y$ axes, respectively. The corresponding force resultants are the axial force $N$ and two bending moments $M_z$ and$M_y$. The section forces and deformations are grouped in the following vectors:

Section force vector$\left\{Q\right\}={\left[\begin{array}{ccc}{M}_z& M_y& N\end{array}\right]}^T$

Section deformation vector$\left\{q\right\}={\left[\begin{array}{ccc}{\chi }_z& {\chi }_y& \epsilon \end{array}\right]}^T$

The incremental section force vector at each integration points is determined based on the incremental element force vector $\left\{\Delta F\right\}$ as

$$\left\{\Delta Q\right\}=\left[B\left(x\right)\right]\left\{\Delta F\right\}$$E24

where $\left[B\left(x\right)\right]$ is the force interpolation function matrix given as

$$\left[B(x)\right]=\left[\begin{array}{cccccc}0& 0& 0& \left(x/L-1\right)& x/L& 0\\ 0& \left(x/L-1\right)& x/L& 0& 0& 0\\ 1& 0& 0& 0& 0& 0\end{array}\right]$$E25

The section deformation vector is determined based on the section force vector as

$$\left\{\Delta q\right\}={\left[k_{\sec }\right]}^{-1}\left\{\Delta Q\right\}$$E26

where $\left[k_{\sec }\right]$ is the section stiffness matrix given as

$$\left[k_{\sec }\right]=\left[\begin{array}{ccc}{\displaystyle \sum _{i=1}^m{E}_i{A}_i{y}_i^2}& {\displaystyle \sum _{i=1}^m{E}_i{A}_i{y}_i{z}_i}& {\displaystyle \sum _{i=1}^m{E}_i{A}_i\left(-y_i\right)}\\ {\displaystyle \sum _{i=1}^m{E}_i{A}_i{y}_i{z}_i}& {\displaystyle \sum _{i=1}^m{E}_i{A}_i{z}_i^2}& {\displaystyle \sum _{i=1}^m{E}_i{A}_i{z}_i}\\ {\displaystyle \sum _{i=1}^m{E}_i{A}_i\left(-y_i\right)}& {\displaystyle \sum _{i=1}^m{E}_i{A}_i{z}_i}& {\displaystyle \sum _{i=1}^m{E}_i{A}_i}\end{array}\right]$$E27

Following the hypothesis that plane sections remain plane and normal to the longitudinal axis, the incremental uniaxial fiber strain vector is computed based on the incremental section deformation vector as

$$\left\{\Delta e\right\}=\left[\Gamma \right]\left\{\Delta q\right\}$$E28

where $\left[\Gamma \right]$ is the linear geometric matrix given as follows

$$\left[\Gamma \right]=\left[\begin{array}{ccc}-y_1& z_1& 1\\ -y_2& z_2& 1\\ \mathrm{...}& \mathrm{...}& \mathrm{...}\\ -y_m& z_m& 1\end{array}\right]$$E29

Once the incremental fiber strain is evaluated, the incremental fiber stress is computed based on the stress-strain relationship of material model. The tangent modulus of each fiber is updated from the incremental fiber stress and incremental fiber strain as

$$E_i=\frac{\Delta {\sigma }_i}{\Delta e_i}$$E30

Eq. (32) leads to updating of the element stiffness matrix $\left[K_e\right]$ in Eq. (22) and section stiffness matrix $\left[k_{\sec }\right]$ in Eq. (29) during the iteration process. Based on the new tangent modulus of Eq. (32), the location of the section centroid is also updated during the incremental load steps to take into account the distribution of section plasticity. The section resisting forces are computed by summation of the axial force and biaxial bending moment contributions of all fibers as

$$\left\{Q_R\right\}=\left\{\begin{array}{c}{M}_z\\ M_y\\ N\end{array}\right\}=\left\{\begin{array}{c}{\displaystyle \sum _{i=1}^m{\sigma }_i}{A}_i\left(-y_i\right)\\ {\displaystyle \sum _{i=1}^m{\sigma }_i}{A}_i{z}_i\\ {\displaystyle \sum _{i=1}^m{\sigma }_i}{A}_i\end{array}\right\}$$E31

2.4. Shear deformation effect

To account for transverse shear deformation effect in a beam-column element, the member force and deformation relationship of beam-column element in Eq. (12) should be modified. The flexibility matrix can be obtained by inversing the flexural stiffness matrix as

$$\left\{\begin{array}{c}\Delta {\theta }_{MA}\\ \Delta {\theta }_{MB}\end{array}\right\}=\left[\begin{array}{cc}\frac{k_{jj}}{k_{ii}{k}_{jj}-k_{ij}^2}& \frac{-k_{ij}}{k_{ii}{k}_{jj}-k_{ij}^2}\\ \frac{-k_{ij}}{k_{ii}{k}_{jj}-k_{ij}^2}& \frac{k_{ii}}{k_{ii}{k}_{jj}-k_{ij}^2}\end{array}\right]\left\{\begin{array}{c}\Delta M_A\\ \Delta M_B\end{array}\right\}$$E32

where $\Delta {\theta }_{MA}$ and $\Delta {\theta }_{MB}$ are the slope of the neutral axis due to bending moment. The flexibility matrix corresponding to shear deformation can be written as

$$\left\{\begin{array}{c}\Delta {\theta }_{SA}\\ \Delta {\theta }_{SB}\end{array}\right\}=\left[\begin{array}{cc}\frac{1}{G{A}_{S}L}& \frac{1}{G{A}_{S}L}\\ \frac{1}{G{A}_{S}L}& \frac{1}{G{A}_{S}L}\end{array}\right]\left\{\begin{array}{c}\Delta M_A\\ \Delta M_B\end{array}\right\}$$E33

where $G{A}_S$ and $L$ are shear stiffness and length of the element, respectively. The total rotations at the two ends $A$ and $B$ are obtained by combining Eqs. (34) and (35) as

$$\left\{\begin{array}{c}\Delta {\theta }_A\\ \Delta {\theta }_B\end{array}\right\}=\left\{\begin{array}{c}\Delta {\theta }_{MA}\\ \Delta {\theta }_{MB}\end{array}\right\}+\left\{\begin{array}{c}\Delta {\theta }_{SA}\\ \Delta {\theta }_{SB}\end{array}\right\}$$E34

The basic force and deformation relationship including shear deformation is derived by inverting the flexibility matrix as

$$\left\{\begin{array}{c}\Delta M_A\\ \Delta M_B\end{array}\right\}=\left[\begin{array}{cc}\frac{k_{ii}{k}_{jj}-k_{ij}^2+k_{ii}{A}_{s}GL}{k_{ii}+k_{jj}+2k_{ij}+A_{s}GL}& \frac{-k_{ii}{k}_{jj}+k_{ij}^2+k_{ij}{A}_{s}GL}{k_{ii}+k_{jj}+2k_{ij}+A_{s}GL}\\ \frac{-k_{ii}{k}_{jj}+k_{ij}^2+k_{ij}{A}_{s}GL}{k_{ii}+k_{jj}+2k_{ij}+A_{s}GL}& \frac{k_{ii}{k}_{jj}-k_{ij}^2+k_{jj}{A}_{s}GL}{k_{ii}+k_{jj}+2k_{ij}+A_{s}GL}\end{array}\right]\left\{\begin{array}{c}\Delta {\theta }_A\\ \Delta {\theta }_B\end{array}\right\}$$E35

The member force and deformation relationship can be extended for three-dimensional beam-column element as

$$\left\{\begin{array}{c}\Delta P\\ \Delta M_{yA}\\ \Delta M_{yB}\\ \Delta M_{zA}\\ \Delta M_{zB}\\ \Delta T\end{array}\right\}=\left[\begin{array}{cccccc}\frac{EA}{L}& 0& 0& 0& 0& 0\\ 0& C_{iiy}& C_{ijy}& 0& 0& 0\\ 0& C_{ijy}& C_{jjy}& 0& 0& 0\\ 0& 0& 0& C_{iiz}& C_{ijz}& 0\\ 0& 0& 0& C_{ijz}& C_{jjz}& 0\\ 0& 0& 0& 0& 0& \frac{GJ}{L}\end{array}\right]\left\{\begin{array}{c}\Delta \delta \\ \Delta {\theta }_{yA}\\ \Delta {\theta }_{yB}\\ \Delta {\theta }_{zA}\\ \Delta {\theta }_{zB}\\ \Delta \varphi \end{array}\right\}$$E36

in which

$$\begin{array}{l}{C}_{iiy}=\frac{k_{iiy}{k}_{jjy}-k_{ijy}^2+k_{iiy}{A}_{sz}GL}{k_{iiy}+k_{jjy}+2k_{ijy}+A_{sz}GL}\text{(a)}\\ C_{ijy}=\frac{-k_{iiy}{k}_{jjy}+k_{ijy}^2+k_{ijy}{A}_{sz}GL}{k_{iiy}+k_{jjy}+2k_{ijy}+A_{sz}GL}\text{(b)}\\ C_{jjy}=\frac{k_{iiy}{k}_{jjy}-k_{ijy}^2+k_{jjy}{A}_{sz}GL}{k_{iiy}+k_{jjy}+2k_{ijy}+A_{sz}GL}\text{(c)}\\ C_{iiz}=\frac{k_{iiz}{k}_{jjz}-k_{ijz}^2+k_{iiz}{A}_{sy}GL}{k_{iiz}+k_{jjz}+2k_{ijz}+A_{sy}GL}\text{(d)}\\ C_{ijz}=\frac{-k_{iiz}{k}_{jjz}+k_{ijz}^2+k_{ijz}{A}_{sy}GL}{k_{iiz}+k_{jjz}+2k_{ijz}+A_{sy}GL}\text{(e)}\\ C_{jjz}=\frac{k_{iiz}{k}_{jjz}-k_{ijz}^2+k_{jjz}{A}_{sy}GL}{k_{iiz}+k_{jjz}+2k_{ijz}+A_{sy}GL}\text{(f)}\end{array}$$E37

where $A_{sy}$ and $A_{sz}$ are the shear areas with respect to $y$ and $z$ axes, respectively.

2.5. Element stiffness matrix

The incremental end forces and displacements used in Eq. (38) are shown in Fig. 10(a). The sign convention for the positive directions of element end forces and displacements of a frame member is shown in Fig. 10(b). By comparing the two figures, the equilibrium and kinematic relationships can be expressed in symbolic form as

$$\begin{array}{l}\left\{f_n\right\}={\left[T\right]}_{6\times 12}^T\left\{F\right\}\text{(a)}\\ \left\{d\right\}={\left[T\right]}_{6\times 12}\left\{d_L\right\}\text{(b)}\end{array}$$E38

where $\left\{f_n\right\}$ and $\left\{d_L\right\}$ are the nodal force and nodal displacement vectors of the element expressed as

$$\begin{array}{l}{\left\{}^{f_n}=\left\{\begin{array}{cccccccccccc}{r}_{n1}& r_{n2}& r_3& r_4& r_5& r_6& r_7& r_8& r_9& r_{10}& r_{11}& r_{12}\end{array}\right\}\text{(a)}\\ {\left\{}^{d_L}=\left\{\begin{array}{cccccccccccc}{d}_1& d_2& d_3& d_4& d_5& d_6& d_7& d_8& d_9& d_{10}& d_{11}& d_{12}\end{array}\right\}\text{(b)}\end{array}$$E39

and $$\left\{F\right\}$$ and $$\left\{d\right\}$$ are the basic member force and displacement vectors given in Eqs. (20) and (21), respectively. ${\left[T\right]}_{6\times 12}$is a transformation matrix written as

$${\left[T\right]}_{6\times 12}=\left[\begin{array}{cccccccccccc}-1& 0& 0& 0& 0& 0& 1& 0& 0& 0& 0& 0\\ 0& 0& -1/L& 0& 1& 0& 0& 0& 1/L& 0& 0& 0\\ 0& 0& -1/L& 0& 0& 0& 0& 0& 1/L& 0& 1& 0\\ 0& 1/L& 0& 0& 0& 1& 0& -1/L& 0& 0& 0& 0\\ 0& 1/L& 0& 0& 0& 0& 0& -1/L& 0& 0& 0& 1\\ 0& 0& 0& 1& 0& 0& 0& 0& 0& -1& 0& 0\end{array}\right]$$E40

Figure 10.

Force and displacement notations

Using the transformation matrix, the nodal force and nodal displacement relationship of element may be written as

$$\left\{f_n\right\}=\left[K_n\right]\left\{d_L\right\}$$E41

where $\left[K_n\right]$ is the element stiffness matrix expressed as

$${[K_n]}_{12\times 12}={[T]}_{6\times 12}^T{[K_e]}_{6\times 6}{[T]}_{6\times 12}$$E42

It should be noted that Eq. (43) is used for the beam-column member in which side-sway is restricted. If the beam-column member is permitted to sway, additional axial and shear forces will be induced in the member. These additional axial and shear forces due to member sway to the member end displacements can be related as

$$\left\{f_s\right\}=\left[K_s\right]\left\{d_L\right\}$$E43

where $\left[K_s\right]$ is the element stiffness matrix due to member sway expressed as

$${\left[K_s\right]}_{12\times 12}=\left[\begin{array}{cc}\left[G_s\right]& -\left[G_s\right]\\ -{\left[G_s\right]}^T& \left[G_s\right]\end{array}\right]$$E44

in which

$$\left[G_s\right]=\left[\begin{array}{cccccc}0& \left(M_{zA}+M_{zB}\right)/L^2& \left(M_{yA}+M_{yB}\right)/L^2& 0& 0& 0\\ \left(M_{zA}+M_{zB}\right)/L^2& P/L& 0& 0& 0& 0\\ \left(M_{yA}+M_{yB}\right)/L^2& 0& P/L& 0& 0& 0\\ 0& 0& 0& 0& 0& 0\\ 0& 0& 0& 0& 0& 0\\ 0& 0& 0& 0& 0& 0\end{array}\right]$$E45

By combining Eqs. (43) and (47), the general force-displacement relationship of beam-column element obtained as

$$\left\{f_L\right\}=\left[K\right]\left\{d_L\right\}$$E46

where

$$\left\{f_L\right\}=\left\{f_n\right\}+\left\{f_s\right\}$$E47

$$\left[K\right]=\left[K_n\right]+\left[K_s\right]$$E48

2.6. Solution algorithm

The generalized displacement control method proposed by Yang and Shieh (1990) appears to be one of the most robust and effective method because of its general numerical stability and efficiency. This method is adopted herein to solve the nonlinear equilibrium equations. The incremental form of the equilibrium equation can be rewritten for the

th iteration of the

th incremental step as

$$\left[K_{j-1}^i\right]\left\{\Delta D_j^i\right\}={\lambda }_j^i\left\{\widehat{P}\right\}+\left\{R_{j-1}^i\right\}$$E49

where $$\left[K_{j-1}^i\right]$$ is the tangent stiffness matrix, $$\left\{\Delta D_j^i\right\}$$is the displacement increment vector, $$\left\{\widehat{P}\right\}$$is the reference load vector, $$\left\{R_{j-1}^i\right\}$$is the unbalanced force vector, and $${\lambda }_j^i$$ is the load increment parameter. According to Batoz and Dhatt (1979), Eq. (51) can be decomposed into the following equations:

$$\left[K_{j-1}^i\right]\left\{\Delta {\widehat{D}}_j^i\right\}=\left\{\widehat{P}\right\}$$E50

$$\left[K_{j-1}^i\right]\left\{\Delta {\overline{D}}_j^i\right\}=\left\{R_{j-1}^i\right\}$$E51

$$\left\{\Delta D_j^i\right\}={\lambda }_j^i\left\{\Delta {\widehat{D}}_j^i\right\}+\left\{\Delta {\overline{D}}_j^i\right\}$$E52

Once the displacement increment vector $$\left\{\Delta D_j^i\right\}$$ is determined, the total displacement vector $$\left\{D_j^i\right\}$$ of the structure at the end of $j$th iteration can be accumulated as

$$\left\{D_j^i\right\}=\left\{D_{j-1}^i\right\}+\left\{\Delta D_j^i\right\}$$E53

The total applied load vector $$\left\{P_j^i\right\}$$ at the $j$th iteration of the $i$th incremental step relates to the reference load vector $$\left\{\widehat{P}\right\}$$ as

$$\left\{P_j^i\right\}={\Lambda }_j^i\left\{\widehat{P}\right\}$$E54

where the load factor $${\Lambda }_j^i$$ can be related to the load increment parameter $${\lambda }_j^i$$ by

$${\Lambda }_j^i={\Lambda }_{j-1}^i+{\lambda }_j^i$$E55

The load increment parameter $${\lambda }_j^i$$ is an unknown. It is determined from a constraint condition. For the first iterative step$\left(j=1\right)$, the load increment parameter $${\lambda }_j^i$$ is determined based on the generalized stiffness parameter $\left(GSP\right)$ as

$${\lambda }_1^i={\lambda }_1^1\sqrt{\left|GSP\right|}$$E56

where $${\lambda }_1^1$$ is an initial value of load increment parameter, and the $$GSP$$ is defined as

$$GSP=\frac{{\left\{\Delta {\widehat{D}}_1^1\right\}}^T\left\{\Delta {\widehat{D}}_1^1\right\}}{{\left\{\Delta {\widehat{D}}_1^{i-1}\right\}}^T\left\{\Delta {\widehat{D}}_1^i\right\}}$$E57

For the iterative step$\left(j\ge 2\right)$, the load increment parameter $${\lambda }_j^i$$ is calculated as

$${\lambda }_j^i=-\frac{{\left\{\Delta {\widehat{D}}_1^{i-1}\right\}}^T\left\{\Delta {\overline{D}}_j^i\right\}}{{\left\{\Delta {\widehat{D}}_1^{i-1}\right\}}^T\left\{\Delta {\widehat{D}}_j^i\right\}}$$E58

where $$\left\{\Delta {\widehat{D}}_1^{i-1}\right\}$$ is the displacement increment generated by the reference load at the first iteration of the previous incremental step; and $$\left\{\Delta {\widehat{D}}_j^i\right\}$$and $$\left\{\Delta {\overline{D}}_j^i\right\}$$ denote the displacement increments generated by the reference load and unbalanced force vectors, respectively, at the $j$th iteration of the $i$th incremental step, as defined in Eqs. (52) and (53).

3. Numerical examples

In this section, three numerical examples are presented to verify the accuracy and efficiency of two proposed analysis methods: (1) the refined plastic hinge method and (2) the fiber method. The predictions of strength and load-displacement relationship are compared with those generated by commercial finite element packages and other existing solutions. The first example is to show how the stability functions capture the $P-\delta$ effect accurately and efficiently. The second one is to show how well the refined plastic hinge model and fiber hinge model predict the strength and behavior of frames. The last one is to demonstrate the capability of two proposed methods in predicting the strength and behavior of a large-scale twenty-story space frame. Five integration points along the length of a member and eighty fibers on the cross-section are used in the fiber model.

3.1. Elastic buckling of columns

The aim of this example is to show the accuracy and efficiency of the stability functions in capturing the elastic buckling loads of columns with different boundary conditions. Fig. 11 shows cantilever and simply supported columns. The section of columns is W8×31. The Young’s modulus and Poisson ratio of the material are $E=200,000$ MPa and$\nu =0.3$, respectively. The buckling load of the columns is obtained using the load-deflection analysis. The geometric imperfection is modeled by equivalent notional lateral loads as shown in Fig. 11.

Fig. 12 shows the load-displacement curves of the columns predicted by the present element and the cubic frame element of SAP2000. Since the present element is based on the stability functions which are derived from the closed-form solution of a beam-column subjected to end forces, it can accurately predict the buckling load of columns with different boundary conditions by using only one element per member. Whereas the cubic frame element of SAP2000, which is based on the cubic interpolation functions, overpredicts the buckling loads by 18% and 16% for the cantilever column and simply supported column, respectively, when the columns are modeled by one element per member. The load-displacement curves shown in Fig. 12 indicate that SAP2000 requires more than five cubic elements per member in modeling to match the results predicted by the present element. This is due to the fact that when the member is divided into many elements, the $P-\delta$ effect is transformed to the $P-\Delta$ effect, and hence, the results of cubic element are close to the obtained results.

Figure 11.

Steel columns

3.2. Two story space frame

A two-story space subjected to combined action of gravity load and lateral load is depicted in Fig. 13 with its geometric dimension. The Young modulus, Poisson ratio, and yield stress of material are $E=19,613$ MPa, $\nu =0.3$, and ${\sigma }_y=98$ MPa, respectively. This frame was previously analyzed by De Souza (2000) using the force-based method with fiber model. De Souza (2000) used one element per member in the modeling. The B23 element of ABAQUS is also employed to model this frame. Each framed member is modeled by one present element. The aim of this example is to demonstrate capability of the present element in capturing the effects of both geometric and material nonlinearities.

The ultimate loads of the frame obtained by different methods are presented in Table 1. The load-displacement responses of the frame are also plotted in Fig. 14. It can be seen that the results of the present element are well compared with those of De Souza (2000) using the force-based method. It should be noted that only one element per member is used in present study and De Souza (2000). The B23 element of ABAQUS overestimates ultimate strength of this frame if each framed member is modeled by less than fifty B23 elements. The difference between B23 element and present element is negligible when more than fifty B32 elements are used, and the ultimate strength and load-displacement curve obtained by ABAQUS and present study are then close each other.

Figure 12.

Load-displacement curves of steel columns

Figure 13.

Two-story space frame

Figure 14.

Load- displacement curves of two-story space frame

Method	Ultimate load (kN)	Difference (%)
De Souza ADDIN EN.CITE <EndNote"/><Cite ExcludeAuth="1""/><Author"/>De Souza</Author"/><Year"/>2000</Year"/><RecNum"/>250</RecNum"/><record"/><rec-number"/>250</rec-number"/><foreign-keys"/><key app="EN" db-id="pxx2ew5ez0dxsne50ahvs0ena0ved5v92v5w""/>250</key"/></foreign-keys"/><ref-type name="Thesis""/>32</ref-type"/><contributors"/><authors"/><author"/>De Souza, RM</author"/></authors"/></contributors"/><titles"/><title"/>Force-based finite element for large displacement inelastic analysis of frames.</title"/></titles"/><volume"/>PhD Dissertation</volume"/><dates"/><year"/>2000</year"/></dates"/><publisher"/>Department of Civil and Environmental Engineering, University of California at Berkeley</publisher"/><work-type"/>PhD Dissertation</work-type"/><urls"/></urls"/></record"/></Cite"/></EndNote"/>(2000)	128.05	-
ABAQUS (5 element/member)	140.26	9.53
ABAQUS (20 element/member)	132.19	3.23
ABAQUS (50 element/member)	130.74	2.10
Present (refined plastic hinge model)	128.50	0.35
Present (fiber model)	128.82	0.60

Table 1.

Comparison of ultimate load of two-story space frame

3.3. Twenty-story space frame

The last example is a large scale twenty-story space steel frame as shown in Fig. 15. The aim of this example is to demonstrate the capability of two proposed methods in predicting the strength and behavior of large-scale structures. A50 steel with yield stress of 344.8 Mpa, Young’s modulus of 200 Gpa, and Poisson’s ratio of 0.3 is used for all sections. The load applied to the structure consists of gravity loads of 4.8 kN/m² and wind loads of 0.96 kN/m² acting in the Y-direction. These loads are converted into concentrated loads applied at the beam-column joints. The obtained results are also compared with those generated by Jiang et al. (2002) using the mixed element method.

Jiang et al. (2002) used both the plastic hinge and spread-of-plasticity elements to model this structure to shorten the computational time because the use of a full spread-of-plasticity analysis is very computationally intensive. When a member modeling by one plastic hinge element detected yielding to occur between the two ends, it was divided into eight spread-of-plasticity elements to accurately capture the inelastic behavior. In this study, each framed member is modeled by only one proposed element. The load-displacement curves of node A at the roof of the frame obtained by the present elements and mixed element of Jiang et al. (2002) are shown in Fig. 16. The ultimate load factor of the frame is also given in Table 2. A very good agreement between the results is seen.

Method	Ultimate load factor	Difference (%)
Jiang et al. (2002)	1.000	-
Present (refined plastic hinge model)	1.021	2.10
Present (fiber model)	1.0002	0.02

Table 2.

Analysis result of twenty-story space frame

Figure 15.

Twenty-story space frame

Figure 16.

Load-displacement curves of twenty-story space frame

4. Conclusion

This chapter has presented two advanced analysis methods for space steel frames. In these methods, the geometric nonlinearities are captured using the stability functions, while the material nonlinearities are considered using the refined plastic hinge model and fiber model. The benefit of using the stability functions is that they require only one element per member, and hence, minimize the modeling and solution time. The advantage of refined plastic hinge model is its simplicity and efficiency. However, it is limited to steel material. Although the fiber model is a little bit time consuming compared to the refined plastic hinge model, it can be used for both steel and concrete or concrete-filled steel tubular structures as shown in the works of Thai & Kim (2011a).