Using Dynamic Analysis to Adjust the Rheological Model of Three Parameters to the Eurocode Creep Criteria

A dynamic analysis of vibration for considering a three-parameter rheological model to fit the same results as predicted for creep by the Eurocode (EN 1992) criteria is performed based on the adjustment of its parameters. The use of a rheological model of three parameters as a valid alternative for real problems brings a huge facility for mathematical implementation and manipulation due the simplic-ity of the solution. For adjustment of the elastics and the viscous parameters, a numerical simulation to calculate the fundamental frequency of an actual reinforcement concrete pole is carried out in comparison with the standard Eurocode criteria. In this determination, the geometry variation, a concentrated force present at the free end of the structural element, and the self-weight of the structure are considered. The physical nonlinearity of the concrete due to the cracks is also considered by reducing of the flexural stiffness, and its viscoelastic behavior is included in the calculation through a temporal modulus of elasticity. In the analysis, the ground was modeled as a set of distributed springs along the foundation length. The frequency over time is then analytically calculated as the critical buckling load for different instants after the structure to be loaded.


Introduction
A column represents a continuous structural member whose vibrations are governed by nonlinear partial differential equations for which exact analytical solutions cannot be found, as pointed out by [1]. Columns constitute continuous systems, and their analysis can be reduced to an analogous system containing a single degree of freedom. The vibration mode is restricted to a configuration previously established by a mathematical function that describes the vibratory movement, and the properties of the system can be expressed as generalized coordinate functions [2]. In his study on the vibration of elastic systems applied this technique considering the function valid throughout the problem domain. However, for real cases, where the properties of the structural elements vary along their length, the formulation developed for calculating the stiffness and mass must be solved by observing the intervals defined in the geometry. In these cases, the integrals obtained can be solved within the limits established for each interval, i.e., the generalized properties can be calculated for each discrete segment of the structure, as defined by its geometry within that segment. A variety of vibration problems using that mathematical concept were solved by [3] who mentioned a previous one [4] where the buckling load is calculated for stepped and tapered columns and where how laborious or even impossible it is to apply it for problems with variable geometry is registered. With the advent of digital computers, these problems passed to be solved by modeling that use discretization technics of the continuum [5].
To analytically define the fundamental frequency for the case modeled in this study, all the elastic stiffness components are considered in the calculation, including the conventional stiffness, which depends on the material behavior; the geometric stiffness, which depends on the normal force acting on the structure; and the soil parcel, which accounts for the soil-structure interaction. It is important to note that the soil-structure interaction cannot be ignored, particularly in the case of a monopile foundation, because it may significantly influence the dynamic behavior of the structure [6].
The structure selected for this study is a slender reinforced concrete (RC) having both full and hollow circular section with variable geometry, for which the natural frequency and the critical buckling load were calculated considering all nonlinearities present in the system. It is important to highlight that nonlinearities play an important role when calculating dynamic proprieties of a system, as well pointed by [7]. In this work, the geometric nonlinearity was taken in consideration by using the geometric stiffness parcel into the total stiffness of the system. The nonlinearity of the material was computed by reducing its flexural stiffness, as similarly done by [8], reflecting the development of cracking in the concrete when bended, which is dependent on the magnitude of the stress. Another kind of material nonlinearity is creep, which occurs due the viscoelastic behavior of the concrete, it being considered in two ways. The first one is the mathematical model for creep predicted by Eurocode 2 (European Standard EN 1992-1-1) [9]. The second one is a threeparameter viscoelastic model whose parameters are adjusted in order to meet the results obtained when using the Eurocode. In this sense, the use of the threeparameter viscoelastic model to represent the creep of concrete brings an enormous facility of employment for actual cases due the reduced number of variables which are manipulated. Indeed, just one of them is necessary because two of the three parameters can be expressed in terms of the modulus of elasticity of the concrete, a data easily calculated for any standard procedure or obtained in laboratory. Figure 1 presents the bar model of a structure in free vibration. Consider the following trigonometric function, taken as valid throughout its domain:

Analytical solution of the structural frequency
where x is the location of the calculation, originating at the base of the cantilever, and L is the length of the column.
That model represents a column under an axial compressive load, N s (x), with either constant or variable properties along its length. These properties include the geometry, elasticity/viscoelasticity, and density. Applied springs of variable stiffness k so (x) act as the lateral soil resistance until the foundation elevation Gr.
The system is under the action of gravitational normal forces, originating from the distributed mass along the length of the column and of a lumped mass at the tip m 0 .
In the case of vibration of a cantilevered column that is clamped at its base and free at its tip, the shape function given in Eq. (1) satisfies the boundary conditions of the problem. The use of Eq. (1) as a shape function for an actual structure with varying geometry has been validated by [10]. This validation involved a comparison with a computational solution derived using computational modeling by finite element method (FEM) and other mathematical expressions.
By applying the principle of virtual work and its derivations, the dynamic properties of the subject system are obtained. The elastic/viscoelastic conventional stiffness is given by where for a segment s of the structure, E s (t) is the viscoelastic modulus of the material with respect to time; I s (x) is the variable moment of inertia of the section along the segment in relation to the considered movement, obtained by interpolation of the previous and following sections and if it is constant, it is simply I s ; k 0s (t) is the temporal term for the stiffness; K 0 (t) is the final conventional stiffness varying over time; and n is the total number of segment intervals given by the structural geometry. In Eq. (2), obviously, t vanishes when the analysis considers a material with purely elastic, time-independent behavior. The geometric stiffness appears as a function of the axial load, including the self-weight contribution and is expressed as where k gs (m 0 ) is the geometric stiffness in segment s, K g (m 0 ) is the total geometric stiffness of the structure with n as defined previously, and N 0 (m 0 ) is the concentrated force at the top, all of which are dependent on the mass m 0 at the tip, given by Further, N j is the normal force from the upper segments, obtained by Then, the total generalized mass is given by considering that where m s x ð Þ is the mass distributed to each segment s, which is obtained by multiplying the cross-sectional area, A s (x), by the density, ρ s , of the material in the respective interval. Therefore, m s x ð Þ is the mass per unit length, and m is the generalized mass of the system owing to the density of the material, with n as previously defined. If the cross section has a constant area over the interval, A s (x) will be just A s ; consequently, the distributed mass will also be constant. Similarly, if the mass m 0 does not vary, all the other parameters that depend on it will also be constant.
One approach for considering the participation of the soil in the vibration of the system is to consider it as a series of vertically distributed springs that act as a restorative force on the system. With k Sos (x) denoting the spring parameter, the effective soil stiffness (as a function of the location x along the length) is generally defined as where the parameter K So is an elastic characteristic consisting of the sum of k Sos (x) along the foundation depth, which depends on the geometry of the foundation D s (x) and the soil parameter S os . Considering the normal force as positive, the total structural stiffness is obtained as Finally, the natural frequency (in Hertz), as a function of the time and the mass at the tip, is calculated according to Eq. (11). The great advantage of using that equation in terms of two independent variables is that it can be employed to calculate the critical load of buckling as well, because all the generalized parameters are expressed as a function of the mass at the top. Details of this analytical procedure can be seen in [11]:

Creep consideration
The creep represents the increase of deformation under constant stress, which occurs in some materials due to its viscoelastic nature. It is essential to consider it in the analysis of slender structural elements, because the stiffness of these members is modified as a function of the rheology of the material. Usually, viscoelasticity is associated with creep of structural elements and can be characterized by models where the immediate elastic deformation is increased by viscous deformation, resulting in a temporal function for deformation. Consequently, the modulus of elasticity must also be provided as a temporal function that provides accurate results under normal levels of stress. Due to the viscous nature of the concrete, even at a constant stress level, the deformation of a structural element tends to increase over time. An increase in strain over time under constant stress is a viscoelastic phenomenon.

Solution of the three-parameter rheological model
Mathematically, viscoelasticity can be represented by a time-dependent function associated with rheological models capable of describing the phenomenon. It is conceptually convenient to consider classic viscoelastic models in which there are only two types of parameters, relating to elasticity and viscosity. Classic viscoelastic models are obtained by arranging springs and dampers, or dashpots, in different configurations. Springs are characterized by elastic moduli and dashpots by viscosity coefficients. The best known of these mechanical models are the Maxwell model, containing a spring in series with a dashpot, and the Kelvin-Voigt model, containing a spring and dashpot in parallel [12]. One model used to represent the viscoelasticity of solids is the three-parameter model, in which the elastic parameter E e is connected to the viscoelastic Kelvin-Voigt model with parameters E v and η, which is a simplification of Group I of the Burgers model, as shown in Figure 2.
The three-parameter model is an appropriate model for describing the viscoelastic nature of many solids [13] and is often used to study the phenomenon in various scientific fields. Adaptation of the Burgers model in different fields of structural analysis can be found in [14][15][16][17][18][19]. The total deformations of the Kelvin-Voigt model are given by ε = ε e + ε v , where ε e is the deformation of the elastic model and ε v is the deformation of the Kelvin-Voigt model. When differentiated with respect to time, the total deformation is obtained as which includes the constitutive equations of the elastic and Kelvin-Voigt models, respectively. Considering the modulus of elasticity for both parts, elastic and viscous, the stress becomes From the previous equations, one derives the following differential equation: where σ = 0 for t < 0 and σ = σ 0 for t > 0, with t representing the time and t = 0 the instant of loading application. As the stress remains constant, the derivative of the stress with respect to time is zero. Applying the previous stress condition, the following ordinary differential equation is found: for which the general solution for t > 0, taking the initial condition ε( From Eq. (16), it is possible to extract the temporal function for the modulus of elasticity of the three-parameter model: Therefore, it is easily seen that for It is important to note that the viscoelastic behavior of the considered material is completely represented by the temporal modulus of elasticity and it can be used for static or dynamic applications. For instance, the previous solution was used in numerical simulations as can be seen in [20,21].
It is possible to transform the parameters of the viscous part to being just a function of the modulus of elasticity of the elastic part, which can easily be calculated by any standard procedure or obtained in the simplest laboratory. Therefore, these parameters can be written as where α is a real positive number and γ brings together a temporal unit.

Model predicted by Eurocode
The method specified in European Standard EN 1992-1-1 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 E c . 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 and creep. All the factors related to the phenomenon, such as loading and environment humidity, 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 f cm (f cm = f ck + 8, f ck in MPa). The creep coefficient φ(t, t 0 ), as defined in Eq. (22), is the product of two factors, namely, φ 0 and β c (t, t 0 ), 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 factor φ 0 defined in Eq. (23) consists of three other factors. The first of them, φ RH (given by Eq. (24)), considers concrete compressive strengths >35 MPa (as in the case that will be seen) and accounts for the effects of the environmental relative humidity RH, the equivalent thickness h 0 of the member which is a function of the cross-sectional area A c , and the external perimeter u e of the member in contact with the environment. The second one, β(f cm ) (Eq. (26)), represents the direct effect of the resistance on φ 0 . The third, β(t 0 ) (Eq. (27)), takes into account the age of the concrete at the beginning of loading, i.e., at t 0 .
The second factor, β c (t, t 0 ) (Eq. (28)), is a function of the coefficient β H (given by Eq. (29) for average concrete compressive strengths upper than 35 MPa), 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. (22), 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 E c (t 0 ) is the modulus of elasticity at the beginning of loading and E c (t 28 ) is the modulus of elasticity 28 days after the commencement of loading.

An application
The case selected for the present study involves calculating the fundamental frequency and the critical buckling load of an actual slender reinforced concrete pole with variable geometry that presents both geometrical and material nonlinearities as shown in Figure 3.
The structure is 46 m high, which includes a 40 m superstructure with a hollow circular section and a 6-m-deep, full circular-type foundation. The moduli of elasticity adopted for the superstructure and foundation are 30.24 and 24.97 GPa, calculated by Eq. (33) considering characteristic resistances (f ck at 28 days after production) of 45 and 20 MPa, respectively: A set of antennas and a platform are installed at the tip of the structure, constituting a concentrated mass of 1098 kg. Cables and a ladder are installed along the entire length, adding a distributed mass to the system of 40 kg/m. The densities of the reinforced concrete were defined as 2600 and 2500 kg/m 3 for the super-and infrastructure, respectively. The physical nonlinearity of the material was computed for the superstructure and the foundation reducing the gross moment of  inertia by a multiplier factor equal to 0.3, allowing the performing of a simplified nonlinear analysis according to Eurocode 2, as presented in [22], but being possible the use of other coefficients as explained by [23].
The foundation is a relatively deep shaft having a bell diameter of 140 cm, bell length of 20 cm, shaft diameter of 80 cm, and shaft length of 580 cm. The lateral soil resistance is represented by an elastic parameter, S os , equal to 2667 kN/m 3 .
The geometric details of the evaluated pole are shown in Figure 4, where g denotes gravitational acceleration; Gr means ground; s represents each structural segment; S, D, and th are the type, the external diameter, and the wall thickness of the section; d b represents the reinforcing bar diameter; n b is the number of reinforcing bars; and c´is the reinforcing cover. The slenderness ratio of the tower structure is approximately 400.
Because this is an RC structure, it is necessary to account for the presence of the reinforcing bars when calculating the moment of inertia, which is accomplished by homogenizing the cross section. Therefore, according to the theorem of parallel axis, the factors, which multiply the nominal moment of inertia of the section in terms of the total moment of inertia of the reinforcing steel, in the homogenized section are appropriately calculated. Studies that assure the occurrence of the transfer of creep to the reinforcement of columns were development by [24,25].

Simulation results
Considering E e is equal to E(t 0 ) = E(t 28 ) and setting α and γ as 3.913 and 10 8.3066 seconds, respectively, and adopting an environmental humidity of 70%, the variation of the fundamental frequency for different instants in the lifetime of the structure can be obtained. The produced results by using the Eurocode model can be observed in Figure 5a, and the result by using the adjusted three-parameter viscoelastic model can be seen in Figure 5b. It is important to mention that these chosen values for α and γ were defined so that simulation leads a good agreement for instants approaching and after 2000 days. Therefore, they were intentionally defined so that the frequency met the same values as given by Eurocode. The choice of these coefficients has been done because the convergence of the deformations occurs at 4000 days, at which time the interest of the structural engineering normally lies, being, however, possible to define other pairs of values for α and γ in the case of a particular objective or even to choose which can match both curves in the whole time interval. Therefore, the mentioned coefficients have been adjusted so that the frequency is equalized by both models considering a precision of six significant digits, as can be highlighted in Table 1. When the modulus of elasticity is calculated by both models, that precision is not reached. Figure 6 shows a comparison between results produced through both models, considering each selected instants of time.
By using the presented dynamic procedure, the critical buckling load is determined when the frequency is zero at any arbitrary time after the structure gets into service. Taking all the previous explanation into consideration and varying the mass at the tip, the force acting at the top also varies according to Eq. (5), as does the frequency of the structure that varies according to Eq. (11). The results obtained for the buckling load for both models can be seen in Figure 7. To obtain it, a short routine of programming has been elaborated considering increments of 0.1 kg to the lumped mass.  Table 1.
Frequencies for both models at selected instants.

Conclusions
• Because of the viscoelastic behavior of the material, the modulus of elasticity presents a variation along the time, reflecting on the structural frequency and critical buckling load.
• A three-parameter viscoelastic model has been adjusted to fit the same results as predicted by Eurocode creep criteria for a specific interval of time. The use of two parameters makes the adjustment process more flexible.
• It is important to stay clear that that adjustment does not lie on the adjustment of the modulus of the elasticity which does not have the same precision when observed for both models; for that a dynamic analysis is important.
• This article demonstrated the possibility of adjustment of a simple model to the standard one and how easy it is used for practical applications to calculate the first natural frequency as the critical buckling load.
• For future works, a programing routine for obtaining a finer adjustment of the curve between the viscoelastic rheological model of three parameters and that of the model for creep as predicted by Eurocode must be developed.
• Further, comparative analyses considering other values of environmental humidity should be also performed.