Theoretical Framework forCreepEffectAnalysis ofAxially Loaded Short CFST Columns under High Stress Levels

Due to its excellent mechanical performances, axially loaded concrete-filled steel circular tube (CFST) columns have been widely used in structural engineering. As an important long-term behaviour of CFST structures, the creep has an obvious nonlinear property under high stress levels, which makes the influence of creep more complicated. In this study, to analyze the impacts of nonlinear creep effect on the behaviour of axially loaded short CFST columns, a complete theoretical framework for coupling analysis of 3D creep effect and material nonlinearity was presented. First, the concrete damaged plasticity model with a uniform constraint (UCCDP) was established to simulate the plasticity and damage evolution of a concrete core. Next, based on the UCCDP, a method of 3D nonlinear creep analysis and a corresponding numerical analysis method were established and implemented in the ABAQUS secondary platform. Finally, by comparing the predicted results with the experimental results, it was observed that the method proposed to predict the creep of axially loaded short CFST columns had satisfactory accuracy.


Introduction
Owing to their merits of high compressive strength, durability, construction, and low cost, concrete-filled steel tubes (CFST) have been widely used in structures [1,2]. Nevertheless, the results obtained by intensive research have demonstrated that creep causes structural deformation and internal force redistribution, which has an important effect on the long-term mechanical behaviour of CFST columns [3]. eoretical and experimental studies have demonstrated that there exists a nonlinear relationship between creep and stress, when the stress of the concrete core exceeds a certain level (0.6f c ), which obviously indicates nonlinear characteristics [4,5]. Fu [6] analyzed the stress levels of a concrete core in super-high-rise buildings, such as the Harbour Building in Xiamen city and the Dacheng Financial & Business Center in Kunming city, and reported that the stress levels have reached 0.69f c . e nonlinear creep effect does not only cause excessive deformation, but also has further impact on dynamic structural characteristics [7]. To analyze the long-term behaviour of CFST columns under high stress levels, it is necessary to establish a nonlinear creep analysis theory.
For the nonlinear creep theory, three different theories have been proposed: one is to establish a rheological model to simulate the creep characteristic directly through the uniaxial strain creep test under high stress level [8]; the second is to establish the nonlinear functional relationship between the creep responses and the stresses based on the uniaxial creep characteristic under low stress levels; the third is to establish the nonlinear creep constitutive equations of concrete under high stress levels by using damage mechanics [9][10][11][12]. For the first two theories which are both based on the creep test under uniaxial stress, it is difficult to apply to the creep analysis under a complex stress state. For the third kind of method, the reason for the nonlinear creep is attributed to the material damage evolution, and the application is also more convenient; however, the results of existing studies have various shortcomings. Li and Qian [9] established the elastic creep damage analysis model by using continuum damage mechanics, but did not consider the influence of the concrete's plastic behaviour. Huang and Xiang [12] proposed a method of analyzing the nonlinear creep effect on concrete by using the classical uniaxial stressstrain relationship (Hognestad model). However, the uniaxial stress-strain curve does not describe the nonlinear behaviour of concrete completely.
Additionally, for CFST columns with axial compression, the concrete core is always under a uniform constraint, i.e., the creep has obvious three-dimensional (3D) characteristics [13,14].
In this study, a new theoretical framework for modeling the 3D nonlinear creep effect of axially loaded CFSTcolumns with high stress levels is proposed on the basis of continuous damage mechanics and concrete elastic-plastic mechanics. e first part of this theoretical framework is the concrete damaged plasticity model under a uniform constraint (UCCDP) to simulate the plasticity and damage evolution of concrete. e second part is the 3D nonlinear creep analysis model and the corresponding numerical analysis method. Finally, using the creep prediction B4 model [15] and the finite element (FE) package ABAQUS, the numerical analysis method was implemented into the ABAQUS secondary platform, and the analytical results were compared with the experimental results in different stress levels.

Constitutive Model for Concrete under
Uniform Constraint e concrete damaged plasticity (CDP) model has been widely used in engineering projects, owing to its ease of use [16,17]. e CDP model mainly includes the Lubliner yield criterion, Drucker-Prager plastic flow criterion, and damage evolution and hardening criterion. However, this model is mainly suitable to the analysis of the elastic-plastic behaviour of concrete structures under low confining stress [17]. e concrete core of circular CFST columns under axial compression is in the triaxial compression state, and the confining pressure has already affected the mechanical behaviour of concrete.
is is a convenient and effective method to modify the CDP model and analyze the elasticplastic behaviour of concrete under uniform confinement. erefore, in this paper, the defects of the Lubliner yield criterion under uniform constraints are discussed, and the modification of the yield and plastic flow criteria is explained.

Defect of the Lubliner Yield Criterion in CFST.
e Lubliner yield criterion is a fundamental part of the concrete constitutive model and is always used to determine whether concrete has entered the plastic stage. e Lubliner yield criterion is expressed as follows: where σ ij is the effective stress; ε pl is the equivalent plastic strain; are the Macaulay brackets and 〈x〉 is taken as (|x| + x)/2; I 1 is the first stress invariant; J 2 represents the second deviatoric stress invariant; σ max is the maximum effective stress eigenvalue; c is the cohesion; α, β, and c are the undetermined parameters. α and β are determined by the uniaxial tensile strength, compressive strength, and biaxial isobaric strength, and c is determined by the tensile meridian (TM) and compressive meridian (CM) [16,17]. However, for concrete in axially compressed CFST columns, the abovementioned method of determining c has unavoidable drawbacks. e expression of TM and CM in Haigh-Westergaard stress space can be expressed as follows: where ξ is the hydrostatic pressure taken as − I 1 /3 and ρ is the von Mises equivalent stress taken as ρ � (3J 2 ) 1/2 . At a given ξ in triaxial compression, the ρ ratio in equations (2) and (3) is defined as K c , K c � (c + 3)/(2c + 3).
e value of K c determines the shape of the trace line in the deviatoric plane, as shown in Figure 1. Because the value of K c is determined, the value of parameter c is also determined, which means that the opening shape of the yield surface is determined by c. Figure 2 shows the von Mises equivalent stress (ρ) versus the negative hydrostatic pressure (− ξ) curves in TM and CM. e angle between the curves and the − ξ axis is defined as φ. us, the slope of the curves can be expressed as follows: e value of the slopes indicates the opening angle of the yield surface. erefore, parameter c influences both the opening shape and the opening angle of the yield surface. However, after selecting an appropriate parameter c to determine the opening shape, the yield surface may not be able to obtain an appropriate opening angle.

Calibration of the Parameter c.
For concrete in axially compressed circular CFSTcolumns, the opening shape of the yield surface should be controlled by the slope of the CM [18]. Because of the confining pressure exerted by the steel tube, the compressive strength of the concrete core is improved as expressed by the following equation: where f cc is the compressive strength of concrete under a uniform constraint, σ l is the confining pressure, and k is the strength enhancement factor.

Advances in Civil Engineering
When the concrete core fails under uniform confinement, I 1 � f cc − 2σ l and J 2 � 1/3(f cc − σ l ) 2 . By substituting these into equation (1), we obtain the following equation: According to the Mohr-Coulomb criterion [17], where φ is the internal friction angle of concrete.
For concrete under a uniform constraint, σ max � − σ l and σ min � − f cc . By substituting σ max , σ min , and equation (8) into equation (7), c can be rewritten as follows: e angle of internal friction can be calculated as follows [19]: e result of comparing the original calibration method with the current calibration method is presented in Figure 3.

Flow Rule.
As part of the constitutive model, the flow criterion is used to determine the plastic deformation of concrete when the concrete enters the plastic state. e widely used Drucker-Prager flow criterion, which is expressed by equation (11), was adopted in this study: where ζ is the eccentricity ratio, always taken as 0.1, σ t0 is the uniaxial tensile strength of concrete, and ψ is the concrete's dilation angle.
To accurately analyze the plastic deformation, we must obtain a reasonable dilation angle (ψ). e physical meaning of ψ represents the plastic deformation rate of the concrete.
us, a larger dilation angle means a larger plastic deformation rate and leads to the steel tube having a stronger constraint effect on concrete [20]. Wang et al. [21] investigated the relationship between the dilation angle and the hoop stress of ultrahigh-performance concrete, but did not consider the effect of material damage evolution, which leads to an obvious underestimation of the components' plastic deformation. Based on this, the calculation method of the dilation angle considering the influence of damage evolution was deduced again. e strain of concrete includes elastic strain and plastic strain, as expressed by the following equation: According to the theory of elastic-plastic mechanics, the plastic strain increment can be expressed as follows: where λ is a nonnegative scalar. e circumferential plastic strain (ε 1 p ), the radial plastic strain (ε 2 p ), and the longitudinal plastic strain (ε 3 p ) can be expressed as follows: where d is the damage factor as a function of equivalent plastic strain and v is the elastic Poisson's ratio of concrete. By substituting equation (11) into equation (13), we obtain the following equations: Advances in Civil Engineering 3 Because ζ is always taken as 0.1,

���������������
(ςσ t0 tan ψ) 2 + 3J 2 can be approximately replaced by ��� 3J 2 . By adding up equations (16) to (18), we obtain Equation (18) − 0.5 × (equation (16) + equation (17)) gives By substituting equation (20) into equation (19), we obtain e positive value of tanѱ represents the volumetric dilation, and its negative value represents the volumetric contraction. If reliable test data are not available, the dilation angle can be calculated by the following equation [19]: where ξ c is the confinement coefficient and is taken as ; A c and A s are the sectional area of the steel tube and concrete core, respectively; f y is the yield strength of the steel tube; and f c ′ is the cylinder compressive strength of the concrete core.

Damage Evolution and Hardening
Rule. By introducing the damage factor d, the relationship between the strain tensor (σ ij ) and the plastic strain (ε ij ) can be expressed as follows [16]: where D 0 ijkl is the concrete's initial elastic matrix and d is the damage index expressed as follows: where dt and dc are the tensile and compressive damage indices of concrete, respectively, as expressed in equations (27) to (29); s t and s c represent the recovery functions of the damage index, which consider the effects of crack closure and the reduction of tensile strength caused by compressive damage, as expressed in equations (30) to (31): For axially loaded short CFST columns, r(σ) equals 0 because 0 > σ1 � σ2 > σ3. Additionally, s t is always taken as 1 and s c is always taken as 0.
e evolution law of the concrete yield surface (hardening criterion) is still driven by the equivalent plastic strain (ε p ij ) as expressed by equation (33).

Processing of 3D Characteristic Creep Effect.
As mentioned above, the creep has an obvious 3D characteristic for CFST columns. However, research on the creep law under multiaxial stress states is ongoing and a prediction model that can be used directly has not yet been developed. erefore, by introducing the concept of creep Poisson's ratio [22], the creep law under uniaxial stress is extended to multiaxial stress states. is can be considered as a simple and feasible method when the experimental data are insufficient.
Neville et al. [23] summarized the results of existing experimental research and reported that the creep Poisson's ratio μ cp,i can be expressed as follows: where μ cp,i is the creep Poisson's ratio corresponding to the principal stress direction and σ i is the principal stress. For convenience in engineering practice, the effect of Poisson's ratio caused by creep is established based on a linear elastic formula. e transformation matrix [A], which considers the influence of the creep Poisson's ratio, can be expressed as follows: For axially compressed CFST columns, owing to 0 > σ1 � σ2 > σ3, equation (34) can be expressed as follows: where q′ is the hoop stress coefficient.

3D Nonlinear Creep Constitutive Model.
e total strain at time t should include three parts: the first part is the elastic strain at loading age τ and the corresponding creep strain; the second part is the elastic strain and the corresponding creep strain caused by a stress increment in the (t − τ) period; the third part is the shrinkage strain (ε S ) and the temperature strain (ε T ). erefore, the total strain can be expressed as follows [24]: where J(t, τ, σ) is the creep compliance function and is expressed as follows: where E(τ) is the elastic modulus of concrete at loading time τ and C(t, τ, σ) is the creep function obtained as follows: where d(τ) is the damage factor related to the concrete stress states, whose value is in the range [0, 1], C(t, τ) must be expressed in the form of a Dirichlet series to avoid integrating the stress histories of the entire time, and a i (τ) and λ i can be determined as follows: (1) High precision is generally already achieved when m is equal to 4 and λ 1 � 1, λ 2 � 0.1, λ 3 � 0.01, and λ 4 � 0.001 [22].
(2) If the creep effect of time t is calculated, then (t − τ) must be divided into n segments (t 1 , t 2 , . . ., t i , . . ., t n t i > t 0 , t n � t), where τ 0 is the curing time and C(t i , τ) is calculated according to the equations given by the corresponding prediction models.
(3) By substituting C(t i , τ) and λ i into equation (41), the relationship between a i (τ) and C(t i , τ) can be expressed as follows: Advances in Civil Engineering 5 (4) Finally, {a}m × 1 is determined using the leastsquares method, which is expressed as follows:

Numerical Analysis Method.
Regarding the creep and shrinkage of the CFST columns, the shrinkage of concrete is restrained by the steel tube, and the concrete core that contacts the steel tube is pulled. is results in tensile strain and a corresponding creep strain. Moreover, the creep deformation of the concrete core caused by the load is partially responsible for the load transfer to the steel tube. Out of the many analysis methods, the method of step-by-step recurrence can handle this complex process easily and can also be combined with FE software. us, this method has been widely used in practical engineering projects [12]. Equation (38) indicates that the creep strain at t n and t n+1 can be expressed as follows: Δσ 0 C t n , t 0 , σ t 0 + · · · + Δσ n− 1 C t n , t n− 1 , σ t n− 1 , Δσ 0 C t n+1 , t 0 , σ t 0 + · · · + Δσ n C t n+1 , t n , σ t n , By substituting equation (48) into equation (47), it can be seen that the creep strain increment between the adjacent time increment steps has a recurrence relationship, which can be expressed as follows: 6 Advances in Civil Engineering Δσ n a i t n ,

Program of Numerical Calculation.
A corresponding numerical calculation program was developed using the FORTRAN programming language, according to the analysis model of the 3D nonlinear creep and shrinkage effects, and the powerful secondary development function of ABAQUS.
Two user subroutines (USDFLD and UEXPAN) and one utility routine (GETVARM) are required to implement creep analysis through the ABAQUS secondary platform. e purpose of USDFLD is to determine the evolution of the concrete elastic modulus with age. First, the elastic modulus is defined as a function of field variables in the main program. Next, the field variables are specified in USDFLD at each time increment step. e function of UEXPAN is to calculate the strain increment caused by the creep and shrinkage at each time increment step. Because shrinkage is independent of stress, the shrinkage strain increment can be calculated according to the formulas provided by different models. Remarkably, when calculating the current creep strain for the creep strain increment using equation (49), it is necessary to call ω(i, n − 1), which is stored during a previous step, owing to the recurrence relation between the adjacent time increment steps. In this case, the solution-dependent state variables (SDVs) are required to store ω(i, n). e utility routine GETVRM is used to extract the variable value of the material points at each moment, including the stress and strain values and the SDVs. Subsequently, they transfer them to the UEXPAN subroutine. e

Program Validation.
Creep tests on short CFST columns under different stress levels have been conducted by Fu [6] and Lin [25], respectively. e specimen details, loadings, and stress levels are listed in Table 1. e mix proportion of the concrete is presented in Table 2, and the material parameters of the concrete core and steel tubes are presented in Tables 3 and 4. ese tests were carried out in laboratories, where the temperature was relatively constant. e strain gauges were located in the middle of the specimens. e equipment used to perform the creep test in these experiments is shown in Figure 5. Loading was achieved through hydraulic jacks in conjunction with the force sensors.
In this analysis, the prediction model B4 was chosen because this model classified creep either as drying creep or as basic creep. For the simulation conducted in this experiment, it was considered that the concrete core did not exchange humidity with the external environment by excluding the drying creep. Moreover, the moisture content decreased via cement hydration. e relative humidity of the concrete core was taken as 85%, according to the actual measurement results reported by Zhang et al. [26]. e relevant parameters of the B4 model are listed in Table 5. e analysis was implemented using the FE package ABAQUS 6.14.
e steel tubes were modeled by C3D8I, which can improve the convergence of nonlinear analysis. e concrete core and the steel-concrete interface were modeled by C3D8 and hard contact, respectively, as shown in Figure 6. In the absence of load-displacement data for the specimen, the dilation angle was calculated according to equation (29). e parameters of UCCDP are listed in Table 6, and the FE models are presented in Figure 6. e comparison between the measured data and the analysis results is presented in  e results revealed that the linear creep theory (L) and nonlinear creep theory (NL) have little difference when the stress level of the concrete core is lower (less than 0.6), and both theories could approximate the experimental results. However, for the specimens with a higher stress level (greater than 0.6), the results of the linear and nonlinear creep theory were obviously different. For the experiment conducted by Fu [6], after 110 days of loading, the result obtained using nonlinear creep theory was approximately 84% more than that obtained using linear creep theory. Moreover, the difference in relation to the experiment conducted by Lin [25] was 72%. It was observed that the calculated results of nonlinear creep theory are in agreement with the measured results. us, the reliability of the proposed model was confirmed.  [24] 0.38-I 108 × 4 × 330 377 0.38 Lin [24] 0.8-I 108 × 4 × 330 793 0.8 Note. N L is the long-term loading load and n c is the load ratio, which is the ratio of the initial stress of concrete core (σ 0 ) to the average cylinder compressive strength (f c (t 0 )),  Note. f c (t 0 ) is the average compressive strength of Φ150 × 300 mm specimen and E c (t 0 ) is the elastic modulus of the specimen.  Figure 5: Creep test apparatus: (a) conducted by Fu [6] and (b) Lin [25].  Figure 6: e FE models: (a) Fu [6] and (b) Lin [24].

Conclusions
is study can be summarized as follows: (1) A new theoretical framework for 3D nonlinear creep effect analysis was proposed to accurately analyze the nonlinear creep effect of short CFSTcolumns with an axial load, and the corresponding numerical analysis method was established. (2) A new damaged plasticity model for uniform confined concrete (UCCDP) was proposed, and the calibration methods of the key parameters were provided. (3) An analysis program was developed using the FORTRAN programming language based on the ABAQUS secondary platform. e results revealed that this method can easily realize nonlinear creep effect analysis for short CFST columns under axial compression with a certain amount of accuracy.
Moreover, it should be noted that the method proposed in this study is mainly applicable to the analysis of the concrete core's nonlinear creep effect under a uniform constraint. For long CFSTcolumns subjected to axial load, or short columns subjected to an eccentric load, the concrete core is no longer under a uniform constraint, and further investigation is required for a method of analyzing the nonlinear creep effect.

Data Availability
e data used to support the findings of this study are available from the corresponding author upon request.