Abstract

This paper investigates different approximation techniques for planar beam column elements in force-based methods. The three fields, introduced in this review, are: curvature-based displacement interpolation (CBDI) used in matrix-based flexibility formulations, linear displacement approximation applied in state space, and higher-order displacement approximation utilized again in state space. Using these three approximation fields, the responses and their accuracies in some systems are compared in examples. Finally, focusing on the accuracy and regarding the performed analyses, it seems that the computational cost is reduced and accuracy of responses is elevated in many engineering problems using the higher-order approximation field in state space.

1. Introduction

Traditionally, displacement-based methods are implemented in structural analyses. These methods have some drawbacks caused by their displacement field interpolations, particularly, in cases of high nonlinearities and nonprismatic beam sections. In this regard, equilibrium-based elements have been recently introduced to conquer the inefficient solution processes in displacement-based approaches.

There are a few researchers who have developed and worked on force-based techniques. Ciampi and Carlesimo [1] developed an equilibrium-based element for structural seismic analyses. Spacone et al. [2, 3] as well as Petrangelti et al. [4] developed the mentioned element to more material nonlinearity effects in concrete elements.

Neuenhofer and Filippou [5] presented a force-based element for geometrically nonlinear analyses of plane frame structures; Taylor et al. [6] presented a mixed finite element method for beam-column problems.

Generally, there are two main aspects of element state determinations, in equilibrium-based methods: the matrix-based approach, used primarily by Spacone et al. [2, 3] and Taucer et al. [7]; the dynamical state space approach recently applied by Simeonov et al. [8] and then Sivaselvan and Reinhorn [9]. The former method is very common and straight forward in a finite element program. On the other hand, the later technique, called State Space Approach (SSA), is based on direct solution of governing differential Algebraic equations (DAEs) of a structural system. Therefore, the basic solution method differs from finite element formulations, which are developed based on structural matrix analyses.

It is noted that the accuracy of force-based elements is directly affected by their interelement displacement field approximation, which is used to evaluate the element nonlinear flexibility matrices. Three main approximation procedures have been proposed till now: displacement interpolation (CBDI) technique, linear field in state space approach, higher-order field in state-space approach. The first two approaches are introduced for Euler-Bernoulli straight beam-column elements; however, the third one is capable of analyzing Timoshenko curved frame elements.

The first technique is primarily developed by Neuenhofer and Filippou [5] for a geometrical nonlinear Euler-Bernoulli force-based beam-column element. Magalhaesde Souza [10] extended this formulation to inelastic spatial frame analyses. Lagrange polynomials are used in their technique to interpolate curvatures at integration points for displacement field approximations. Because of the fact that rotations are assumed small in this procedure, large deformation problems need more elements in this technique.

The second one was at first proposed by Sivaselvan and Reinhorn [9]. According to this technique, a linear interpolation is used to evaluate the effects of large internal rotations. Although building structural systems without high nonlinearities are effectively analyzed by their approach, no reliable responses are obtained in very large displacements and deformations.

The third approach is presented by Jafari et al. [11] and offers a higher-order displacement field approximation to be applied to the governing differential algebraic equations upon Lagrange polynomials in state space. The procedure is applied to planar elastic Timoshenko curved beams considering large displacements and interelement rotations.

This paper reviews firstly the basic concepts of a general Timoshenko force-based planar curved element. Afterwards, three mentioned approximation procedures are explained. A brief comparison is performed using different approximation procedures applied for equilibrium-based elements till now. Finally, the numerical examples compare the results obtained by these different approximation techniques.

2. State Equations

An initially curved frame is considered without any internal loading. The end displacement and force vector components of this element are shown in Figure 1(a) in the global planar coordinate system. The adopted local coordinate system is shown in Figure 1(b) based on the corotational transformation technique. Three independent local degrees of freedom, demonstrated in Figure 1(b), are𝐩=𝑃1𝑃2𝑃3𝑇=𝑁𝑀𝐼𝑀𝐽𝑇𝜹=𝐷1𝐷2𝐷3𝑇=𝑢𝐼𝜃𝐼𝜃𝐽𝑇,(2.1) where 𝑁 = axial load; 𝑀𝐼 and 𝑀𝐽 = two corotational moments; 𝑢𝐼 = axial displacement; 𝜃𝐼 and 𝜃𝐽 = two end rotations.

In order to evaluate the element flexibility matrix, the equilibrium and compatibility equations should be presented for curved beams. The Reissner’s kinematical exact stress resultant theory [12] is applied for curved beams to obtain equilibrium equations for geometrical nonlinear beam with end forces by Jafari et al. [13, 14]. Hence, we have the following

(a) Equilibrium Equations
One has 𝐬𝑠0=𝑁𝑠0𝑀𝑠0𝑉𝑠0𝑠=𝚪0𝑃1𝑃2𝑃3𝑠=𝚪0𝐩,(2.2) where 𝜃𝑠𝚪(𝑠)=cos0𝜂𝜃𝑠sin0𝜃𝑠sin0𝜃𝑠/𝑙𝜉/𝑙1cos0𝜃𝑠/𝑙sin0𝜃𝑠/𝑙𝜉/𝑙cos0/𝑙(2.3) and  𝑙=𝐿+𝐷1, shown in Figure 1(b).
Here, a special point is assumed with local coordinates of (𝑥,𝑦) on the element centre line moving to (𝜉,𝜂) in the deformed configuration.
𝜃 is angle between beam cross-section and 𝑦-axis;
𝐬(𝑠0) is vector of arbitrary section force components, Figure 1(c); 𝑁(𝑠0)is axial force; 𝑀(𝑠0)is moment; 𝑉(𝑠0)is shear.

(b) Compatibility Equations
Compatibility equations in state space are:
Finally, the matrix form of compatibility equation is: ̇𝜹=𝑠0𝚪𝐓̇𝐝𝑠0𝑑𝑠0,(2.4) where 𝐝(𝑠0)= vector of generalized strains as 𝐝𝑠0=𝜀0𝑠0𝑠𝜃0𝜃0𝑠0𝛾𝑠0𝑇.(2.5)

(c) Element State Equations
The element state equation is 𝑠0𝚪𝑠0𝑇𝚽𝑠0𝚪𝑠0𝑑𝑠0̇̇𝐩=𝜹𝑠0𝚪𝑠0𝑇𝚽𝑠0̇𝚪𝑠0𝑑𝑠0𝐩,(2.6) where 𝚽𝑠0=𝐸𝑇𝐴000𝐸𝑇𝐼000𝐺𝑇𝐴𝑠1,(2.7) where Φ(𝑠0) is section flexibility matrix; 𝐸𝑇 is axial material tangent modulus; 𝐺𝑇 is shear material tangent modulus; 𝐴 is section area; 𝐼 is inertia; 𝐴𝑠is shear area.

3. Displacement Field Approximations in Force-Based Methods

An interelement displacement field is required to evaluate the flexibility matrices in (2.6) matrix for geometrical nonlinear state. Besides, the element accuracy is extensively related to displacements approximation procedures. So, the existed techniques are explained hereafter.

3.1. CBDI

Neuenhofer and Filippou developed a geometrically force-based element in matrix-based formulation. The flexibility matrix in their formulas is a simple form of (2.6) and requires a field approximation; it consists of displacement terms. It is noted that a straight Euler- Bernoulli beam is considered in this formulation. They named their technique, curvature-based displacement interpolation, “CBDI.”

Here, transverse displacements are required in the flexibility matrix, and rotations are assumed small. So, CBDI technique is used only for evaluating the transverse displacement field. This field is obtained using curvatures at integration points for planar case by 𝐯=𝐋𝜿𝐋=𝐿212𝑋21𝑋116𝑋31𝑋11𝑋NG(NG+1)NG1+1𝑋112𝑋2NG𝑋NG16𝑋3NG𝑋NG1𝑋NG(NG+1)NG+1NG𝑋NG𝐆1,𝐆=1𝑋1𝑋21𝑋NG111𝑋NG𝑋2NG𝑋NG1NG,(3.1) where 𝐯 is transverse displacement field vector; 𝐋 is influence matrix including polynomials; 𝐆 is so called Vandermnode matrix; 𝑋𝑖 is coordinate of integration point 𝑖 in the local system; 𝑖 is changes from 1 to NG; NGis Number of integration points; 𝜿 is curvature vector at integration points.

Although this approximation performs well in many analyses, it lacks in accuracy where large internal rotations happen in an analysis. This is because of the fact that the rotations have been neglected in this formulation.

Jafari et al. [14] extended CBDI to curvature/shearing-based displacement interpolation (CSBDI) to capture shear effects in force-based elements.

3.2. Linear Displacement Field Approximation in State Space

As it is seen in (2.6), the displacement field includes 𝜉,𝜂, and 𝜃 for the element. Using a rate form equation based on the strain displacement relationship (Huddleston [15]), a linear approximation procedure was proposed by Sivaselvan and Reinhorn [9] in their analyses.

As mentioned in matrix-based formulation, shear effects as well as initial curvature effects are not considered in this method. However, its extension to shear deformable curved element is very straightforward as described in Section 2. The strain-displacement relations presented by Huddleston [15] are 𝑑𝜉=𝑑𝑥1+𝜀0(𝑥)cos(𝜃(𝑥)),𝑑𝜂=𝑑𝑥1+𝜀0(𝑥)sin(𝜃(𝑥)).(3.2) In this procedure, the deformations and local coordinates are monitored at NG integration points in terms of NG parameters as 𝜃𝑗𝑥=𝜽𝑗,𝜀𝑗=𝜺0𝑥𝑗,𝜉𝑖𝑥=𝝃𝑗,𝜂𝑖𝑥=𝜼𝑗,(3.3) where 𝜽(𝑥𝑗) is the 𝑗th element of rotation derivative vector; 𝜺𝟎(𝑥𝑗) is the 𝑗th element of axial strain vector; 𝝃(𝑥𝑗) is the 𝑗th element of 𝜉 local coordinates vector; 𝜼(𝑥𝑗) is the 𝑗th of 𝜂 local coordinates vector.

Here, an implicit first-order method is used to integrate (3.2) as initial value problems (IVP), for two times from two ends. Eventually, a weighted average is taken for approximation. For this purpose the formulas are mentioned bellow, starting from one end of the element 𝜃𝑖+1=𝜃𝑖+𝑋𝑖+1𝑋𝑖𝜃𝑖+1,𝜉𝑖+1=𝜉𝑖+𝑋𝑖+1𝑋𝑖1+𝜀𝑖+1𝜃cos𝑖+1,𝜂𝑖+1=𝜂𝑖+𝑋𝑖+1𝑋𝑖1+𝜀𝑖+1𝜃sin𝑖+1.(3.4)

This time, the equations are integrated as IVPs starting from the other end of element. Finally, displacement fields are evaluated at every integration points using weighted average of two solutions.

Due to the fact that this formulation counts for large internal rotations, its accuracy in analyses with large deformations is well elevated in comparison with pervious formulations. Therefore, a linear field in approximation can result in acceptable responses just in building structures where not very high geometric nonlinearities happen. In order to improve this accuracy and obtain reliable results, a large number of integration points and/or more elements are required. Consequently, the computational cost is considerably increased.

In addition, using Huddleston strain-displacement relationship leads to nonlinear rotational terms, and consequently more accurate results are predicted. However, the accuracy is obviously weakened applying this technique; that is, how to approximate the field is of high importance. The deficiency of linear approximation is confirmed in many examples in comparison with CBDI and CSBDI (Jafari et al. [14]).

3.3. Higher-Order Displacement Field Approximation in State Space

As mentioned previously, the displacement field components 𝜉, 𝜂, and 𝜃 need to be evaluated at integration points. Here, a higher-order displacement approximation method is adjusted using Lagrange polynomials.

The method proposed by Jafari et al. [11] is mainly composed of two techniques, linear approximation and CBDI. In this formulation, again the deformations are monitored at NG integration points in terms of NG parameters as 𝜃𝑗𝑠=𝜽𝑗,𝛾𝑗𝑠=𝜸𝑗,𝜀𝑗=𝜺0𝑠𝑗,𝜉𝑖𝑠=𝝃𝑗,𝜂𝑖𝑠=𝜼𝑗,(3.5) where 𝑠𝑗 is local coordinates of integration points on curved beam; 𝜸(𝑠𝑗) is the 𝑗th element of shear strain vector.

According to above assumption, the bellow formulas can be written as 𝜀(𝑠)=NG𝑗=1𝑙𝑗(𝑠)𝜀𝑗,𝜃(𝑠)=NG𝑗=1𝑙𝑗(𝑠)𝜃𝑗,𝛾(𝑠)=NG𝑗=1𝑙𝑗(𝑠)𝛾𝑗,(3.6) where, 𝑙𝑗(𝑠)is Lagrangian polynomials, obtained from (3.7).

We have 𝑙𝑗(𝑠)=NG𝑖=1,𝑖𝑗𝑠𝑠0𝑖NG𝑖=1,𝑖𝑗𝑠0𝑗𝑠0𝑖,(3.7) thus, 𝑙𝑗𝑠𝑖=𝛿𝑖𝑗.(3.8)

Regarding the necessity of integrating over functions, the more applicable form of Lagrange polynomials is used for closed form integration as 𝑙1𝑙(𝑠)2𝑙(𝑠)NG=1𝑠𝑠(𝑠)NG1𝐆1.(3.9) This matrix is determined only one time as the number of integration points remains constant during analysis. Starting from one end of the element, the incremental form of internal displacement field components can be stated as𝜃𝑖+1=𝜃𝑖+Δ𝜃𝑖,𝑖+1,(3.10)𝜉𝑖+1=𝜉𝑖+Δ𝜉𝑖,𝑖+1,(3.11)𝜂𝑖+1=𝜂𝑖+Δ𝜂𝑖,𝑖+1.(3.12) The incremental values are written, using a higher-order and more accurate integration, asΔ𝜃𝑖𝑗=𝑆0𝑠𝑠22𝑠23𝑠NGNG𝐆1||||||𝜽𝑠=𝑠𝑗𝑠=𝑠𝑖,(3.13)Δ𝜉𝑖𝑗=𝑠𝑗𝑠𝑖1+𝜀0𝑠0𝜃𝑠cos0𝜃𝑠𝛾sin0𝑑𝑠0,(3.14)Δ𝜂𝑖𝑗=𝑠𝑗𝑠𝑖1+𝜀0𝑠0𝜃𝑠sin0𝑠+𝛾0𝜃𝑠cos0𝑑𝑠0.(3.15)

These integrations are achieved based on the local Gauss or Gauss-Lobatto procedure performed between 𝑠=𝑠𝑖 and 𝑠=𝑠𝑗.

In this way, the displacement field is obtained from one end of the element. The same procedure is repeated from other end performing the operations as (3.10) and (3.15). Finally, the weighted averages of these two approximations are evaluated, and three displacement field components are obtained. In this study, the weights are considered the same for 𝜉, 𝜂, and 𝜃, in every integration points.

As it is demonstrated next, the results from this approximation technique are more accurate in comparison with those of other procedures. Besides, the shear and curvilinear effects are considered in the analyses.

It is seen from different analyses that there is no significant change in computational cost using this method, in comparison with the previous ones. Regarding the former, Lagrangian polynomials are used in both CBDI method and the one presented here. However, large deformations, including shear and curvilinear effects, are considered by this new method, and, consequently, more accurate results are obtained.

In linear approximation, despite considering large deformations for Euler- Bernoulli beam, a large number of integration points are needed in order to meet accurate results. The computational cost is increased significantly as the consequence of this phenomenon.

4. Numerical Examples

4.1. Williams Toggle Frame (Figure 2)

This example is to confirm the accuracy of the results obtained by this new displacement field approximation in comparison with linear field one. Toggle frame has analyzed elastically by several researchers for surveying geometrical nonlinear formulations; the frame and its properties are shown in Figure 3.

In Figure 4, the equilibrium curve is evaluated using one element per half of the structure with four Lobatto points inside. The results, presented by the formulation with higher-order approximation, are more reliable and more efficient comparing with linear method in force-based technique. It should be noted that ten Lobatto points per element are used in linear approximation.

The results of CBDI technique are approximately similar to those of the proposed method as large internal rotations are not included in the problem. As it is seen in Figure 3, the accuracy is not significantly changed by increasing the number of internal integration points; therefore, four Lobatto points seem enough.

4.2. Postbuckling Analysis of a Clamped-Hinged Circular Arch

This well-known example presents the geometric nonlinear analysis of the pre- and Postbuckling deformation of a circular arch, hinged at one end and clamped at the other end, under a vertical force applied at the apex. This problem was suggested by da Deppo and Schmidt [16] and has been analyzed by many researchers. (Simo and Vu-Quoc [17]; Ibrahimbegovic, and Frey [18]).

Material and geometric data for the arch are shown in Figure 4. As responses show, the problem includes a sharp negative slope after the first buckling load is reached. The use of higher-order approximation field as well as continuous form of constraint equation results in robust responses as are demonstrated in Figure 5 for eight elements and five Lobatto points.

The use of arc length solution technique assists to trace the complete equilibrium path which includes a steep downward slope. For this problem, ten elements with five internal Lobatto points result in accurate responses.

5. Conclusions

Three main approximation fields are reviewed: CBDI in matrix-based formulations, linear approximation in state space, and Lagrange polynomial-based displacement field approximation for curved elements in state space. The later field is basically resulted from two former methods. While fewer elements and integration points are used in the latter method, the results gained for high nonlinear problems are more accurate and reliable. Besides, the efficiency of presented method is not affected, comparing to linear approximation procedure, because of the fact that applying lower number of integration points results in reliably accurate responses. The approximating polynomials are not generally high ordered; therefore, nonphysical wriggling is avoided. However, it seems that using an approximation like Bsplines may enhance the efficiency even better.