Analytical determination of the vibration frequencies and buckling loads of slender reinforced concrete towers

This study focused on improving the design of slender structures with reinforced concrete (RC) telecommunication towers as the main application. Analytical procedure based on Rayleigh’s method to compute the first natural vibration frequency and the critical buckling load was development. All the nonlinearities present in the system were considered, in addition to the soil–structure interaction and the variation of the geometric properties along the length of the structure. The geometric nonlinearity and imperfections of the tower structure were computed as functions of the axial load using a geometric stiffness matrix. Further, the material nonlinearity was accounted for by reducing the flexural stiffness. As concrete structures exhibit viscoelasticity, creep was calculated using the Eurocode 2 model. The soil– structure interaction was modeled as a set of distributed springs. To validate the proposed method, the first frequency and critical buckling load were compared with those yielded by FEM simulations. The frequency results were in good agreement with those of the FEM simulations, indicating that the proposed method is sufficiently accurate for use in engineering design applications and easy to implement. On the other hand, the buckling load results obtained using the proposed method and FEM differed significantly, motivating further investigation.


INTRODUCTION
Telecommunication towers are slender beam-column structures that are frequently subjected to dynamic loading due to wind excitation.The wind speed is usually analyzed by considering two components: the mean wind speed and the floating wind speed (Eurocode Standard, 2010).The mean wind speed induces static loads, whereas the floating wind speed is responsible for the vibrations of such structures.According to Silva et al. (2013) and Silva and Brasil (2006), the dynamic response depends mainly on the first vibration frequency and the height of the structure.Hence, precise prediction of the first resonance frequency is necessary for dynamic analysis of such structures.As, for these cases, the frequencies are especially very low, suggesting proximity to the losing of stability, it is necessary, additionally, to verify the structural safe in relation to the vertical capacity of loading, understood it as the limit value for installation of antennas and platforms.Considering both aspects, the main objectives of this study were to develop a novel analytical expression for computing the first vibration frequency and modal shape of slender reinforced concrete (RC) structures, to compute the critical buckling load considering geometrical and physical nonlinearities, and to compare the obtained values with those given by the finite element method (FEM), what is going to allow analysis of similar structural systems when exposed to the environmental wind.
In relation to the determination of the lateral forces produced by the wind, Vieira, Brito and Loredo-Souza (2018) point that it is an extreme complex problem because of the large number of variables involved.Amiri and Bucher (2017), affirmed that the design codes for wind loading provide useful instructions for engineers, but wind tunnel testing, which is generally employed for their evaluation, has several limitations.Consequently, field measurements can be very useful for verifying wind tunnel test data, allowing more realistic and reliable analysis.On the other hand, the comparison of numerical simulations and practical field results provides useful information for computing modal wind loads properly.Lateral forces due to wind produces lateral displacement on a structure.This configures an important aspect of the structural analysis for telecommunication systems.In terms of the serviceability state limit constraints to deflections at position of antennas are imposed in order to avoid discontinuity on communication service, transmission/reception of signal of cellular mobile phone.These deflections obviously depend on the lateral forces due to the wind, which in turn depend on the structural frequency (flexible or rigid body).Therein lies the importance of a correct calculation, which must consider all nonlinear aspects present in the problem.When calculating a slender concrete element must be observed that its frequency change on time because of the viscoelastic behavior of the material.Under wind action, if flexible, a structure vibrates around its balance position and presents displacement that can drive to the loss of signal.In this context, the permanent challenge for the mobile telecommunication operators is the establishment of a distribution system with a network configuration that be able to operate continuously to keep the transmission of the radio signal, ensuring local and long-distance coverage.
Taking into consideration the complexity characteristic of a phenomena of environmental nature, Durbey and Hansen (1996) presented several methods for computing the wind loads on structures.In addition, they stated that dynamic wind excitation can cause structures to resonate.Silva and Brasil (2006) and Simiu and Scanlan (1996) showed that for structures with low first frequencies, i.e., less than 1 Hz, the dynamic effects of wind speed cannot be ignored, and the classical static analysis does not properly capture the structural behavior.They also noted that the fundamental frequency is the frequency that contributes the most to the dynamic behavior of the structure.This is general principle, which is applied either for a telecommunication tower as a building.Lo et al. (2017) who performed vibration tests of high-rise buildings mention that someone of factors have been investigated for many researchers in relation to wind action over these structures are their modal frequency and mode shapes.In this sense Pozos-Estrada and Hong (2015) point out that the value that most of the structures are sensitive to dynamic effect of the wind, especially tall buildings, falls within 0.05 to 0.5 (Hz).In that context is important to have clear that a telecommunication tower is a flexible and very slender structural element.Kwak and Kim (2004) affirmed that the failure of an RC column having a high slenderness ratio occurs well within the cross-section interaction diagram because of a pronounced second-order effect.In the case of columns with a slenderness ratio that does not exceed 100, failure occurs rather close to the material strength.The wind-induced collapse of three tall hyperbolic RC cooling towers in Ferrybridge, UK, in 1965, was investigated by Jia (2013), who concluded that buckling, in the classical sense of the term, was the most likely cause of the accident.
The compressive capacity of a column can be quantified in several ways, one of which is based on the Euler buckling load (Uzny, 2011).The elastic buckling load of a pin-ended column was first determined by Euler in 1774 (El-Sawy, Sweedan and Martini, 2009).However, he found it difficult to consider the effects of the self-weight of a column when determining its buckling response.Although this problem was finally solved by Greenhill in 1881 (Timoshenko and Gere, 1961), consideration of the column self-weight remains a subject of intense and extensive mathematical discussion (Lubbers et al., 2017), and several contributions have been made to the field since then.In the present study, the influence of the self-weight was considered when calculating the critical buckling load, which was computed by imposing a null value for the frequency.This dynamic approach was previously adopted by Deng et al. (2017) to calculate the stability of a multi-span beam used to represent a fluid-conveying pipe.In the present study, three components were considered to describe the elastic behavior of the structural system: the structure material and its nonlinear behavior, the geometric nonlinearity, and the soil.
In relation to the material nonlinearity of the concrete Liu et al. (2016) developed a mathematical model to investigate the influence of cracks on the dynamic behavior of a cantilever beam.For this purpose, they developed a model with local flexibility at the crack tips to analyze the beam vibrations.Challamel et al. (2005) performed several numerical simulations to study the effects of creep and the strain rate on the failure stresses of quasi-brittle materials.
In the present study, the effects of the material nonlinearities were captured by reducing the flexural stiffness, while the creep effects were computed in accordance with the relevant Eurocode 2 model ( 2004).
The geometric nonlinearity whether produced by the eccentricity of the load applied to a column or a constructive imperfection can be considered as reducing the stiffness of the structure in terms of its geometric stiffness, as recommended by Wilson and Habibullah (1987), because it is a function of existent axial forces.The effect of the axial load on the fundamental vibration frequency was further examined by Wahrhaftig and Brasil (2016), who studied a cantilever beam with large initial displacement.The geometric stiffness is independent of the elastic properties of the material, as explained by Noh (2006).Houmat (2015) applied the concept of geometric stiffness to the free vibration of laminates with variable stiffness along the structure domain.Further, Challamel and Hellesland (2013) performed buckling and post-buckling analysis of softening columns composed of quasi-brittle materials, specifically, RC columns, by applying eccentric vertical loads.They used continuum damage mechanics to include the influence of micro cracking on the material.Zuo, Bi, and Hao (2018) noted that the soil-structure interaction could not be ignored, because, particularly in the case of a monopile foundation, there is an important interaction with and a significant influence on the dynamic behavior of the structure.Amini et al. (2018) investigated this phenomenon associated with the soilstructure interaction.They asserted that a major challenge when modeling structures under wind loading is proper design of the foundations.
Owing to the large number of factors that influence the structures analyzed in this study, whose behavior is governed by partial non-linear differential equations, the exact response cannot be computed analytically, affirm Awrejcewicz et al. (2014).Challamel et al. (2015a) noted that, in numerical approaches, unsolvable mathematical equations of continuous systems are often converted into finite numbers of variables associated with discrete equivalent systems to reduce the complexity of the problem.According to Challamel et al. (2015b), the possibility of solving complex continuous systems using a computer makes the discretization approach particularly advantageous over continuous methods.Alternatively, a continuous system can be reduced to a single-degree-of-freedom (SDOF) system.In this approach, the system is restricted to have only one shape mode, and the system properties can be expressed as functions of a single generalized coordinate.This technique was proposed by Lord Rayleigh (1877) in his work on elastic system vibrations.Most actual structures are complex systems with properties that vary along their lengths.In such cases, the integrals obtained using Rayleigh's method must be solved within the limits established for each interval, i.e., generalized properties can be obtained for each discrete segment of the structure.Whereas Rayleigh's method involves only scalar numbers, in the FEM, it is necessary to deal with matrices with the order given by the number of degrees of freedom considered in the mathematical model.

Overview
Creep is the gradual increase of structural deformation with time, and it is a typical occurrence in concrete structures owing to the viscoelastic nature of concrete.It is important to consider creep in the ultimate limiting state when verifying the stability of a slender member under compression, because the stiffness of such a member varies with the rheology of the material.Charpin and Sanahuja (2017) noted that understanding concrete behavior under extended sustained loading had emerged as an important topic since the study performed by Granger (1995) Granger (1995) further observed that when stress remains constant over time, deformation can be calculated classically by replacing the Young's modulus with a deferred or temporal modulus of elasticity.Boshoff and van Zijl (2007) also used a temporal modulus and proposed an expression for the long-term material stiffness by simply adjusting the Young's modulus using a creep coefficient.Hamed and Bradford (2010) performed a similar analysis in their mathematical investigation of creep in RC beams strengthened with fiber-reinforced polymers.Han et al. (2017a) experimentally investigated creep in a concrete-filled steel tube.In addition, Irfan-ul-Hassan et al. (2016) established the relationship between the modulus of elasticity and time and quantified the Young's modulus by performing an experimental examination of cylindrical specimens.
To confirm the occurrence of creep in concrete structures subjected to tensile and compressive loads, Rossi, Tailhan, and Le Maou ( 2013) conducted an experimental study that revealed the apparent equivalence of compressive and tensile creep under curing conditions.Ranaivomanana, Multon, and Turatsinze (2013) also performed an experimental investigation using different levels of stress, and they concluded that nonlinear stress behavior occurs under compression that exceeds 30%-50% of the concrete compressive strength.In the present study, the compressive stress was varied from 0.14% of the compressive strength at the tower tip to 3.62% of that at the base.Assmann and Reinhardt (2014) used concrete specimens that included polymers to experimentally investigate creep under tension.Further, Charpin and Sanahuja (2017) calculated the Poisson's ratio of concrete using classical integral expressions of the linear viscoelastic behavior of concrete, thereby assuming isotropic linear viscoelastic creep behavior in a concrete structure.
Mathematically, creep can be characterized by an increase in the initial elastic deformation through viscous deformation, resulting in deformation as a function of time.Consequently, the modulus of elasticity must also be represented as a temporal function.Theoretically, creep strain can be calculated under arbitrary load conditions by using the superposition principle.However, it is difficult to distinguish between the effects of creep and those of elastic strain during initial loading, as reported by Su et al. (2017).This difficulty arises from the fact that it is not possible to distinguish between the beginning of immediate deformation and long-term deformation, as immediate deformation already includes a certain portion of the creep deformation.Qian and Kawashima (2016) investigated the behavior of concrete specimens and observed that the viscosity increased when the applied stress was relatively low, causing the material to yield and viscoelastic solid-like behavior to appear.Rheological models for evaluating creep were developed by Sellier et al. (2016), who performed a series of computational simulations in an FEM environment.A viscoelasticplastic model that superposes the strains obtained using the Maxwell, Kelvin, Bingham, and Saint-Venant models was proposed by Han et al. (2017b) to describe early-age creep in concrete.This model essentially represents a modified version of a system that maintains the conventional initial viscoelastic behavior of concrete.In the field of structural analysis, Pascon (2018) employed the well-known three-parameter model, normally used to represent the behavior of a solid, to numerically evaluate the time-dependent mechanical behavior and large deformation of functionally graded visco-hyperelastic materials.A similar strategy of using a set of a spring and a damper for representing the rheological behavior of a material was adopted by Milašinovic´ and Landovic (2018) for studying viscoelastoplastic vibrations in the analysis of low cycle fatigue in internally damped inelastic frame structures.A more complex rheological model was recently developed and experimentally validated by Delsaute, Torrenti, and Staquet (2017).Awruch and Gomes (2015) suggested that the enduring effects of loading a structure originate from the timedependent deformations that occur during the first few months of loading but last for most of the service life of the structure.This characteristic sometimes results in displacement of the same order as that of the instantaneous deformation induced by a service load.Mehta, Monteiro and Filho (1994) noted that the slow deformation of a material depends on parameters such as the loading direction, speed, and duration; material moisture content; ambient temperature; and relative humidity.The rheological properties of concrete that are essential to understanding creep behavior were recently discussed in detail by Jiao et al. (2017).
The effects of creep on the dynamic behavior of structures have been investigated by several researchers.On the basis of CEB-FIP 1990 (Comité Européen du Béton (CEB)-Fédération Internationale du Béton (FIB) or International Federation for Structural Concrete), which includes the recommendations of the European Standard (2004) (EN 1992-1-1 -Eurocode 2: Design of concrete structures -Part 1: General rules and rules for buildings) Yang, Chen, Zhang, and Fan (2016) investigated the creep induced by the response of a prestressed concrete railway bridge using a temporal function and found that the phenomenon and its associated deformation rate are influenced by temperature.Temperature is particularly crucial for the operation of high-speed trains because irregularities caused by temperature-induced creep and deformation of the railway concrete affect the wheel unloading rate and vertical acceleration of the train, affirms Wang (2014).Suspension bridges constitute a different type of slender structure, and their creep behavior accounts for a significant portion of the required structural analysis.Arco et al. (2001) observed that the creep response of a suspension bridge and its analysis were complicated by the combined effects of the geometric nonlinearity, pretension, and time-variant material behavior.

Creep as addressed in Eurocode 2
As noted by Hołowaty (2015), the creep model in the European Standard EN 1992 -1-1 (Eurocode 2, 2004) is based on the CEB-FIP Model Code (1990).In addition, the method specified in Eurocode 2 for incorporating creep into structural analysis considers the effects of the creep behavior and its variation with time.Eurocode 2 provides hypothetical and model limitations for creep calculation, wherein the creep coefficient φ is predicted as a function of the tangent modulus of elasticity Ec, which can be considered as the secant modulus of elasticity Ecm + 5%.The creep deformation of concrete is computed by multiplying the immediate deformation by the creep coefficient.The total concrete deformation at time t, under constant temperature, can be obtained as the sum of the terms that represent the immediate deformation, creep, and shrinkage, as defined by Neville (1997).The method specified in EN 1992-1-1 ( 2004) makes no distinction between slow deformation and fast-reversible, slow-irreversible, or slow-reversible deformation.All the factors of the concrete, load, and environment are calculated under the assumption that they remain constant over the considered time interval, affording a specific result for the creep coefficient φ.This coefficient is then directly introduced into the slow deformation equation and used as input data for various procedures.The basic equations for determining the creep coefficient of concrete over time are based on the average compressive strength fcm (fcm = fck + 8, fck in MPa).The creep coefficient φ(t, t0), as defined in Eq. ( 1), is the product of two macrofactors, namely φ0 and βc(t, t0), which respectively characterize the effects of the rheological properties of the concrete under environmental conditions and the evolution of creep with time after loading of the structure.
The first macrofactor φ0 defined in Eq. ( 2) consists of three microfactors.The first microfactor φRH (given by Eq. ( 3) for average concrete compressive strengths ≤ 35 MPa and by Eq.( 4) for greater than 35 MPa) accounts for the effects of the environmental relative humidity RH, equivalent thickness h0 of the member (a function of the cross-sectional area Ac), and external perimeter ue of the member in contact with the environment.The second microfactor β(fcm) (Eq.( 6)) represents the direct effect of the resistance on φ0.The third microfactor β(t0) (Eq.( 7)) represents the age of the concrete at the beginning of loading, i.e., at t0.
The second macrofactor βc(t, t0) (Eq.( 8)) is a function of the coefficient βH (given by Eq. ( 9) for average concrete compressive strengths ≤ 35 MPa and by Eq.( 10) for upper), and it is used to regulate the combined effects of the relative humidity and the equivalent member thickness.The percolation path of the adsorbed water in a robust section of concrete is so large that the effects of creep due to differential moisture are less important for slimmer sections.
Thus, the creep coefficient can be obtained using Eq. ( 1), and the temporal function that describes the deformation in accordance with EN 1992-1-1 can be expressed as Based on the above equations, the modulus of elasticity with respect to time can be expressed as where Ec(t0) is the modulus of elasticity at the beginning of loading and Ec(t28) is the modulus of elasticity 28 days after the commencement of loading.

Description of the practical problem
The practical problem considered in the present study involved the calculation of the first natural vibration frequency and the critical buckling load of an actual slender RC tower (or pole), with its inherent geometrical and material nonlinearities (see Figure 1).Such RC towers are commonly used to support antenna systems for cellular telephony services.The geometrical details of the considered tower are shown in Figure 1(b).The geometric division of the structure and the cross-sections along the tower are shown in Figure 2, where s denotes the segment, S denotes the section type, g is the gravitational acceleration, Gr denotes the ground, db is the reinforcement bar diameter, nb is the number of reinforcement bars, and c´ is the reinforcement covering.The tower has a height of 46 m and a hollow circular cross-section.The external diameter and wall thickness of each structural segment of the tower are denoted as D and th (mm), respectively.Figure 2 shows the cross-section along the length of the tower.The slenderness ratio of the tower is 334.A set of antennae and a platform are installed at the top of the structure, contributing 1097.76 kg to the total mass of the structure.Cables and a ladder are installed along the length of the structure, yielding a distributed mass of 40 kg/m.The centrifuged RC pole and foundation shaft have densities of 2600 kg/m 3 and 2500 kg/m 3 , respectively.The foundation of the tower is composed of a relatively deep shaft with a bell diameter and length of 1400 mm and 200 mm, respectively, and a shaft diameter and length of 800 mm and 5800 mm, respectively.In this study, the lateral soil resistance was represented by an elasticity parameter with a magnitude of 2668.93 kN/m 3 , corresponding to the total lateral stiffness.The concrete resistance was considered to be 45 MPa and 20 MPa in the pole and foundation, respectively.The physical nonlinearity of the RC was computed in accordance with the recommendations of the NBR 6118:2014 (2014) -Design of concrete structures -Procedure, from the Brazilian Association for Standardization (ABNT), suggesting a 50% reduction of the gross moment of inertia.Note that Brasil et al. (2007) and Silva et al. (2013) investigated the same problem under similar conditions and confirmed that the effective flexural stiffness was no less than 40% of the total stiffness when the structure was subjected to the designed gross moment of inertia.The details of the tower and installed equipment are summarized in Table 1.As the focus of the present study was the RC structure, it was necessary to account for the contribution of the reinforcement bars to the moment of inertia of the entire cross-section, which was accomplished by homogenization.A circular-ring cross-section was assumed, with external diameter D, wall thickness th (for the pole, the foundation does not have a thickness th because it is a full cross-section), and reinforcement covering c´ (= 25 mm).The modulus of elasticity of the steel reinforcement (Est) is 205 GPa.The number of bars (nb) distributed along the periphery of the cross-section is 20, and each bar has a diameter (db) of 12.7 mm, as stated in Figure 3  The variable Rbi describes the position of the barycenter of each bar bi relative to the barycenter of the crosssection.Because all the bars have the same radius Rbi, the following expression applies: where c´ is the covering and dbi is the diameter of the ith bar.Because θi is an independent variable that takes values between 0 and 2π, the distance between the center of each bar and the axis of the cross-section can be expressed as The spacing sp between the centers of adjacent bars is given by where nb denotes the number of longitudinal bars.The angular phase shift Δθ between two adjacent bars is thus The moment of inertia of each bar about the cross-section barycenter, Isi, can be computed using the following expression: The homogenized moment of inertia of the reinforcement steel, Ihs, can be expressed as The moment of inertia of the concrete in the cross-section (without neglecting the steel domain) for the pole and foundation are, respectively, given by ( ) while the total homogenized moment of inertia of the cross-section is given by Itot = Icon + Ihs.Therefore, the factor Fhs, which is used to multiply the nominal moment of inertia Ihs, computed using the total moment of inertia of the reinforcement steel, can be written as Therefore, the values of the factor Fhs, which multiplies the nominal moment of inertia for each section from 1 to 5, were Fh1 = 1.0199,Fh2 = 1.0568,Fh3 = 1.0811,Fh4 = 1.0671, and Fh5 = 1.0859, respectively.

Parameters of the practical problem
The parameters employed in the analytical and computational process conducted during this study are summarized in Table 2.The modulus of elasticity of the concrete was based on Eurocode 2. In Figure 2, the height of each segment s of the structure is defined from the bottom to the top in the following sequence: s1 extends from 0 to L1, s2 extends from L1 to L2, s3 extends from L2 to L3, s4 extends from L3 to L4, and s5 extends from L4 to L5 (= L).Points L1, L2, L3, L4, and L5 are at elevations of 0.2, 6.0, 12.0, 19.0, and 46.0 m, respectively.
For the base of the foundation, D1 = 1400 mm, x D L Thus, the area and moment of inertia of each segment are given by ( ) x A L and ( ) , respectively.The external diameter and thickness of the initial segment of the superstructure, D3 and th3, are 800 mm and 150 mm, respectively.Hence, the internal diameter, area, and moment of inertia of the segment are given by = , respectively.Therefore, the internal diameter, area, and moment of inertia of the variable section are given by ( ) , respectively.

Overview of the mathematical development of the SDOF model
As observed by Clough and Penzien (1993), the approximation of an infinite-degree-of-freedom system, such as the slim RC tower considered in the present study, requires the use of an SDOF system to formulate the equations of motion for the continuous problem.This SDOF system (a basic cantilever frame element with undamped free movement, as shown in Figure 4) was assumed to consist of a bar possessing either constant or variable properties, such as geometry, elasticity (or even viscoelasticity), and density, represented by Is(x), Es(t), and ( ) s m x , respectively, along its length, defined between ordinates Ls and Ls-1, as shown in Figure 2. A set of springs of variable stiffness kSos(x) was laterally attached to the system to model the soil-structure interactions.In addition to its self-weight and distributed mass, the system supported a mass m0 at its tip.Ns(x) denotes the variable normal generalized force, while u(t) and q(t) denote the vertical and lateral displacements, respectively, of the tip section.The essential properties of the system are its flexural stiffness and mass per unit length.It was further assumed that the motion of the system did not alter the direction of the normal force.To approximate the motion of the system in Figure 4 using an SDOF model, it was assumed that the system deflected in a single manner (or shape) that was adequately chosen to represent the vibration mode considered.The shape function describing this deflection is defined as ϕ(x), where x is an independent variable representing the location along the height of the tower and the amplitude of the motion is represented by the generalized coordinate q(t).This generalized coordinate was conveniently adopted as the displacement of the tip of the tower.The shape function was then obtained by calculating the following dimensionless ratio: The equation of motion of this generalized system could be conveniently formulated using Hamilton's principle.Hence, the external virtual work done by the forces of inertia (WEIN) could be expressed as The internal virtual work done by the flexural deformation (WIFL(t)) is given by and the work done by the deformation of the foundation springs (WIKm) is given by To calculate the work done by the axially directed normal force, whose direction and amplitude remained unchanged during the response, it was necessary to evaluate the vertical component of the motion of the tower tip (u(t)).An infinitesimal element dl of the curved bar was used to determine this axial displacement.The shortening of the axis owing to this axial displacement could be expressed as Taking the binomial expansion and considering the higher-order terms to be smaller than the first-order term, the initial series can be reduced.Thus, Eq. ( 27) can be written in a more compact form as Therefore, the total displacement u(t) along the column can be obtained by integrating Eq. ( 28): Thus, the potential energy of the axial load Ns(x) can be derived from and the external and internal virtual work obeys the following relations: where Latin American Journal of Solids and Structures, 2019, 16(5), e196 15/31 In Eq.( 32) and ( 33), ( ) ′ v x and ( ) ′′ v x represent the first-and second-order derivatives, respectively, of the virtual work done with respect to the independent variable of the lateral displacement, whereas ( ) v t  and ( ) v t  represent the first-and second-order derivatives, respectively, of the virtual work done with respect to the time elapsed during the lateral displacement.The real and virtual displacements and their derivatives, expressed as functions of generalized coordinates, as well as the shape functions used to represent the considered vibration modes, can be calculated using the following set of equations: x q t v x t x q t v x t x q t q t t q t v x t x q t v x t x q t v x t x q t q t t v x t x q t v x t x q t v x t x q Substituting the appropriate terms from Eq. ( 34) into Eq.( 32) yields Because the variation δq(t) can assume any value, this forces the term in brackets to become zero.Therefore, it can be rewritten as where the generalized mass M of the system, including the mass at the tip, is given by ( ) ( ) .
The generalized conventional stiffness K0(t) is given by whereas the geometric stiffness Kg becomes where N(x) denotes a function of the normal force, and the elastic stiffness of the spring (KSo) is given by ( ) ( ) .
The analytical formulation described above is similar to that of Rayleigh [15], who assumed that an infinitedegree-of-freedom system could be replaced by an SDOF system to estimate the frequency and buckling load.For example, Simão (2017) used the principles set forth by Rayleigh to evaluate the critical buckling loads of columns under the influence of shear.Leissa (2005) claimed that the precision achieved by applying this method depends solely on the form of the function used to represent the free vibration mode; if the exact shape is employed, the exact corresponding frequency is generated.Scheble et al. (2003) asserted that the classical Rayleigh-Ritz method requires a priori admissibility of the shape function to satisfy the underlying variational principle.However, it is difficult to satisfy this requirement for a complex structure.Scheble et al. (2003) also noted that although diverse approaches have been used to determine the natural frequencies and normal vibration modes of large structures, they generally involve strategies for dividing structures into assemblages of finite components or substructures.
The underlying concept of both Rayleigh's method and the analytical method proposed in this paper is a simplified form of the principle of energy conservation in mechanical systems.Of some relevance here is the suggestion by Temple and Bickley (1933) that the fundamental principles developed by Rayleigh are applicable to both linear and nonlinear continuous systems with finite or infinite degrees of freedom.The advantage of the method proposed in this paper is its applicability to the determination of the fundamental vibration period and analysis of the stability of elastic systems with sufficient precision for engineering applications.

Mathematical modeling of a practical problem
The system is considered to be under the normal gravitational forces induced by the distributed masses (the selfweight, ladder, and cables) along the length of the tower, as well as the lumped masses at the tip of the tower (the antennae and platform).In addition, the following well-known trigonometric function (Timoshenko and Gere, 1961;Challamel and Hellesland, 2013) is assumed to be valid throughout the domain of the structure: where x represents the location of the calculation, measured upwards from the base along the tower height, and L is the total length (or height) of the column.To obtain an analytical solution to the present problem, the trigonometric function in Eq. ( 41) was considered to be a function of x only, thereby effectively reducing the cantilever column to an SDOF system.Shang et al. (2017) and Pramod et al. (2017) employed trigonometric equations in their analysis to define a shape function for the free vibration of planar elements, having the first one worked in the context of the Unit Partition Method.Ding, Tan, and Dowell (2017); Masjedi and Maheri (2017); and Pradhan and Chakraverty (2017) also used trigonometric equations to account for beam vibrations under different operating conditions.In the case of the vibration of a cantilever column (fixed at its base but free at its tip), the shape function in Eq. ( 41) satisfies the boundary conditions, i.e., when x = L (at the free extremity of the column), ϕ(L) = 1, and when x = 0 (at the base), ϕ(0) = 0.The use of Eq. ( 41) as the shape function for an actual structure with varying geometry was validated by Wahrhaftig (2017), by comparing the results obtained using this function with those obtained using the FEM and other mathematical functions.
Therefore, as described in Section 4.1, the dynamic properties of the system can be computed.The conventional elastic/viscoelastic stiffness can be expressed as shown in Eq. ( 38), where, for a given segment s of the structure, Es(t) is the time-dependent viscoelastic modulus of the material; Is(x) is the variable moment of inertia of a given cross-section along the segment in relation to the considered movement, obtained by interpolating the already homogenized previous and following cross-sections (if the moment of inertia is constant, it is simply Is); k0s(t) is a temporal stiffness term; K0(t) is the final conventional time-varying stiffness; and n is the total number of segments based on the geometry of the structure.In Eq. ( 38), the variable t obviously vanishes for a material with purely elastic timeindependent behavior.
The geometric stiffness, as shown in Eq. ( 39), can be expressed as a function of the axial load, including the selfweight, as follows: where Kg(m0) is the total geometric stiffness of the structure with n segments, N0(m0) is the concentrated force at the top of the structure as a function of the mass m0 at the tip, and Nj is the normal force developed in the upper segments, given by x is the mass per unit length and m is the generalized mass of the system based on the density of the material.If the cross-sectional area of the cantilever column is constant over a given interval, As(x) can be simplified to As, and the distributed mass also becomes constant.Similarly, if m0 does not vary, all the parameters that depend on it become constant.
One way of considering the effects of the soil on the vibration analysis of the system involves regarding the soil as a series of vertically distributed springs that generate a restorative force.If kSos(x) denotes the spring parameter, the effective soil stiffness (as a function of x along the length of the column) can be generally defined as where KSo is an elastic characteristic given by the sum of kSos(x) along the foundation depth Ds(x) (dependent on the geometry of the foundation) and the soil parameter Sops.In this study, Sops was considered to be constant for each layer of soil, which is not necessarily the case.The variables s and n in Eq. ( 46) are as previously defined.
Considering the normal force to be positive, the overall structural stiffness can be obtained by Finally, the natural frequency can be expressed as a function of time and the tip mass as follows: It is important to observe that the conventional stiffness fundamentally depends on the elastic/viscoelastic/plastic behavior of the material, whereas the geometric stiffness is predominantly a function of the normal force acting on the system, including the self-weight of the structural element.All problems are more or less nonlinear from a geometric viewpoint.Hence, Eq. ( 47), which includes Eq. ( 42), automatically captures the level of nonlinearity of the problem in its solution.
The mathematical procedure described above was used in this study to determine not only the frequency of the structure but also the critical buckling load.The relevant criterion for determining the critical buckling load of the dynamic structure was established by assuming zero frequency at the instant at which the structure lost its stiffness.The load was thus expressed as a function of the lumped mass at the free extremity of the column, N (=N(m0), Eq. ( 43)).To introduce creep into the analysis, it was necessary to use a model that accurately represented the viscoelastic behavior of the concrete.After introducing creep, the frequency and critical buckling load become functions of time, owing to the variation of the modulus of elasticity of the structure with time.Therefore, the frequency can be expressed in terms of the elapsed time and the mass at the top of the cantilever column (Eq.( 48)).
For a column with constant elastic/viscoelastic properties along its length, a unique moment of inertia, and a uniform material density, all the previous integrals can be taken over a single interval covering the entire height of the column.Thus, the final equation for calculating the first undamped resonance frequency (Hz) can be expressed in the form of Eq. ( 49), which can be obtained by substituting Eqs. ( 44) and (47) into Eq.( 48).Note that the soil-structure interaction is not included in this formulation.All the previous observations regarding creep remain valid, including the correspondence between the critical buckling load and zero frequency, which occurs at an arbitrary time.
Latin American Journal of Solids and Structures, 2019, 16(5), e196 18/31 ( ) ( , ) where L is the length of the column, g is the acceleration due to gravity, E(t)´ is the viscoelastic modulus of the material, and I´ is the moment of inertia, which must be homogenized and reduced from its full value to account for the presence of steel reinforcements in the RC and the developed cracks within the structure.In Eq. ( 49), the conventional portion of the column stiffness, including the modulus of viscoelasticity, is to the left of the addition sign in the numerator, whereas the geometric portion of the stiffness is to the right.In this equation, the force concentrated at the free end of the column (mass m0) is considered in addition to the self-weight through the parameter ′ m (mass distributed per unit length), and both are multiplied by the acceleration due to gravity.
Based on Eq. ( 49), both the geometric and material nonlinear analyses, which involve reducing the flexural bending and considering the effects of creep, can be performed in a single mathematical operation, which also accounts for the effects of the normal force, without requiring complex computational tools and/or iterative calculations.Further details on the development of this specific analytical procedure have been provided in a previous report (Wahrhaftig et al., 2013).

Creep calculation
To determine the effects of creep, the initial zero-load condition at time t0 was assumed to extend over 28 days from the time of production under standard conditions, while the humidity rose to 70%.The parameters listed in Table 2 were used for the analysis.The area of the concrete was calculated using the average of the external and internal diameters of the superstructure (tower), i.e., D = 733.3mm and d = 446.7 mm, respectively; thus, the perimeter in contact with the environment (πD) was 2300 mm.

Calculation of generalized mass
Hereafter, Roman numeral subscripts are used to denote the tower segments.The generalized mass can be calculated by taking the following integrals: Latin American Journal of Solids and Structures, 2019, 16(5), e196 20/31 Summing the above integrals yields the generalized distributed mass: The total generalized mass of the structure is given by Eq. ( 52), which is a modification of Eq. ( 44) that accounts for mass variation.For m0 in Table 2, the total generalized mass of the structure considered in the present study is 7848.06kg. (52)

Calculation of generalized geometric stiffness
To compute the generalized geometric stiffness, it was necessary to determine the normal forces within the parts defined by the geometry.From top to bottom along the length of the structure, the normal forces are as follows: where g is the gravitational acceleration.Summing these six terms yields the total normal force N: The geometric stiffness can then be calculated as follows: [ ] Latin American Journal of Solids and Structures, 2019, 16(5), e196 21/31 [ ] The generalized geometric stiffness Kg(m0) of the structure was obtained by summing these five terms: Based on Eq. ( 56) and the value of m0 in Table 2, Kg(m0) was computed to be 2.631 kN/m for the present structure.

Calculation of generalized conventional stiffness
The conventional elastic stiffness components can be expressed as The summation of the above components yields the generalized elastic stiffness: Based on Eq. ( 58), the initial generalized elastic stiffness of the present structure was found to be 9.471 kN/m.Note that no creep occurs in the foundation; thus, the modulus of elasticity (Ef) is not a temporal function.However, E(t) represents the material behavior of the tower, and it can be described using Eq. ( 13), implying that Ec(t0) = Ec(t28) = Ec.

Calculation of the generalized soil-spring stiffness
For the specific case under consideration, the elastic soil parameter can be set to Sop = 2668.93kN/m 3 (Sop = Sop1 = Sop2).The distributed spring stiffnesses in the first and second segments can then be expressed as ( ) ( ) = Finally, the generalized total stiffness of the structure was obtained as the algebraic sum of the different contributors, i.e., a function of two variables in accordance with Eq. ( 47).At t0 and for the value of m0 in Table 2, the generalized total stiffness was found to be 7.363 kN/m.

Frequency calculation
The frequency of the first vibration mode was computed numerically using Eq. ( 48), which is a function of the mass at the top and time and considers creep.Based on the parameter values in Table 2, the frequency of the first mode without considering creep was found to be 0.160 Hz.The corresponding value computed by performing linear analysis without considering the geometric stiffness was 0.185 Hz.By contrast, the frequency varies with time when creep is included, as shown in Figure 6.At 4000 days after the commencement of loading, the frequency decreased toward an asymptotic value of 0.147 Hz, illustrating the need for structural stability verification.Discarding the geometric stiffness, that frequency is 0.173 Hz.
To apply Eq. ( 49) directly in an analysis using E´, I´, and ′ m as the unique modulus of elasticity, moment of inertia, and unit mass of the structural system, respectively, it is necessary to employ a criterion to weight the geometric properties of the structure.One possible weighting scheme is shown in Eq. ( 60), where A1 and I1 denote the area and moment of inertia, respectively, calculated using the average diameter of the body structure, and A2 and I2 denote the area and moment of inertia, respectively, calculated using the average diameter of the foundation.ρ, ρ1, E and E1 are the density and modulus of elasticity of the pole and foundation.
, , , , Considering the mass of the structure and foundation shaft, the weighted parameters for the present case are as follows: area A´ = 0.289 m 2 , moment of inertia I´ = 0.0138 m 4 , homogenizing factor Fhsave = 1.062, modulus of elasticity E´ = 18615.81MPa, and density ρ´ = 2586.957kg/m 3 .Therefore, the geometric and material nonlinear frequency was calculated using the weighted parameters to be approximately 0.120 Hz at zero time and 0.101 Hz after 400 days in operation.It is also possible to consider creep in this calculation by using the temporal modulus of elasticity, as estimated using Eq. ( 13), and simply considering E´ as Ec(t0) = Ec(t28) to obtain E´(t) for Eq. ( 49).

Calculation of critical buckling load
In addition to the changing frequency of the structure, the variation of the mass m0 at the top of the tower causes the force acting at the top to vary (in accordance with Eq. ( 61)).Thus, the critical buckling load (Nbuck) at zero frequency can be defined as The variation of the critical buckling load with the mass m0 at the top of the tower is shown in Figure 7; thus, the critical buckling load of the structure can be computed at any stage during its service life.

FINITE ELEMENT ANALYSIS
It should be noted that, while the analytical solution presented in the previous Section provides a single functional form for the entire problem domain, the FEM formulation establishes interpolation functions that are restricted to the domain of each finite element.Regarding the mathematical procedure, FEM analysis represents a geometric nonlinear formulation based on a consideration of the geometric stiffness.According to Bathe and Wilson (1973), in terms of modal analysis, the relevant eigenvalues and eigenvectors can be obtained by solving the following secular equation: where (M) is the mass matrix and (K) is the stiffness matrix, which includes the geometric stiffness term for nonlinear cases, formulated similarly to Eq. ( 47).In the FEM environment, [ω 2 ] represents the eigenvalues and [Φ] represents the eigenvectors.The spring matrix is a 6 × 6 symmetric matrix of the spring coefficients kij, including all the translational and rotational degrees of freedom of a bar element.The components of this matrix are nodal springs.The mathematical development of the FEM system, as described above, is based on modal analysis performed by reducing the stiffness of the system using geometric stiffness matrix components.The importance of the geometric stiffness to the FEM solution for cases similar to that of the present study was highlighted by Armand and Lin(1993), who employed a full stiffness matrix to derive a dynamic solution for the frequency of a 30-m-tall tower.
The buckling analysis requires the solution of another generalized eigenvalue problem.The equation governing the buckling of an assemblage of structural elements can be expressed as where [λ] is the diagonal matrix of the eigenvalues and [ψ] is the matrix of the corresponding eigenvectors (mode shapes).The stiffness matrices (K0), (Km), and (Kg) are as previously stated.
To assess the accuracy of the proposed equations for the frequency and critical buckling load, their results were compared with those given by an FEM computational model.The considered structure was modeled using frame elements with constant or variable cross-sections, as appropriate.The forces presented in Table 1 were applied to the model with the corresponding masses.The spring factor Sop was assigned to the foundation frame element as a linearly distributed parameter.Only a lateral spring was used to model the foundation system, and the spring stiffness matrix was thus considered to be null at all positions corresponding to the axial and rotational degrees of freedom, while the values extrapolated from the interior of each element joint were considered to be active components of the matrix.Figure 8 shows the three-dimensional model constructed using SAP2000 (2015) structural analysis software, as well as the discretization of the structure in Figure 1 using 45 frame elements.The rightmost panel in the figure shows the cross-sectional details of the structure, where S indicates the type of section.In nonlinear problems, discretization significantly influences the results under dynamic conditions, with deviations of up to 21% observed between the results obtained by linear and nonlinear analyses (Wahrhaftig, 2013).In the present study, the natural frequency of the tower was calculated by performing nonlinear static FEM analysis, with the effects of the normal forces taken into consideration using the equations provided above.This method is regarded as a geometric nonlinear approach (SAP2000, 2015).Note that the global stiffness matrices of the structure were obtained automatically by the analysis software.The first vibration mode frequency and buckling load were computed by using the FEM and assuming an isotropic homogeneous material using the parameters employed in the analytical investigation and a Poisson's ratio of 0.2.Note that the interpolation functions commonly used in FEM equations are towers Alexandre de M. Eq. ( 49) does not have a practical form, and the frequency is overestimated.The soil stiffness accounted for roughly 12.36% of the total stiffness of the structure.The natural frequency of the structure 4000 days after installation (0.147 Hz) was 8.56% lower than that at the time of installation.This result implies that, when the mean lateral loading due to wind is calculated principally as a function of the first natural frequency and the corresponding vibration mode, any frequency reduction would be critical to the dynamic behavior of the structure.The critical buckling load was estimated using the analytical process as a compressive force of 307.687 kN at time zero, corresponding to the first natural frequency.The FEM-predicted critical buckling load under identical operating conditions was 263.602 kN.For lightly loaded telecommunication towers, the risk of buckling may be insignificant.Nevertheless, verification is necessary owing to the extreme slenderness of the structures, associated with the growing demand for expansion of the mobile cellular telephone system.An alternative means of meeting this demand is the sharing of existing structures to apply new loads.For other slender structures subjected to heavy loads, such as elevated water storage reservoirs, this verification process is mandatory.
The comparison of the FEM results with those of Rayleigh's method indicated that the frequency given by the former was initially 3.01% lower than that yielded by the latter, whereas for buckling, the FEM results were 14.33% lower too.The difference was not as large for the frequency, and it is considered to be acceptable for engineering practice; however, for buckling, the difference was found to be significant.Considering the values at the final time, the difference was computed to be 5.94% for frequency and 21.45% for buckling.Although the frequency difference is larger in this case, it is still considered to be reasonable, as there is usually a safety margin of around 40% (EN 1992(EN -1-1, 2004)), which includes model errors.Silva et al. (2013) studied 90 RC towers similar to the one investigated in this study and developed a procedure to compute the dynamic wind effects on such a structure.In one approach, they proposed that the dynamic magnification due to wind loading is linearly proportional to the first natural vibration frequency; thus, a decrement in the frequency is proportional to an increment in the dynamic load.For buckling, this difference is extremely large and may become infeasible for engineering purposes, as the FEM values are better.Similarly, considering the mathematical approach, the FEM involves many more vectors in the function vector space compared to Rayleigh's method, which requires only one function for the entire domain.These functions are the interpolation functions.Consequently, the FEM tends to be more precise, but this statement cannot be generalized.Simultaneously, a correction factor could be defined, as a function of the structure age, to fit both methods.These could be topics of future research, including tests and numerical analysis.A structural analysis further refined the grid of the FEM, considering 51 elements (seven elements more), and the results were similar to those presented here (with differences of 0.2% for the buckling load and 0.1% for the frequency), implying that the adopted discretization is adequate.Wahrhaftig and Brasil (2017) performed an experimental investigation of the wind loads on slender towers by using low frequencies such as those employed in this study.As can be seen in Figure 10, they captured the effects of the frequencies in the acceleration measured on site.It was observed that the first frequency was the most important one in the behavior of the structure, being responsible for more than 90% of the acceleration of the structure.The second frequency had a small influence, and higher frequencies did not influence the dynamic loads at all.Considering the root sum square (RSS) to combine the responses of several modes, the first mode is responsible for more than 99% of the dynamic effects of the wind loading.Therefore, the computation of the first vibration frequency, in addition to the first mode, is crucial in this type of design.It is important to state that this kind of structures, used as telecommunication towers, has the function of transmitting electromagnetic waves from the tower to the cell phones and from one tower to other.In this case, the displacements and rotations of a given section may cause loss of signal and consequently loss of the communication between towers.So, during the design process it is necessary to introduce constraints limiting the rotation of the cross-sections to avoid loss of signal.

CONCLUSIONS
In this study, a new way of application of an analytical procedure for computing the first natural vibration frequencies of RC slender towers was investigated.The method was applied to a 46-m-high RC telecommunication tower, and the results were compared with those obtained using the FEM.The buckling critical load was computed as well, by considering the geometrical and physical nonlinearities in addition to the soil-structure interaction.To account for these factors, the effects of cracking, as a cause of stiffness reduction, and material creep were included.The properties of the cross-sections changed as functions of the height of the structure above the ground.All of the analysis performed in this study was conducted using only one shape function, which included clear boundary conditions.The solution developed in this study is based on Rayleigh's method conducted by adopting a single shape function and considering the boundary conditions.The main contributing of this method is that the first natural vibration frequency can be calculated accurately for relatively complex problems without sophisticated tools that would require considerable modeling time and expensive software.An important application of the proposed mathematical procedure is the calculation of wind-induced forces.For the considered high-rise towers, the main loading comprises the wind-induced forces.This type of slender cantilever structure, when excited by wind, usually vibrates mainly in the first mode, and the equations developed in this study are thus useful for practical engineering applications.The proposed technique can also be applied to similar structures with slightly different geometries and mass distributions.
Excellent agreement was observed between the first vibration frequencies obtained using the proposed method and FEM analysis.Specifically, the difference between the results of the proposed method and FEM analysis was only 3.0% initially and 6.3% finally, considering creep.As this type of structure is designed with a safety margin of 40%, which also includes the model errors, these differences are very reasonable for practical engineering purposes.It was noted in the present study that the effect of creep reduces the first frequency by 9.3%, which indicates that it is important to consider this effect on the material.
However, the buckling load computed using this analytical method was 13.9% larger than that given by the FEM.A correction factor function could be used to fit the results and minimize the error between the different formulations.The difference between the results could be related to the fact that many more interpolation functions are used in FEM analysis than in Rayleigh's method, in which only one function is used for the entire domain of the problem.However, confirming this statement is not a straightforward task.Future analysis and specific studies will be required to obtain an appropriate conclusion.Therefore, further investigation, including the use of experimental methods, is necessary to complement this study.

Figure 1 .
Figure 1.Case study of mobile communication tower.

Figure 2 .
Figure 2. Segments and cross-sections along the tower.
(a).Each reinforcement bar bi assumes a given position i in the cross-section domain, defined by Rbi and θi, as shown in Figure3(b).

Figure 3 .
Figure 3. Cross-section adopted for the RC structure

Table 2 .
Parameters used for the analysis and computations.Parameter Value Compressive strength of tower concrete fcm = 45 MPa Elastic modulus of tower concrete at 28 days Ec = 38097.35MPa Reduced elastic modulus of tower concrete 0.5 Ec = 19048.67MPa Compressive strength of foundation concrete fcm = 20 MPa Elastic modulus of foundation concrete at 28 days Ef = 31460.05MPa Reduced elastic modulus of foundation concrete 0.5 Ef = 15730.02MPa Density of structure concrete ρ = 2600 kg/m 3 Density of foundation concrete ρ1 = 2500 kg/m 3 Lumped mass at tower top m0 = 1097.76kg Additional distributed mass per unit height e m = 40 kg/m Elastic parameter of soil Sop = 2668.93kN/m 3 Gravitational acceleration g = 9.806650 m/s 2

π
The variation of the segment diameter between the base and the shaft is represented by the linear relationship ( ) Similarly, D5, th5, d5, A5, and I5,    which denote the external diameter, thickness, internal diameter, area, and moment of inertia of the last segment, respectively, are D5 = 700 mm, th5 = 130 mm, segments are connected by the segment s4 of variable cross-section.This begins with a diameter and thickness equal Latin American Journal of Solids and Structures, 2019, 16(5), e19612/31 to those of the third section and ends with those of the fifth section (i.e., th4 = 150 mm).Its area and moment of inertia can be expressed as
m0 is the lumped mass at the tip of the element joint and ( ) s m x is the mass per unit length.The total generalized mass is thus given by the distributed mass in each segment s.Therefore, ( ) s m

Figure 5 .
Figure 5. Variations of the structural deformation and modulus of elasticity with time.
respectively.Using the following equation, the generalized spring stiffness KSo was computed to be 1.123 kN/m:

Figure 6 .
Figure 6.Variation of the fundamental natural frequency of the structure with time.

Figure 7 .
Figure 7. Variation of the critical buckling load with time.

Figure 8 .
Figure 8. FEM model of the tower structure.

Table 1 .
Physical details of the tower structure and installed equipment.