Investigation into the nonlinear structural behavior of tapered axially functionally graded material beams utilizing absolute nodal coordinate formulations

Abdur Rahman Shaukat; Bing Yan; Feilong Nie; Tengfei Wang; Jia Wang

doi:10.1515/secm-2025-0063

Article Open Access

Investigation into the nonlinear structural behavior of tapered axially functionally graded material beams utilizing absolute nodal coordinate formulations

Abdur Rahman Shaukat , Bing Yan , Feilong Nie , Tengfei Wang and Jia Wang

Published/Copyright: October 6, 2025

Published by

Become an author with De Gruyter Brill

Submit Manuscript Author Information

From the journal Science and Engineering of Composite Materials Volume 32 Issue 1

Abstract

This study investigates the nonlinear structural behavior of axially functionally graded material beams composed of steel and aluminum through a systematic computational approach. Three absolute nodal coordinate formulation (ANCF)-based approaches were applied: the first employs a planar ANCF shear deformable beam element, the second utilizes an enhanced continuum mechanics-based ANCF element, and the third leverages a higher order ANCF beam element. The research methodology follows a rigorous stepwise analysis: first, the accuracy of the ANCF beam elements was verified through a convergence study; second, a large deformation analysis was carried out focusing on the effect of material nonlinearity; third, comprehensive buckling analyses were conducted, including both nonlinear postbuckling evaluation using Crisfield’s arc-length method and nonlinear eigenvalue buckling estimation via a dichotomy scheme, identifying the complex nonlinear behavior and demonstrating the variation in critical loads depending on material gradation profiles, respectively. The overall results demonstrate the ANCF framework’s effectiveness in capturing complex geometric and material nonlinearities. This study highlights the mechanical tradeoffs between steel dominated regions and aluminum rich zones, providing insights for the design of graded structural components in advanced engineering applications. It was also found that the ANCF-based higher order beam element provides more accurate and reliable deflection curves than other beam models.

Keywords: nonlinear buckling; nonlinear postbuckling; large deformation analysis; axially functionally graded material; absolute nodal coordinate formulation

1 Introduction

Functionally graded materials (FGMs) are advanced composites with spatially varying composition and properties, enabling tailored performance for specific applications [1]. Unlike conventional composites with abrupt interfaces, FGMs exhibit continuous gradients in material properties [2]. Their unique thermal, mechanical, electrical, and chemical characteristics make them suitable for diverse fields, including biomedical devices, sensors, and electronics [3]. Recent advances in computational modeling techniques and new material synthesis have further expanded the potential applications of FGMs [4]. Notably, the use of generalized thermoelastic theory in conjunction with finite element discretization has proven effective in modeling FGM behavior [5,6]. FGMs were first conceptualized in 1984 to address thermal protection challenges in aerospace applications [7]. On the other hand, non-prismatic members offer a unique combination of weight reduction and strength enhancement, meeting both architectural esthetics and civil engineering requirements [8]. These members have been applied in large span bridges, communication towers, and wind turbines. They have found applications from structural to the aerospace industry [9,10]. They provide esthetic appeal and engineering benefits, such as structural optimization, weight reduction, material savings, load distribution, and aerodynamic efficiency [11].

Recent decades have seen significant progress in large deformation and buckling analysis of tapered beams, including both homogeneous and non-homogeneous materials. Early work by Eisenberger [12] established explicit stiffness matrices for non-prismatic members, incorporating transverse shear effects. Subsequent developments by Baker [13] applied Bernoulli–Euler theory with weighted residual methods to solve large deflection problems in cantilevered non-prismatic beams under arbitrarily distributed loading. Furthermore, Brojan et al. [14] extended nonlinear analysis to non-prismatic cantilevers using Ludwick’s constitutive law, incorporating material nonlinearity under end moments. For buckling analysis, Singh and Li [15] developed an efficient low-dimensional model for functionally graded columns with variable cross-sections using equivalent piece-wise constant geometric and material properties. A comprehensive study was conducted by Kang and Li [16] comparing small and large deformation theories for FGM cantilever beams under end moments. Their results demonstrated less than 5% discrepancy in deflection for rotations below 45°, which increases to over 10% beyond this threshold. In a separate approach, Attarnejad [17] developed exact displacement functions for non-prismatic beams using flexibility methods, employing Euler-beam assumptions for large deflection analysis. Shahba and Rajasekaran [18] enhanced convergence rates in stability analysis of axially FGM (AFGM) tapered Euler–Bernoulli beams using the differential transform element method (DTEM). In subsequent work [19], they developed exact shape functions for static analysis of these beams through basic displacement function (BDF) derived via the unit-dummy-load method. In his research, Nguyen [20] developed a finite element method (FEM)-based co-rotational approach to analyze the large displacement response of tapered axially functionally graded material (AFGM) beams. Meanwhile, Niknam et al. [21] investigated the effects of thermal and mechanical loading on tapered FGM beam bending using both analytical and generalized differential quadrature methods. The investigation was carried out by Sun et al. [22] on the buckling behavior of non-prismatic FGM columns using Timoshenko beam theory and the initial value method, determining critical buckling loads and validating the approach. Separately, the static response of axially graded FGM beams was analyzed by Ghayesh and Farokhi using coupled field equations discretized via the Galerkin technique [23]. Additionally, the large transverse deformation of tapered circular FGM cantilever beams was studied by Horibe and Mori using the Runge–Kutta method [24]. A power series-based approach was developed for linear buckling analysis of AFGM Timoshenko beams with variable cross-sections, demonstrating its effectiveness in a study by Soltani and Asgarian [25]. In another approach, the finite difference method was employed by Soltani and Asgarian to predict linear buckling loads of tapered FGM Timoshenko beams demonstrating better accuracy [26].

The absolute nodal coordinate formulation (ANCF), a non-incremental nonlinear FEM, is promising for soft robotics due to its ability to model nonlinear materials. The ANCF accurately represents complex geometries using absolute coordinates and global slopes, which make it suitable for continuum mechanics-based analysis in flexible and soft robotic systems [27]. The shift from a rigid link to deformable arms marks a significant change in robotic design. Meanwhile, the innovations in smart materials, computational modeling of nonlinear deformations, and bio-inspired control strategies provide exciting opportunities in the field. This research aims to investigate a new computational tool to study the deformability characteristics of structural components with non-uniform cross section and graded material compositions. In this study, the nonlinear structural behavior of tapered AFGM cantilever beams using ANCF-based approaches was investigated. Nonlinear structural properties of AFGM beams, whose material distribution along the longitudinal axis is governed by a power law, with steel and aluminum as the constitutive materials, were explored. Large deformation analyses were conducted utilizing the full Newton–Raphson algorithm, while the nonlinear post-buckling responses were obtained through Crisfield’s arc-length method. Finally, nonlinear eigenvalue buckling loads using a dichotomy method based scheme was estimated. Table 1 summarizes the existing studies on AFGM beams, highlighting various analytical and numerical methods used by different researchers.

Table 1

AFGM beams found in the literature review

Year	Authors	Beam	Loading	Method
1991	Eisenberger [12]	Non-prismatic members	Transverse shear	The stiffness matrices including the shear effect are derived explicitly
1993	Baker [13]	Non-prismatic cantilever	Transverse arbitrarily distributed loads	The governing equations were derived using Bernoulli–Euler beam theory and solved via the weighted residual method
2009	Brojan et al. [14]	Non-prismatic cantilever	End moment loading	This study analyzed a materially nonlinear beam governed by the Ludwick constitutive law, deriving exact moment-curvature relations for rectangular cross-sections
2009	Singh and Li [15]	FGM non-uniform column	Compressive load	This study analyzed FGM columns with variable cross-sections using a low-dimensional mathematical model with piecewise discretization of material and geometric characteristics
2010	Kang and Li [16]	FGM cantilever beam	End moment loading	The non-linear cantilever FGM beam with end moment was investigated; the explicit formulas for deflection and rotations were derived
2010	Attarnejad [17]	Non-prismatic beams	Transverse loading	The authors implemented the flexibility method by employing fundamental displacement functions to satisfy element compatibility equations
2012	Shahba and Rajasekaran [18]	AFGM beam	Compressive loading	The Euler–Bernoulli assumption with DTEM was utilized to improve the convergence rate
2013	Shahba et al. [19]	AFGM tapered beam	Compressive loading	The BDF via the unit-dummy method was employed to obtain the exact shape functions
2013	Nguyen [20]	AFGM tapered beam	Transverse loading	The FEM-based co-rotational beam elements was used
2014	Niknam et al. [21]	Tapered FGM beam	Thermal and mechanical transverse loadings	A closed-form solution was derived from the governing equations, with solution methodologies employing both the Galerkin method and generalized differential quadrature method
2016	Sun et al. [22]	Varying cross-section column	Axial and transverse loadings	The critical tip force and axial load were determined using Timoshenko beam theory combined with the initial value method
2018	Ghayesh and Farokhi [23]	AFGM tapered beam	Transverse loading	Hamilton’s principle formulated the governing equations, which were then discretized via the Galerkin modal approach
2018	Horibe and Mori [24]	AFGM tapered cantilever beam	Transverse loading	The potential energy method was employed, with the Runge–Kutta method used for numerical solution
2018	Soltani and Asgarian [25]	AFGM tapered beam	Axial compression	Stiffness matrices for linear stability analysis were derived using power series expansion
2020	Soltani et al. [26]	AFGM tapered beam	Compressive load	Timoshenko beam theory with small displacement assumption was used; while, the finite difference method was employed to discretize the problem in order to obtain linear buckling loads

2 ANCF-based approaches

Shabana introduced a non-incremental nonlinear finite element formulation for studying the dynamics of flexible multibody systems, termed the ANCF [28]. Unlike conventional FEM, which rely on rotation angles as nodal coordinates, the ANCF employs absolute nodal positions and global slopes, thereby circumventing the need for infinitesimal or finite rotation assumptions. This approach ensures an exact representation of rigid body motion while maintaining a constant mass matrix and eliminating Coriolis and centrifugal forces in the system’s equations of motion [29]. The ANCF framework has been extended to a wide range of structural elements – such as beams, solids, and plates – while incorporating diverse shape functions, thereby significantly expanding the ANCF element library [30,31]. Although originally developed for dynamic analysis, the ANCF has also proven effective in nonlinear structural analysis. This section provides a brief overview of the three ANCF approaches employed in the present study.

2.1 Omar–Shabana beam element (OmSh)

The shear deformable planar beam element based on ANCF was proposed by Omar and Shabana. In this element, elastic forces were derived using the principles of general continuum mechanics [32]. In this study, this ANCF element is referred to as “OmSh” for convenience. The position vector is expressed as a function of element nodal coordinates, r = S e , where e represents the element nodal coordinates vector as given below:

(1) e = r i T r x i T r y i T r j T r x j T r y j T T .

In this expression, r represents the vector of global position between nodes i and j . While r x = ∂ r / ∂ x and r y = ∂ r / ∂ y represent gradient vectors. In this method, the shape function matrix is expressed as follows:

(2) S = s 1 I s 2 I … s 6 I .

In this equation, I is a 2 × 2 identity matrix and s i ( i = 1 , 2 , … , 6 ) are shape functions defined as

(3) s 1 = 1 − 3 ξ 2 + 2 ξ 3 , s 2 = l ( ξ − 2 ξ 2 + ξ 3 ) , s 3 = l ( η − ξ η ) s 4 = 3 ξ 2 − 2 ξ 3 , s 5 = l ( − ξ 2 + ξ 3 ) , s 6 = l ξ η .

In the above expressions, x and y represent the local coordinates, while l denotes the length of element. The coordinates in the non-dimensional form ξ = x / l and η = y / l were used, where 0 ≤ x ≤ l and − h / 2 ≤ y ≤ h / 2 . As illustrated in Figure 1, the matrix of position vector gradients is expressed as

(4) J = ∂ r ∂ X = ∂ r ∂ x ∂ X ∂ x − 1 = J e J 0 − 1 .

Figure 1

ANCF beam element.

In the Jacobian matrix, J e = ∂ r / ∂ x and J 0 = ∂ X / ∂ x , while the straight configuration axis is given by x . The reference configuration was expressed in terms of X = S ( x ) e 0 , where e 0 is the corresponding nodal vector. The Green–Lagrange strain tensor is expressed as

(5) ε m = 1 2 ( J T J − I ) .

The term I denotes a 2 × 2 identity matrix, and the strain is expressed in Voigt notation as:

(6) ε v = ε x x ε y y 2 ε x y T .

The terms ε x x , ε y y , and ε x y represent the axial, transverse, and shear components of the strain tensor, respectively. For the linear elastic isotropic material, the strain energy is expressed as

(7) U = 1 2 ∫ V ε v T E m ε v ∣ J 0 ∣ d V .

In the above equation E m represents the elastic coefficient matrix. Using the material law for the homogeneous isotropic material, the elastic coefficient matrix under plane stress assumption is expressed as

(8) E m = E 1 − γ 2 1 γ 0 γ 1 0 0 0 ( 1 − γ ) / 2 .

In the above equation E denotes Young’s modulus and γ represents Poisson’s ratio. The differentiation of the elastic energy by the nodal vector gives the elastic force as:

(9) Q s = ∂ U ∂ e = ∫ V ∂ ε v ∂ e T E m ε v ∣ J 0 ∣ d V .

This expression relates the volume V of the straight configuration to the volume V 0 of the initially curved (referenced) configuration straightforwardly by the simple expression d V 0 = det ( J 0 ) d V . Hence, the integration process is simple and straightforward. The straight configuration is directly used to find the deformation for the initially curved configuration. The initially curved stress-free configurations was utilized for nonlinear structural analysis. This is one of the benefits of using ANCF. A detailed discussion is provided in the study of Shabana [33].

2.2 Higher order beam element (HOBE)

In addition to the strain-split method, Patel and Shabana proposed a planar HOBE with 16 degrees of freedom as an alternative to eliminate locking [34]. The HOBE has a curvature vector r y y , and the element nodal vector is expressed as

(10) e = r i T r x i T r y i T r y y i T r j T r x j T r y j T r y y j T T .

In this representation r denotes the position vector. The terms r x = ∂ r / ∂ x , r y = ∂ r / ∂ y , and r y y = ∂ 2 r / ∂ y 2 represent the gradient vectors in the beam element. The shape function matrix is given by

(11) S = s 1 I s 2 I … s 8 I ,

where I represents the 2 × 2 identity matrix, and shape functions s i ( i = 1 , 2 , … , 8 ) are defined as

(12) s 1 = 1 − 3 ξ 2 + 2 ξ 3 , s 2 = l ( ξ − 2 ξ 2 + ξ 3 ) , s 3 = l ( η − ξ η ) , s 4 = l 2 ( η 2 − ξ η 2 ) / 2 s 5 = 3 ξ 2 − 2 ξ 3 , s 6 = l ( − ξ 2 + ξ 3 ) , s 7 = l ξ η , s 8 = l 2 ( ξ η 2 ) / 2 .

The variables x , y , l , ξ , and η are the same as shown previously. For the HOBE, the strain energy is expressed as U = 1 2 ∫ V ε v T E m ε v det ( J 0 ) d V , and the corresponding elastic force is obtained by derived by differentiating the strain energy Q s = ∂ U ∂ e = ∫ V ∂ ε v ∂ e T E m ε v det ( J 0 ) d V .

2.3 Enhanced continuum mechanics formulation (OS-EnCM)

Due to the coupling that exists between the axial strain and transverse strain in the standard beam element proposed by Omar and Shabana, it exhibits a Poisson locking phenomenon. Gerstmayr et al. [35] introduced the enhanced continuum mechanics (EnCM) method to address the locking issue. In this method, the strain energy was calculated in two distinct parts. The Poisson ratio was removed in the transverse direction, while it was incorporated solely along the beam centreline. Consequently, this method was utilized with the “OmSh” beam element, hence referred to as “OS-EnCM” in this study. The expression for the elastic energy is given by

(13) U = 1 2 b ∫ − h / 2 h / 2 ∫ 0 l ε T D 0 ε d x d y + 1 2 b h ∫ 0 l ε T D v ε d x .

In the above expression, b represents the thickness of the element. As the first component does not include the Poisson effect, the elastic coefficient matrix is expressed as

(14) D 0 = E 0 0 0 E 0 0 0 G k s .

In this equation, G = E / 2 ( 1 + ν ) denotes the shear modulus. Conversely, the elastic coefficient matrix is established as

(15) D v = E ν 1 − ν 2 ν 1 0 1 ν 0 0 0 0 .

This component incorporates the Poisson effect in the EnCM beam element.

3 Computational strategies

The nonlinear large deformation analysis was conducted using the full Newton–Raphson algorithm. However, postbuckling analysis was performed utilizing Crisfield’s adoption of arc-length method given by [36]

(16) Q in ( e ) − λ ext q ext = 0 ,

where the external force is expressed in terms of λ ext denoting the scalar loading parameter and q ext represents the constant directional vector of the external force F ext = λ ext q ext . The internal force Q in ( e ) is expressed as a nonlinear function of nodal coordinates. The change, as compared to the full Newton–Raphson algorithm, is an addition of the scalar loading parameter λ ext [37]. Consequently, the dimensionality of the displacement field transforms N into an N + 1 load-displacement field. The constraint equation by Crisfield is expressed in the following form:

(17) Δ u T Δ u = Δ l 2 .

The arc-length parameter Δ l was not changed while tracing the complete equilibrium path and Δ u denotes the incremental displacement vector. At each iteration step, the solution is obtained by tracing the arc with the help of the arc-length parameter Δ l until it converges to the static equilibrium. Nevertheless, the additional constraint introduced causes the arc follow the overall cylindrical curve. Consequently, it is also referred to as the “Cylindrical arc-length method.” The detailed implementation of this method can be found in the study of Shaukat et al. [38]. The constant arc-length parameter was used as Δ l = 0.2 . At some instances where the algorithm finds difficulty during convergence it was allowed to halve the parameter.

Wang et al. utilized the energy criterion for the nonlinear eigenvalue buckling analysis. This method was based on the concept of the relative minimum of the system’s total potential energy. They employed a dichotomy-based nonlinear iterative method to determine the critical load of the column. The critical point is identified when the second-order derivative of the strain energy equals zero:

(18) δ 2 U elastic = 0 .

This corresponds to the subsequent equation, whereby the determinant of the system’s tangent stiffness matrix under loading becomes singular:

(19) det ( K t ) = 0 .

In the dichotomy-based method, eigenvalues of the corresponding tangent stiffness matrix were monitored by using the nonlinear iterations. The applied compressive load was increased until the solution converged to the zero eigenvalue. This method was applied to find the critical buckling loads of the tapered AFGM columns. The details can be found in the study of Wang and Wang [39].

4 Numerical examples

In this section, the AFGM beam composed of steel and aluminium as constituent materials was studied. In the study of nonlinear structural analysis, the large deformation analysis was considered as the benchmark. In the first part, the ANCF elements used in this article were verified through the convergence study. Then, large deformation analysis was presented. In the second part, the nonlinear postbuckling analysis was performed. Finally, the nonlinear buckling analysis was conducted.

4.1 Convergence

In order to establish the convergence of the ANCF approaches utilized in this study, a large deformation analysis of tapered AFGM cantilever beams under a constant transverse load was conducted.

The tapered beam was considered to be made of an AFGM composed of aluminum and steel. The material properties considered were Young’s modulus of elasticity, for steel E s = 2.1 × 10 11 N/m 2 and for aluminum E a = 7 × 10 10 N/m 2 , with Poisson’s ratio ν = 0.3 assumed for both. The geometrical properties included the length of beam L = 0.5 m , initial cross-sectional height h i = 0.01 m , and width b = 0.01 m . The final cross-sectional height was defined by the taper parameter as h f = ( 1 − β ) h i , with the taper parameter β = 0.5 . Figure 2 demonstrates the geometrical configuration and loading position on the tapered AFGM beam. The gradual change of the Young’s modulus in the x-axis was governed by the power law distribution, expressed as follows:

(20) E ( x ) = ( E a − E s ) ( x / L ) m + E s ,

where m describes the distribution profile which can be linear, quadratic, cubic, or non-zero and non-negative fraction less than one. In this study, m = 0.2 , 1 , and 5 was selected. At the free end, the vertical downward force P = 1,000 N was employed. The deformations in the vertical and horizontal dimensions at the tip end were recorded as detailed in Tables 2–4. Consequently, the convergence of the three ANCF approaches, i.e., OmSh, OS-EnCM, and HOBE, was assessed through the graphs of load vs number of elements, as illustrated in Figures 3–5.

Figure 2

Initial configuration for large deformation analysis.

Table 2

Large deformation at the free end of the tapered AFGM beam using the ANCF HOBE

Elements	m = 0.2	m = 0.2	m = 1	m = 1	m = 5	m = 5
Elements	Ux	uy	ux	uy	ux	uy
10	−0.170824	−0.316196	−0.122067	−0.272899	−0.096744	−0.252799
20	−0.174853	−0.319311	−0.124821	−0.275357	−0.098856	−0.255367
30	−0.175191	−0.319612	−0.125055	−0.275599	−0.099052	−0.255636
40	−0.175289	−0.319713	−0.125125	−0.275684	−0.099116	−0.255734
50	−0.175335	−0.319765	−0.125158	−0.275728	−0.099147	−0.255785
60	−0.175361	−0.319797	−0.125178	−0.275756	−0.099167	−0.255818
70	−0.175378	−0.319818	−0.125190	−0.275775	−0.099180	−0.255840
80	−0.175389	−0.319833	−0.125199	−0.275788	−0.099189	−0.255856
90	−0.175398	−0.319844	−0.125206	−0.275799	−0.099196	−0.255868
100	−0.175404	−0.319853	−0.125211	−0.275807	−0.099201	−0.255878

Table 4

Large deformation at the free end of the tapered AFGM beam using the ANCF OS-EnCM

Elements	m = 0.2	m = 0.2	m = 1	m = 1	m = 5	m = 5
Elements	Ux	uy	ux	uy	ux	uy
10	−0.171342	−0.316922	−0.122499	−0.273619	−0.097205	−0.253625
20	−0.175109	−0.319700	−0.125046	−0.275752	−0.099106	−0.255828
30	−0.175355	−0.319869	−0.125206	−0.275867	−0.099220	−0.255950
40	−0.175410	−0.319906	−0.125239	−0.275889	−0.099243	−0.255973
50	−0.175430	−0.319920	−0.125250	−0.275895	−0.099251	−0.255980
60	−0.175440	−0.319927	−0.125255	−0.275898	−0.099255	−0.255983
70	−0.175446	−0.319931	−0.125258	−0.275899	−0.099257	−0.255985
80	−0.175449	−0.319933	−0.125260	−0.275900	−0.099258	−0.255986
90	−0.175452	−0.319935	−0.125261	−0.275900	−0.099258	−0.255987
100	−0.175453	−0.319937	−0.125262	−0.275901	−0.099259	−0.255987

Figure 3

Convergence of the AFGM beam for power law index m = 0.2.

Figure 4

Convergence of the AFGM beam for power law index m = 1.

Figure 5

Convergence of AFGM beam for power law index m = 5.

The various deflected configurations of tapered cantilever beams subjected to transverse loading are depicted in Figure 6. Table 2 presents the data for the HOBE elements, Table 3 for the OmSh element, and Table 4 provides the numerical values for OS-EnCM elements employed for the analysis.

Figure 6

AFGM cantilever beam under large deformation.

Table 3

Large deformation at the free end of the tapered AFGM beam using the ANCF OmSh

Elements	m = 0.2	m = 0.2	m = 1	m = 1	m = 5	m = 5
Elements	Ux	uy	ux	uy	ux	uy
10	−0.161559	−0.309282	−0.113037	−0.264094	−0.088004	−0.242590
20	−0.164812	−0.311713	−0.115177	−0.265942	−0.089512	−0.244427
30	−0.165029	−0.311863	−0.115315	−0.266043	−0.089604	−0.244529
40	−0.165079	−0.311897	−0.115344	−0.266062	−0.089623	−0.244549
50	−0.165097	−0.311910	−0.115355	−0.266067	−0.089630	−0.244555
60	−0.165106	−0.311917	−0.115360	−0.266070	−0.089633	−0.244558
70	−0.165112	−0.311920	−0.115362	−0.266071	−0.089634	−0.244559
80	−0.165115	−0.311923	−0.115364	−0.266072	−0.089635	−0.244560
90	−0.165117	−0.311925	−0.115365	−0.266072	−0.089636	−0.244561
100	−0.165119	−0.311926	−0.115365	−0.266073	−0.089636	−0.244561

In the convergence study, the OmSh approach exhibited locking behavior, resulting in lower predicted vertical deflections compared to the HOBE and OS-EnCM approaches. The OS-EnCM effectively alleviated locking, while the HOBE approach performed well, yielding accurate deflection predictions. A convergence study revealed that a mesh configuration of 50 elements ensured convergence up to four decimal points for all three ANCF-based approaches. Consequently, a 50-element mesh was used throughout the study to maintain consistency.

4.2 Large deformation analysis

This section investigates the effect of material nonlinearity during the large deflection of the tapered AFGM beam. An increasing transverse load was employed at the tip end of the tapered cantilever AFGM beam utilizing the parameters used previously. The complete nonlinear vertical and horizontal deflections in response to the increasing transverse load were obtained by using Crisfield’s arc-length algorithm. The resulting beam deflections at various load levels are depicted in Figure 7.

Figure 7

Different AFGM beam configurations under increasing transverse loading.

A comparison of the vertical and horizontal deflections reveals that decreasing the material power law index to 0.2 has the most pronounced effect on the deflection path as compared to m = 1 and 5. Figures 8–10 show the normalized displacement vs normalized load curves for the AFGM beam with power law index m = 5, m = 1, and m = 0.2, respectively. However, the ANCF-based approaches with 50 elements were compared with a reference solution and were found to be in good agreement [20].

Figure 8

Large deformation under increasing transverse loading for power law index m = 5.

Figure 9

Large deformation under increasing transverse loading for power law index m = 1.

Figure 10

Large deformation under increasing transverse loading for power law index m = 0.2.

4.3 Nonlinear postbuckling analysis

The nonlinear postbuckling behavior of the AFGM beam utilizing the parameters from the previous example was investigated. The eccentric compressive force at the top corner of the free end of the considered beam was used. In Figure 11, the geometrical configuration and loading position on the tapered AFGM column can be observed. Crisfield’s arc-length method was employed to trace the postbuckling response.

Figure 11

Initial configuration for nonlinear postbuckling analysis.

The obtained AFGM beam deformations are shown in Figure 12.

Figure 12

AFGM beam configurations under eccentric compressive load.

To study the effect of material inhomogeneity on the AFGM beam, the power law index m = 0.2 , 1 , and 5 were employed in postbuckling analysis. The results were obtained for the complete nonlinear equilibrium behavior using 50 elements. The load vs vertical and horizontal deflections are illustrated in Figures 13–15. It was observed that for m = 0.2, the postbuckling behavior manifested earlier, and as the power law index increased, the stiffness of the AFGM beam also increased. Considerably, for the power law index m = 5, the onset of postbuckling behavior occurs at a relatively higher value. The results obtained were compared with results in the literature and were found to be in good agreement [20].

Figure 13

Postbuckling analysis under eccentric compressive loading for power law index m = 5.

Figure 14

Postbuckling analysis under eccentric compressive loading for power law index m = 1.

Figure 15

Postbuckling analysis under eccentric compressive loading for power law index m = 0.2.

4.4 Buckling analysis

The buckling analysis was first conducted on a homogeneous tapered cantilever beam to validate the ANCF approaches applied in this study.

The material properties considered were Young’s modulus of elasticity E = 2 × 10 11 N/m 2 and Poisson’s ratio ν = 0.3 . The geometric properties of the beam included a total length of beam L = 0.028867 m and a cross-sectional width b = 0.01 m . The tapered cross-sectional height was obtained using h f = ( 1 − β ) h i , where the initial cross section height h i = 0.01 m and taper parameter β = 0.1 ≤ 0.9 were used. The geometrical configuration and loading position on the tapered AFGM column are shown in Figure 16. For comparison, this example follows the configuration presented in the study of Soltani and Asgarian [25]. The buckling load results obtained using three ANCF approaches are provided in Table 5. The OmSh element slightly over-predicted the critical loads. While OS-EnCM slightly under-predicted it. On the other hand, the HOBE results agree well with the reference solution, as shown in Table 5.

Figure 16

Initial configuration for buckling analysis.

Table 5

Comparison of homogeneous beam buckling loads

Taper	Pnorm	Pcr	HOBE	OmSh	EnCM	BEAM188 n = 200
0.1	2.088	417,600	412999.48	442825.59	406050.65	424,870
0.2	1.884	376,800	373547.67	400668.03	367282.75	383,103
0.3	1.676	335,200	333296.04	357634.96	327736.29	340,777
0.4	1.465	293,000	292123.94	313589.03	287286.46	297,780
0.5	1.25	250,000	249860	268336.95	245757.48	253,835
0.6	1.029	205,800	206245.42	221591.96	202888.31	208,981
0.7	0.8	160,000	160862.7	172896.69	158262.26	162,507
0.8	0.56	112,000	112963.88	121442.26	111141.37	113,766
0.9	0.301	60,200	60981.17	65548.82	59990.27	61,167

In the second part, buckling analysis of the non-homogeneous tapered cantilever beam made of an AFGM composed of aluminum and steel was considered. The material properties considered were Young’s modulus for steel E s = 2 × 10 11 N/m 2 and for aluminum E a = 7 × 10 10 N/m 2 , where the Poisson’s ratio ν = 0.3 was considered as per the reference solution. The geometrical properties of the beam including the length L = 0.028867 m , initial cross-sectional height h i = 0.01 m , and width b = 0.01 m were assumed. The final cross-sectional height was determined using the taper parameter h f = ( 1 − β ) h i , where the taper parameter varies between β = 0.1 ≤ 0.9 . The gradually varying Young’s modulus in the axial direction followed the power law distribution as:

(21) E ( x ) = ( E a − E s ) ( x / L ) m + E s .

The buckling load estimation for the three ANCF-based approaches of the tapered AFGM beam with different taper parameters was conducted. The results obtained for the power law index m = 1 are shown in Table 6. As illustrated in Figure 17, the buckling loads decrease as the taper ratio increases.

Table 6

Comparison of buckling loads for the AFGM beam with power law index m = 1

Taper	Pnorm	Pcr	HOBE	OmSh	EnCM	BEAM188 n = 200
0	1.7101	342,020	343009.64	364622.77	335165.54	351,149
0.1	1.5338	306,760	308287.71	327746.79	301199.08	314,844
0.2	1.3574	271,480	273354.87	290646.81	267042.43	278,517
0.3	1.1808	236,160	238237.52	253347.09	232717.87	242,190
0.4	1.0043	200,860	202975.22	215886.5	198261.06	205,896
0.5	0.828	165,600	167629.75	178328.74	163730.42	169,687
0.6	0.6524	130,480	132303.72	140782.09	129225.25	133,658
0.7	0.4783	95,660	97178.68	103441.02	94924.09	97,979
0.8	0.3076	61,520	62608.32	66685.58	61176.13	62,993
0.9	0.1538	30,760	29406.43	31372.47	28769.78	29,483

Figure 17

Buckling loads for the homogeneous tapered cantilever beam.

The comparison of buckling loads for the three ANCF-based approaches for m = 2 was conducted, and the results are shown in Table 7. In Figure 18, the decreasing buckling load can be observed as the taper ratio increases and the capacity to bear compressive load diminishes.

Table 7

Comparison of buckling loads of FGM beam with m = 2

Taper	Pnorm	Pcr	HOBE	OmSh	EnCM	BEAM188 n = 200
0	1.9815	396,300	397172.01	422695.64	388506.84	406,702
0.1	1.7847	356,940	358781.1	381764.77	350818.69	366,350
0.2	1.5868	317,360	319827.46	340244.31	312607.03	325,678
0.3	1.3871	277,420	280287.49	298107.36	273846.76	284,659
0.4	1.1855	237,100	240142.55	255333.9	234518.78	243,271
0.5	0.982	196,400	199388.91	211923.34	194621.35	201,503
0.6	0.7766	155,320	158063.86	167924.09	154196.9	159,385
0.7	0.5699	113,980	116310.82	123506.64	113398.93	117,046
0.8	0.3655	73,100	74563.27	79160.93	72673.35	74,895
0.9	0.1684	33,680	34159.99	36314.8	33329.86	34,208

Figure 18

Buckling loads for the AFGM beam for power law index m = 1.

The results for the AFGM column with power law index m = 3 for the three ANCF-based approaches employed are shown in Table 8. In Figure 19, the same trend was observed as in previous section, the critical load decreases with the increase in taper ratio.

Table 8

Buckling loads of the FGM beam with m = 3

Taper	Pnorm	Pcr	HOBE	OmSh	EnCM	BEAM188 n = 200
0	2.0966	419,320	419358.32	447137.48	410847.91	430,058
0.1	1.8955	379,100	380305.24	405383.81	372412.35	388,788
0.2	1.6919	338,380	340525.34	362861.28	333291.04	347,052
0.3	1.4857	297,140	299949.79	319493.92	293414.07	304,781
0.4	1.2764	255,280	258496.19	275194.42	252700.69	261,894
0.5	1.0635	212,700	216072.2	229869.11	211063.59	218,299
0.6	0.8465	169,300	172589.68	183438.48	168426.94	173,907
0.7	0.6271	125,420	128021.65	135906.77	124789.18	128,689
0.8	0.4028	80,560	82611.99	87597.56	80437.81	82,873
0.9	0.2082	41,640	37729.96	40023.63	36750.92	37,741

Figure 19

Buckling loads for the AFGM beam with power law index m = 2.

From Figures 18–20, it was found that OmSh tends to slightly over-predict the critical buckling load, whereas OS-EnCM and HOBE approaches produce results that closely match well with the reference solution [25].

Figure 20

Buckling loads for the AFGM beam with power law index m = 3.

5 Conclusions

This study comprehensively examined the nonlinear structural behavior of tapered AFGM cantilever beams, composed of steel and aluminum, utilizing three distinct ANCF-based approaches. A series of numerical examples involving homogeneous and non-homogeneous AFGM beams with non-uniform cross-sections were meticulously analyzed. The results obtained from the ANCF approaches were rigorously compared with existing numerical solutions to underscore the efficacy of the ANCF method. To ensure the reliability of the findings, a convergence study was initially conducted to validate the numerical stability and accuracy of each ANCF-based approach. Following that, the nonlinear deflection response under transverse loading was evaluated, providing insights into the nonlinear structural behavior under transverse loading conditions. Detailed simulations of the complete postbuckling equilibrium path for AFGM columns with varying power law indices were performed, offering a comprehensive understanding of the postbuckling behavior. Additionally, critical buckling loads were estimated across different taper parameters, further establishing the ability of the ANCF method in predicting stability characteristics of the AFGM columns. The following conclusions were drawn from these analyses:

The OmSh element exhibited slight locking behavior while predicting large deformations under transverse loading, leading to deviations from the reference solution. The OS-EnCM-based approach demonstrated accuracy with results closely aligning with benchmark solution. The HOBE provided the most accurate predictions, achieving good agreement with the reference solution, thereby validating its effectiveness in capturing geometrically nonlinear responses.
The variations in the power-law index “m” influenced the nonlinear postbuckling response. For m = 5, the HOBE results matched well with the reference solution, while OmSh over-predicted and OS-EnCM under-predicted the postbuckling equilibrium path. Similar trends were observed for m = 1 and 0.2, though the discrepancies diminished with decreasing m, highlighting the sensitivity of postbuckling behavior to material gradation. Increasing m generally increases the stiffness of the structure but reduces its post-buckling flexibility. This trade-off is important to consider when selecting the appropriate value of m for a specific application. Based on our analysis, we believe that m = 1 offers a good compromise between stiffness and post-buckling flexibility.
The HOBE formulation consistently produced buckling load estimates in good agreement with the existing literature for power-law indices m = 1, 2, 3. OmSh exhibited a slight over-prediction tendency, whereas OS-EnCM marginally under-predicted critical loads, emphasizing the need for careful selection of the ANCF approach based on desired accuracy.

Among the examined ANCF-based approaches, the HOBE emerges as the most reliable and accurate method for conducting nonlinear structural analyses of tapered AFGM beams. This is particularly true for applications that demand high precision in predicting large deformations, postbuckling behavior, and critical buckling loads. Although the OmSh and OS-EnCM methods are also viable alternatives, their propensity to either slightly over- or under-predict structural responses indicates that their use should be carefully considered based on the specific accuracy requirements of the problem at hand. These findings provide valuable insights for engineers and researchers working with FGMs. It offers insights for selecting an appropriate ANCF-based approach in the design and analysis of advanced AFGM structures subjected to different loading conditions. Future work could explore the extension of these ANCF approaches to more complex geometries, dynamic loading conditions, and other nonlinear material compositions.

Acknowledgments

The authors acknowledge the Shanghai Zhenhua Heavy Industries Co., Ltd. (ZPMC) for their commitment toward scientific research.

Funding information: This research was funded by the postdoctoral fellowship under grant number RP2024045736/9908000471 for the project “Research on the buckling problems of box girder structures of port cranes.”
Author contributions: All authors have accepted responsibility for the entire content of this manuscript and consented to its submission to the journal, reviewed all the results and approved the final version of the manuscript. Conceptualization, methodology, and writing – Abdur Rahman Shaukat; software – help and validation – Jia Wang; formal analysis and investigation – Tengfei Wang; resources and data curation – Feilong Nie; supervision, project administration, and funding acquisition – Bing Yan.
Conflict of interest: The authors state no conflict of interest.
Supplementary materials: None.
Data availability statement: Data can be provided on reasonable request to the corresponding author.

References

[1] Naebe M, Shirvanimoghaddam K. Functionally graded materials: A review of fabrication and properties. Applied materials today. Vol. 5. Netherlands: Elsevier Ltd; 2016. p. 223–45.10.1016/j.apmt.2016.10.001Search in Google Scholar

[2] Boggarapu V, Gujjala R, Ojha S, Acharya S, Venkateswara Babu P, Chowdary S, et al. State of the art in functionally graded materials. Composite structures. Vol. 262. England: Elsevier Ltd; 2021.10.1016/j.compstruct.2021.113596Search in Google Scholar

[3] Li Y, Feng Z, Hao L, Huang L, Xin C, Wang Y, et al. A review on functionally graded materials and structures via additive manufacturing: From multi-scale design to versatile functional properties. Advanced materials technologies. Vol. 5. United States: Wiley-Blackwell; 2020.10.1002/admt.201900981Search in Google Scholar

[4] Kumar P, Sharma SK, Singh RKR. Recent trends and future outlooks in manufacturing methods and applications of FGM: a comprehensive review. Materials and manufacturing processes. Vol. 38. United States: Taylor and Francis Ltd; 2023. p. 1033–67.10.1080/10426914.2022.2075892Search in Google Scholar

[5] Abbas IA. Generalized thermoelastic interaction in functional graded material with fractional order three-phase lag heat transfer. J Cent South Univ. 2015;22(5):1606–13.10.1007/s11771-015-2677-5Search in Google Scholar

[6] Abbas IA, Othman MIA. Plane waves in generalized thermo-microstretch elastic solid with thermal relaxation using finite element method. Int J Thermophys. 2012;33(12):2407–23.10.1007/s10765-012-1340-8Search in Google Scholar

[7] Koizumi M. FGM activities in Japan. Compos Part B Eng. 1997;28(1–2):1–4.10.1016/S1359-8368(96)00016-9Search in Google Scholar

[8] Mercuri V, Balduzzi G, Asprone D, Auricchio F. Structural analysis of non-prismatic beams: Critical issues, accurate stress recovery, and analytical definition of the Finite Element (FE) stiffness matrix. Eng Struct. 2020;213:110252.10.1016/j.engstruct.2020.110252Search in Google Scholar

[9] Ramakrishnan V, Rajasekaran S, Krishnamoorthy R. Limit design for concrete bridges. In American Concrete Institute, ACI Special Publication. United States: American Concrete Institute; 1971. p. 552–81.Search in Google Scholar

[10] Hasan Z, Hammoudeh H, Atmeh G. Design and analysis of a smart composite wing. In ASME 2011 International Mechanical Engineering Congress and Exposition, IMECE 2011. United States: American Society of Mechanical Engineers (ASME); 2011. p. 213–22.10.1115/IMECE2011-64802Search in Google Scholar

[11] Jeong SK, Yoon SJ, Cho SK. Elastic local buckling of transversely isotropic thin-walled compression members with varying thickness. In Key engineering materials. Switzerland: Trans Tech Publications Ltd; 2004. p. 627–32.10.4028/www.scientific.net/KEM.261-263.627Search in Google Scholar

[12] Eisenberger M. Stiffness matrices for non-prismatic members including transverse shear. Comput Struct. 1991;40(4):831–45.10.1016/0045-7949(91)90312-ASearch in Google Scholar

[13] Baker G. On the large deflections of non-prismatic cantilevers with a finite depth. Comput Struct. 1993;46(2):365–70.10.1016/0045-7949(93)90201-NSearch in Google Scholar

[14] Brojan M, Videnic T, Kosel F. Large deflections of nonlinearly elastic non-prismatic cantilever beams made from materials obeying the generalized Ludwick constitutive law. Meccanica. 2009;44(6):733–9.10.1007/s11012-009-9209-zSearch in Google Scholar

[15] Singh KV, Li G. Buckling of functionally graded and elastically restrained non-uniform columns. Compos Part B Eng. 2009;40(5):393–403.10.1016/j.compositesb.2009.03.001Search in Google Scholar

[16] Kang YA, Li XF. Large deflections of a non-linear cantilever functionally graded beam. J Reinf Plast Compos. 2010;29(12):1761–74.10.1177/0731684409103340Search in Google Scholar

[17] Attarnejad R. Basic displacement functions in analysis of nonprismatic beams. Eng Comput (Swansea, Wales). 2010;27(6):733–45.10.1108/02644401011062117Search in Google Scholar

[18] Shahba A, Rajasekaran S. Free vibration and stability of tapered Euler-Bernoulli beams made of axially functionally graded materials. Appl Math Model. 2012;36(7):3094–111.10.1016/j.apm.2011.09.073Search in Google Scholar

[19] Shahba A, Attarnejad R, Hajilar S. A mechanical-based solution for axially functionally graded tapered euler-bernoulli beams. Mech Adv Mater Struct. 2013;20(8):696–707.10.1080/15376494.2011.640971Search in Google Scholar

[20] Nguyen DK. Large displacement response of tapered cantilever beams made of axially functionally graded material. Compos Part B Eng. 2013;55:298–305.10.1016/j.compositesb.2013.06.024Search in Google Scholar

[21] Niknam H, Fallah A, Aghdam MM. Nonlinear bending of functionally graded tapered beams subjected to thermal and mechanical loading. Int J Non Linear Mech. 2014;65:141–7.10.1016/j.ijnonlinmec.2014.05.011Search in Google Scholar

[22] Sun DL, Li XF, Wang CY. Buckling of standing tapered Timoshenko columns with varying flexural rigidity under combined loadings. Int J Struct Stab Dyn. 2016;16(6):1550017.10.1142/S0219455415500170Search in Google Scholar

[23] Ghayesh MH, Farokhi H. Bending and vibration analyses of coupled axially functionally graded tapered beams. Nonlinear Dyn. 2018;91(1):17–28.10.1007/s11071-017-3783-8Search in Google Scholar

[24] Horibe T, Mori K. Large deflections of tapered cantilever beams made of axially functionally graded material. Mech Eng J. 2018;5(1):17-00268.10.1299/mej.17-00268Search in Google Scholar

[25] Soltani M, Asgarian B. Finite element formulation for linear stability analysis of axially functionally graded nonprismatic timoshenko beam. Int J Struct Stab Dyn. 2019;19(2):1950002.10.1142/S0219455419500020Search in Google Scholar

[26] Soltani M, Asgarian B, Jafarzadeh F. Finite difference method for buckling analysis of tapered Timoshenko beam made of functionally graded material. AUT J Civ Eng. 2020;4(1):91–102.Search in Google Scholar

[27] Shabana AA. Continuum-based geometry/analysis approach for flexible and soft robotic systems. Soft Robot. 2018;5(5):613–21.10.1089/soro.2018.0007Search in Google Scholar PubMed

[28] Shabana AA. Definition of the slopes and the finite element absolute nodal coordinate formulation. Multibody Syst Dyn. 1997;1(3):339–48.10.1023/A:1009740800463Search in Google Scholar

[29] Shabana AA. Definition of ANCF finite elements. J Comput Nonlinear Dyn. 2015;10(5):054506.10.1115/1.4030369Search in Google Scholar

[30] Yoo WS, Dmitrochenko O, Pogorelov DY. Review of finite elements using absolute nodal coordinates for large-deformation problems and matching physical experiments. In International Design Engineering Technical Conferences & Computers and Information in Engineering Conference September 24–28, 2005, Long Beach, California USA; 2005. p. 1–10.Search in Google Scholar

[31] Gerstmayr J, Sugiyama H, Mikkola A. Review on the absolute nodal coordinate formulation for large deformation analysis of multibody systems. J Comput Nonlinear Dyn. 2013;8(3):031016.10.1115/1.4023487Search in Google Scholar

[32] Omar MA, Shabana AA. A two-dimensional shear deformable beam for large rotation and deformation problems. J Sound Vib. 2001;243(3):565–76.10.1006/jsvi.2000.3416Search in Google Scholar

[33] Shabana AA. Computational continuum mechanics. 3rd edn. United States: John Wiley & Sons; 2017.10.1002/9781119293248Search in Google Scholar

[34] Patel M, Shabana AA. Locking alleviation in the large displacement analysis of beam elements: the strain split method. Acta Mech. 2018;229(7):2923–46.10.1007/s00707-018-2131-5Search in Google Scholar

[35] Gerstmayr J, Matikainen MK, Mikkola AM. A geometrically exact beam element based on the absolute nodal coordinate formulation. Multibody Syst Dyn. 2008;20(4):359–84.10.1007/s11044-008-9125-3Search in Google Scholar

[36] Crisfield MA. A fast incremental/iterative solution procedure that handles “Snap-through.”. Comput Struct. 1981;13(1–3):55–62.10.1016/0045-7949(81)90108-5Search in Google Scholar

[37] Borst R, de Crisfield MA, Remmers JJC, Verhoosel CV. Non-linear finite element analysis of solids and structures. John Wiley & Sons; 2012.Search in Google Scholar

[38] Shaukat AR, Lan P, Wang J, Wang T. In-plane nonlinear postbuckling analysis of circular arches using absolute nodal coordinate formulation with arc-length method. Proceedings of the Institution of Mechanical Engineers Part K Journal of Multi-body Dynamics. vol. 253, United Kingdom. 2021, p. 297–311.10.1177/1464419320971412Search in Google Scholar

[39] Wang J, Wang T. Buckling analysis of beam structure with absolute nodal coordinate formulation. Proc Inst Mech Eng Part C J Mech Eng Sci. 2021;235(9):1585–92.10.1177/0954406220947117Search in Google Scholar

Received: 2025-05-17

Revised: 2025-06-23

Accepted: 2025-07-11

Published Online: 2025-10-06

This work is licensed under the Creative Commons Attribution 4.0 International License.

Articles in the same Issue

https://doi.org/10.1515/secm-2025-0063

Keywords for this article

nonlinear buckling; nonlinear postbuckling; large deformation analysis; axially functionally graded material; absolute nodal coordinate formulation

Creative Commons

BY 4.0