Abstract
Buckling of nonprismatic single columns with arbitrary boundary conditions resting on a nonuniform elastic foundation may be considered as the most generalized treatment of the subject. The buckling differential equation for such columns is extremely difficult to solve analytically. Thus, the authors propose a numerical approach by discretizing the column into a finite number of segments. Each segment has constants (modulus of elasticity), (moment of inertia), and (subgrade stiffness). Next, an exact analytical solution is derived for each prismatic segment resting on uniform elastic foundation. These segments are then assembled in a matrix from which the critical buckling load is obtained. The derived formulation accounts for different end boundary conditions. Validation is performed by benchmarking the present results against analytical solutions found in the literature, showing excellent agreement. After validation, more examples are solved to illustrate the power and flexibility of the proposed method. Overall, the proposed method provides reasonable results, and the examples solved demonstrate the versatility of the developed approach and some of its many possible applications.
1. Introduction
By the beginning of the twentieth century, many researchers got into studying the specific cases of bucking using various methods such as continuous and lumped matrix analysis, Finite Element Method (FEM), and Boundary Element Method (BEM) [1]. Dinnik reported the exact solution of simply supported columns with monomial variation in stiffness and axial load in his paper published in 1932 [2]. Simply supported tapered columns were studied and the analytical solutions were presented by Gere and Carter in 1962 [3]. In 1970, Gallagher and Lee solved an axially loaded column with a variable flexural stiffness using Finite Element Method [2]. Additionally, Elishakoff and Bert in 1988 succeeded at improving Rayleigh’s method and obtaining an approximate solution for columns with variable stiffness [2]. A year later, Eisenberger presented the exact buckling solution for columns with variable cross-sections and variable axial load using any polynomial variation and using different boundary conditions [2]. In 1986, Ermopoulos published the solution for buckling of tapered bars under stepped axial loads [4]. Analysis on the buckling of uniform columns with different supports was performed by Iyengar in his book which was published in 1988 [5]. Elastic stability of elastically supported columns under the effect of distributed forces was studied by Lee and Kuo in 1991 [5]. In the same year, Arbabei and Li presented the solution for the buckling of nonprismatic elastic columns [6]. A more recent solution for the exact buckling of constrained stepped columns was presented by Yang and Park in 2003 [7]. Coşkun and Atay analyzed the critical buckling load for elastic columns of constant and variable cross-sections using variational iteration method in 2009 [6]. This variational iteration method produces an approximate solution for the presented problem. In 2012, Huang and Li published their solution for analytically determining the exact critical buckling loads of nonuniform columns [8]. Lee et al. published their calculations of natural frequencies and buckling loads for columns with intermediate multiple elastic springs in 2002 [9]. Atanackovic and Novakovic worked on Lagrange problem, in which they determined the optimal shape of an elastic column on elastic foundation using an area varying equation to represent the column shape [10]. Using Pontryagin’s maximum principle, Levy presented the optimal shape of simply supported columns on elastic foundations in 1990 [10]. Using Galerkin’s method, Lacarbonara solved for buckling and postbuckling of nonuniform nonlinearly elastic rods [11]. To the authors’ best knowledge, this is the first paper that provides a numerical treatment combining the buckling of nonprismatic columns resting on nonuniform elastic foundation with arbitrary end boundary conditions.
2. Formulation
2.1. The Buckling Equation
The buckling differential equation for the column shown in Figure 1 can be derived using the bifurcation method resulting in the following equation:
The boundary conditions (BCs) are as follows:
(1) Shear
(2) Momentwhere is the column elastic modulus function, is the column moment of inertia function, is the column axial load, is the elastic foundation spring stiffness function, is the column lateral deflection function, is the distance from the bottom of the column, is the bottom lateral discrete spring stiffness, is the top lateral discrete spring stiffness, is the bottom rotational discrete spring stiffness, is the top rotational discrete spring stiffness, is the length of the column, is the column shear force function, and is the column moment function.
Such a differential equation can prove to be extremely difficult to solve analytically for arbitrary , , and functions. The authors are not aware of any analytical mathematical method that can solve this differential equation. Thus, the authors decided to solve the differential equation semianalytically by discretizing the column into a finite number of segments. Each segment has constants , , and , as shown in Figure 2. , , and values for each segment were taken to be the average of each function within that segment of the column. The averages were calculated using the general equation for averaging functions within an interval:
An exact solution was then obtained for each segment using the classical buckling equation of columns on an elastic foundation. The domain of for each segment was taken to be :where is number of segments, is segment length, and is segment number.
2.2. Analytical Solution of the Segment Buckling Equation
The general solution to the segment differential equation with arbitrary boundary conditions is complicated and requires the consideration of four different cases. The solution to each case is shown below.
Case 1 (). The segment differential equation becomesLet , which can be substituted into (6) resulting inwhereThe deflection function is then obtained as follows:
Case 2 (). Let , which can be substituted into (5) resulting inThe deflection function is then found to be
Case 3 (). Let , which can be substituted into (5) resulting inwhereThe deflection function is then obtained to be as follows
Case 4 (). Let , which can be substituted into (5) resulting in
Let implying which results in the following two equations with two unknowns, :
From (17) we get
Substituting into (16) we obtain
We know that , but it is obvious that . Therefore, since . This leaves the following: .
Using (17) we can then obtain as follows:
Therefore, is found to bewhere
The deflection function is then obtained to bewhere are the constant coefficients for each segment.
After solving the segment differential equation, the following compatibility conditions were established at each node between the segments.
(1) Displacement
(2) Slope
(3) Moment
(4) Shear
Finally, in order to obtain the coefficient matrix (P-matrix), the following boundary conditions were established.
(1) Shear
(2) Moment
Solving the following system of equations would then render the buckling load:where is the -Matrix and = .
2.3. The -Matrix
The structure of the -matrix is as follows.
(1) One Segment
(2) Two Segments
(3) Three Segments
For the case of an assembly with more than three segments, the structure of the -matrix would be the same as the three segments case but with additional intermediate segment compatibility submatrices that are shifted diagonally downwards. The buckling load would then be the load that produces a singular -matrix.
The submatrices inside the -matrix are obtained by applying the compatibility and boundary conditions and are dependent on the case that governs each segment, as explained earlier. These submatrices were found to be as shown below.
(a) Intermediate Segments
Case 1 ().
Case 2 ().
Case 3 ().
Case 4 ().
(b) First Segment
Case 1 ().
Case 2 ().
Case 3 ().
Case 4 ().
(c) Last Segment
Case 1 ().
Case 2 ().
Case 3 ().
Case 4 ().
2.4. The Mode Shape
The mode shape can be obtained by substituting the buckling load into the -matrix and setting one of the constant coefficients to be equal to 1. The original system looks like the following:
After letting , the following system of equations is obtained.
The constant coefficients can then be calculated by solving the system as follows:where is the reduced -matrix, , and .
After substituting the constants back into their respective equations, one can obtain the deflection functions and plot the mode shape.
3. Sensitivity Analysis
As the proposed approach solves the differential equation for buckling numerically by discretizing the column into a finite number of segments, a sensitivity analysis was conducted to determine the optimum number of segments necessary to obtain accurate results. The sensitivity analysis was conducted for a simply supported column with both cross-section and soil stiffness varying to account for all approximations involved in the solution. The following functions were used: and . The first 3 buckling loads and their associated mode shapes were determined. Table 1 shows the results for 10, 20, 50, 100, and 200 segments.
The results, for the discretization segments above, were reasonably close; however, the results required a high segment count to stabilize. In the next section, it will be shown that this variation is not significant as the obtained results with a segment count of 20 provided excellent agreement with the results in the literature. Figure 3 shows the first three mode shapes, which did not vary with the change in the number of segments.
4. Validation of the Proposed Approach
In order to validate the proposed approach, a few cases were analyzed and compared with available solutions from the literature. For the first case, a column resting on elastic foundation is analyzed. A fixed-pinned column was evaluated with the following parameters: , , and [10]. In this case, a single element was used in the analysis as no variation in any of the parameters occurs along the span. The critical load was determined to be 28.307, which matches the load obtained by Wang and his coauthors [12]. Figure 4 shows the first three mode shapes for this case and their associated buckling loads.
The following cases are those for a nonprismatic column, as derived by Gere and Carter [13]. For Case 2, nonprismatic fixed-pinned column, as derived by Gere and Carter [13], was analyzed. The parameters were as follows: , , and , and the number of elements was taken to be 20. The buckling load obtained was 221.031, compared to 222.2 obtained by Gere’s approach. Figure 5 shows the mode shape obtained for this case.
In Case 3, the variable end conditions implementation was validated. This was accomplished by analyzing cases with different rotational spring values at the ends and comparing the results with the alignment charts. This was done for both the sway permitted and sway not permitted cases. The column analyzed had the following parameters: . For the pinned-pinned reference case, the Euler buckling load was obtained to be 9.8696. The results obtained are shown in Table 2.
The results obtained showed excellent agreement with those obtained from the alignment charts. This validates the derivation of variable end supports and its implementation in the computer program.
Finally, in Cases 4 and 5, the approach was validated for cases with nonprismatic sections [2] and variable soil stiffness [14], respectively. In Case 4, columns with different boundary conditions were analyzed as shown in Table 3. The parameters for these columns were as follows: , , and . The results obtained are also shown in Table 3.
For Case 5, a column resting on a variable elastic foundation was analyzed. The column had the following parameters: N·m2, m, and . The critical load was obtained to be 148514.92 N, which is very close to what was obtained by Eisenberger and Clastornik [14].
5. Examples
5.1. Example 1: Column on Variable (Linear and Nonlinear) Elastic Foundation
After validating the proposed formulation in the previous section, the case of nonuniform soil pressure was investigated. A simply supported column was analyzed with 20 segments. The column’s parameters were as follows: and , and, similar to what Gere and Carter [13] did for the moment of inertia in nonprismatic columns, the soil stiffness was taken to be as follows:where
The soil pressure profiles are shown in Figure 6.
The first three buckling loads were determined for each case as shown in Table 4. It is observed that the buckling loads experienced an increase with the increase in the exponent . This is reasonable because, as shown in Figure 6, the area under the soil pressure distribution increases as power increases.
The mode shapes are shown in Figure 7. The first mode shape shifts to the left as the exponent increases, which is expected as the relative soil stiffness is less in that area.
(a)
(b)
(c)
5.2. Example 2: Crack Propagation
In this example, a simply supported column was analyzed. As the maximum moment in simply supported structural elements develops at the midspan, it was assumed that crack initiation occurred at the center of the column. The crack was modeled by a change in the cross-section of column, making it nonprismatic, as shown in Figure 8. This change was attained by reducing the moment of inertia . The propagation of the cracks in the column was modeled by increasing the cracked length . The following constants were used for the columns in this example: , , , and . The matrix of the columns to be modeled in this example was generated by varying the crack length to column length ratio and moment of inertia of the cracked section to gross moment of inertia ratio. The varying parameters from which the combinations were generated are shown in Table 5.
The analysis was performed, and the results for the columns are shown in Table 6. For each moment of inertia ratio, the first column shows the buckling load obtained, while the second column shows the buckling load normalized against Euler’s buckling load. Euler’s buckling load is obtained for the uncracked prismatic column, which corresponds to a crack length ratio of 0 (row 1). The compression action in the column forces the crack to close; however, this causes an imperfection at that location which encourages the initiation of buckling from it.
In order to further analyze the obtained results, Figures 9 and 10 were plotted. Figure 9 plots the normalized buckling loads against the cracked span ratios. A reduction in the buckling load as the cracked span ratio increased was observed. This is logical as the propagation of the cracks reduces the overall moment of inertia of the column and thus decreases the load required to buckle the column. Each series of data points was fitted with a second-degree polynomial function. These regression equations provided excellent correlation .
Figure 10 plots the normalized buckling loads against the cracked moment of inertia ratios. It was noted that the buckling load decreased as the cracked moment of inertia ratio decreased. This is rational because the smaller the ratio is, the larger the crack is and, thus, the smaller the overall moment of inertia becomes. An interesting trend observed was that the curve was bilinear for small cracked span ratios and approached the standard linear behavior as the span ratio increased.
5.3. Example 3: Column on Linearly Variable Elastic Foundation
In this example, a simply supported column was analyzed under different distributions of linear soil pressure. Details on the analyzed soil stiffness functions are shown in Table 7 and Figure 11. These functions provide a constant area under the soil stiffness profile and only differ in the distribution of that area. For this example the length was taken to be 10.
In order to facilitate the analysis, a few parameters will be defined. The slope of the soil stiffness profile is calculated by taking the derivative as follows:
Another parameter is Kappa , which relates the uniform soil stiffness to the stiffness of the column according to the following equation:
For this example, columns with values that correspond to the following ratios were analyzed: 0.1, 1, 10, 100, and 200. The modulus of elasticity and moment of inertia were varied based on this ratio. The obtained buckling loads were normalized against Euler’s buckling load and are shown in Table 8.
Additionally, the data points were plotted in Figures 12 and 13. Figure 12 plots the normalized buckling load versus . It is observed that as the value of κ increased, the normalized buckling load increased. Also, at low values, which indicates that the EI is relatively high compared to , the soil stiffness distribution did not have a significant effect on the buckling load. On the contrary, columns resting on relatively stiff soil experienced a significant variation in the buckling load based on the soil stiffness distribution (about 51%). Data in Table 8 show that as the slope increases, the buckling load decreases. Figure 13 plots the normalized buckling load versus . Lower results are represented by straight lines, which indicate again that the distribution of soil stiffness did not significantly affect the buckling load for these columns. As increased, the decrease in the critical load became more evident with the increase in .
6. Conclusion
In this paper, an approach to determine buckling loads and their associated mode shapes semianalytically for nonprismatic columns with variable boundary conditions resting on nonuniform elastic foundation was formulated. An exact solution for the segment buckling equation was derived under all possible cases that govern the differential equation. Coupling the solutions of the various segments resulted in a system of equations that can be used to numerically solve for the buckling load and plot the mode shape. Sensitivity analysis and benchmarking verification of the derived solution is performed next. Finally, three examples addressing the comparison of linear and nonlinear elastic foundation, crack propagation using nonprismatic column analysis, and investigating the effect of changing soil stiffness distribution on the buckling load are examined.
Conflicts of Interest
The authors declare that there are no conflicts of interest regarding the publication of this paper.
Acknowledgments
Publication of this article was funded in part by the Kansas State University Open Access Publishing Fund.