Soil-structure interaction with time-dependent behaviour of both concrete and soil

Abstract

As the consideration of soil-structure interaction is increasingly being incorporated into design practice, there is a need for solutions that consider the most relevant aspects of the behaviour of both the concrete and the soil, some of which are time-dependent. The present work introduces a model that suits both the structure and the foundation – including the supporting soil – allowing the coupled analysis of hardening on cure, creep, shrinkage and cracking of concrete, and consolidation of soil. Both the structure and the foundation are modeled as one-dimensional finite elements. The time-dependent behaviour of the concrete and soil is modeled using Kelvin chains. For the structural elements, a mixed finite element formulation is used. Validation tests were conducted in the proposed modeling, comparing its results with available experimental and analytical data. The study of a reinforced concrete continuous beam supported by foundations on consolidating clay considering the time-dependent behaviour of the materials showed considerable changes in the effects of the soil-structure interaction.

Keywords
Soil/structure interaction; Time dependence; Concrete structures; Consolidation; Finite element method

Graphical Abstract

Notation

$A_s$ cross-sectional area

$\mathbf{a}$ element nodal displacement vector

${\mathbf{a}}_{\mathit{r}}$ rotation matrix

$B$ foundation width

$B_s$ first moment of the cross-sectional area

$\mathbf{b}$ internal forces interpolation function

${\mathbf{b}}_p$ particular internal forces interpolation function

$C_c$ compression index

$c_v$ coefficient of consolidation

$\mathbf{D}$ displacement vector

$D_{\mu }$ parameter of μth aging Kelvin unit

$D_{\mu ,k+1/2}$mid-step stiffness modulus

$E_c$ elastic modulus of consolidating soil

$E_i$ asymptotic modulus (age-independent)

${\stackrel{-}{E}}_k$ incremental modulus in step number $k$

$E_0$ initial value of the elastic modulus

$E_{\infty }$ final value of the elastic modulus

$E_s$ elastic modulus of steel

$\mathbf{e}$ cross-sectional strain vector

${\mathbf{e}}^{\mathit{\text{'}}\mathit{\text{'}}}$ creep and shrinkage cross-sectional strain vector

$e_{r\mu }$ dimensionless strain rate in unit μ

$e_o$ initial void ratio

${\mathbf{e}}_1$ non-cracked cross-sectional strain vector

${\mathbf{e}}_2$ cracked cross-sectional strain vector

$F_x$ nodal force in x direction

$F_y$ nodal force in y direction

$\mathbf{F}$ applied nodal force vector

${\mathbf{F}}_R$ resisting nodal force vector

$\mathbf{f}$ element loading vector

$f_{cm}$ mean concrete compressive strength

$H$ Heaviside step function

$\mathbf{H}$ transformation matrix for eliminating rigid body movements

$h$ thickness of the compressible layer

$h_d$ drainage distance in consolidation

$I_s$ second moment of the cross-sectional area

$i$ index for the iterative counter

$j$ creep function

$\mathbf{K}$ structural tangent stiffness matrix

$k$ index representing the time step

${\mathbf{k}}_{el}$ element tangent stiffness matrix without rigid-body modes

${\stackrel{-}{\mathbf{k}}}_s$ transformed cross-sectional stiffness matrix

${\stackrel{-}{\mathbf{k}}}_{1,s}$ non-cracked transformed cross-sectional stiffness matrix

${\stackrel{-}{\mathbf{k}}}_{2,s}$ cracked transformed cross-sectional stiffness matrix

$L$ length of foundation element

$l$ length of beam element

$m$ number of Kelvin units in a chain

$M$ bending moment, moment

$\widehat{M}$ bending moment due to strain

$M_p$ bending moment particular solution for non-zero $p_y$

$M_r$ cracking moment

$M_z$ nodal moment in z direction

$N$ normal force, axial force

$\widehat{N}$ normal force due to strain

$N_p$ normal force particular solution for non-zero $n_x$

$n_x$ loading per unit length of beam in x direction

$p$ excess pore-pressure

$p_y$ loading per unit length of beam in y direction

$q$ applied pressure

$\mathbf{q}$ element nodal internal forces vector

$\mathbf{R}$ residual force vector

$\widehat{\mathbf{s}}$ cross-sectional internal forces vector due to strain

${\widehat{\mathbf{s}}}^{\text{'}\text{'}}$ constrained cross-sectional internal force vector due to creep and shrinkage

${\mathbf{s}}^{\mathbf{p}}$ cross-sectional internal forces vector due to loading

$T$ consolidation time factor

$t$ time, current age of concrete

$t^{\text{'}}$ age at loading in a creep test

$t_k$ time of kth stress jump or kth time in a numerical solution

$U$ average degree of consolidation

$u$ axial displacement

$v$ vertical displacement

$W$ quadrature weight and length

$\mathbf{w}$ nodal displacements without rigid body movement/element deformation vector

$w_c$ consolidation settlement

$x_r$ position of the cracked moment

$y$ vertical distance from the fiber to the cross-sectional reference axis

$\beta$ coefficient for the influence of loading time on mean strain

${\beta }_{\mu }$μth auxiliary constant used by exponential algorithm

$\gamma$ slope of stress diagram

${\gamma }_s$ unit weight of the soil

$\mathrm{\Delta }{\epsilon }_k^{\text{'}\text{'}}$ strain increment due to creep in step number $k$

$\mathrm{\Delta }{\sigma }_{v0}^{\text{'}}$ average stress increase in a soil layer

$\epsilon$ strain

${\epsilon }_o$ strain at the cross-sectional reference point

${\epsilon }_{cr}$ creep strain

${\epsilon }_i$ instantaneous strain

${\epsilon }_s$ steel strain

${\epsilon }_{sh}$ shrinkage strain

${\epsilon }_T$ thermal strain

$\zeta$ interpolation coefficient

$\eta$ index referring to the value of the variable at the quadrature point

${\lambda }_{\mu }$μth auxiliary constant used by exponential algorithm

$\mu$ index referring to the Kelvin unit in a chain

$\sigma$ normal stress

${\sigma }_{v0}^{\text{'}}$ initial effective vertical stress

$\tau$ retardation time

$\psi$ cross-sectional curvature

1 INTRODUCTION

The interaction between the structure and its foundation is increasingly being considered and incorporated into design practice. Structural projects based on rigid supports are still made for small buildings or in cases where the foundations are certainly rigid (not settling). In structures founded on footings or rafts, or even on groups of piles that interact with each other (i.e. subject to group effect), the analysis including soil-structure interaction leads to different and more realistic results. In addition to changing the pattern of settlements, making them more uniform, it is possible to assess the redistribution of forces in the structure.

Foundation and structural standards around the world are beginning to include recommendations for interaction analyses both in working conditions and at failure, although the former is more usual. For the analysis, adequate models must be used, ideally in a coupled solution. Although being recognized as important factors in the behaviour of reinforced concrete structures, concrete shrinkage, hardening (aging of the elastic modulus), cracking and creep, as well as the presence of reinforcement are hardly included in soil-structure interaction analyses.

The main parameter to assess the implications of an interaction analysis is the relative stiffness between the structure and the foundation (Lee and Brown, 1972; Gusmão, 1990). Several studies were carried out to analyze the influence of different factors affecting this stiffness, such as three-dimensional effect (Majid and Cunnell, 1976), number of floors (Goschy, 1978; Lopes and Gusmão, 1991), construction sequence (Brown and Yu, 1986), soil nonlinearity (Jardine et al., 1986), soil constitutive model (Noorzaei et al., 1995), interface elements (Dalili et al., 2014), and even the presence of wall columns (Dalili et al., 2011) and brick panels (Hora, 2014).

Various studies considered the time-dependent soil behaviour (Viladkar et al., 1993), the viscoelastic behaviour of the structure (Gonçalves et al., 2016; Rosa et al., 2018) or both (Santa Maria et al., 1999; Fajman and Sejnoha, 2007). Nevertheless, no study could be found that also considered the effects of concrete hardening, cracking and tension-stiffening (increase in stiffness of a concrete section due to the transmission of stresses from rebars to the concrete between adjacent cracks) and the presence of reinforcement, all of which considerably influence the response of the structure (Espion and Halleux, 1990).

The present work presents a model that can be applied to both the structure and the foundation – including the soil – allowing a coupled analysis. The model allows the consideration, for the concrete, of hardening, shrinkage, creep, reinforcement, cracking and tension-stiffening. As for the soil, it enables the consideration of consolidation (optionally also secondary consolidation). The modeling uses, for both the structure and the foundation, one-dimensional finite elements (truss or beam elements). Kelvin chains are used to model the time-dependent behaviour of the materials, namely creep and consolidation. A mixed finite element formulation is used for the structural elements.

Validation tests were carried out after the models were programmed in a code referenced here by the aceonym SoilStruc, comparing its results with available experimental and analytical data. Then, a hypothetical case is studied. The case consists of a continuous beam supported by foundations on consolidating clay. Although simple, this case allows the examination the effects of the various phenomena on the behaviour of structure, alone and in combination.

This study shows that, although concrete hardening and reinforcement increase the stiffness of the structure, creep and cracking reduce it more significantly, causing both the uniformity of settlements and the redistribution of forces – expected in a soil-structure interaction analysis – to be less significant.

2 PROPOSED MODELING FOR TIME-DEPENDET BEHAVIOUR OF REINFORCED CONCRETE BEAMS

In buildings, continuous (hyperstatic) beams are the main elements responsible for the redistribution of forces between columns when differential settlements occur. Therefore, it is crucial that its behaviour be adequately modeled for a correct force-displacement relationship of the structure.

Although many studies on soil-structure interaction consider the behaviour of structural elements as linear elastic (e.g. Brown and Yu, 1986; Masih, 1993; Viladkar et al., 1993; Arapakou and Papadopoulos, 2020), it is known that reinforced concrete beams tend to crack and develop a non-linear response (Cosenza and Greco, 1990). In addition, due to shrinkage on cure and creep of the concrete, the force-displacement relationship varies over time (Ghali, 1986), increasing the complexity of the analysis.

Building structures are usually modeled using one-dimensional elements. Due to the slenderness of these elements, the strains due to shear can be neglected and the Euler-Bernoulli hypothesis can be assumed. In this way, all strains are null, except the axial strain which is given by:

$\epsilon =u^{\mathrm{\text{'}}}-y{v}^{\mathrm{\text{'}}\mathrm{\text{'}}}={\epsilon }_o+y\psi ={\mathbf{a}}_s\mathbf{e}$ (1)

where ${\epsilon }_o$ is the strain at the cross-sectional reference point, $y$ is the vertical distance from the fiber to the cross-sectional reference axis, $\psi$ is the cross-sectional curvature, ${\mathbf{a}}_s$ and $\mathbf{e}$ are defined as:

${\mathbf{a}}_s=\left\{\begin{array}{cc}1& y\end{array}\right\}$ (2)

$\mathbf{e}={\left\{\begin{array}{cc}{\epsilon }_o& \psi \end{array}\right\}}^{\mathrm{T}}$ (3)

Besides the classical parameters, such as the elastic modulus of the material, geometric properties of the cross-section and element length, the phenomena that most affect beam strain are concrete shrinkage, hardening, creep and tension-stiffening (Espion and Halleux, 1990). Such phenomena can be modeled at different levels. While the first three can be modeled at the material level, the latter has to be modeled at the cross-section level. With the relationship between internal forces and cross-sectional strains, it is sufficient to integrate the strains to obtain displacements at the extremities of the element, disregarding rigid body movements.

2.1 Constitutive laws of the materials

• (a) Concrete

For a uniaxial condition, under constant stresses $\sigma \left(t_0\right)$ applied at $t_0$, the total strain of the concrete $\epsilon \left(t\right)$ at time $t$ can be given by the sum of four parts:

$\epsilon \left(t\right)={\epsilon }_i\left(t_0\right)+{\epsilon }_{cr}\left(t\right)+{\epsilon }_{sh}\left(t\right)+{\epsilon }_T\left(t\right)$ (4)

where ${\epsilon }_i\left(t_0\right)$ is the instantaneous strain due to loading, ${\epsilon }_{cr}\left(t\right)$ is the creep strain, ${\epsilon }_{sh}\left(t\right)$ is the shrinkage strain, and ${\epsilon }_T\left(t\right)$ is the thermal strain.

For the in-service behaviour, the relation $\sigma \left(t_0\right)\mathrm{}\mathrm{x}\mathrm{}{\epsilon }_{cr}\left(t\right)$ can be considered linear. Furthermore, assuming there is no temperature variation, (Eq. 4) can be written as:

$\epsilon \left(t\right)-{\epsilon }_{sh}\left(t\right)=\sigma \left(t_0\right)J(t,t_0)$ (5)

where $J\left(t,t_0\right)$ is the creep function.

For a continuous loading function, the strain generated by the stress history is given by:

$\epsilon \left(t\right)-{\epsilon }_{sh}\left(t\right)={\int }_{t_0}^{t}J\left(t,t^{\mathrm{\text{'}}}\right)\mathrm{d}\sigma \left(t^{\mathrm{\text{'}}}\right)$ (6)

The step-by-step method was selected to model creep as it is applicable for evaluating any stress history. In this method, creep can be represented either in integral or in differential form. The second representation was adopted because it is computationally more efficient as it is not necessary to store the entire stress history, but only a fixed number of variables that are updated at each step (Bazant and Jirásek, 2018).

The creep function can either be obtained in tests or taken from standards and recommendations. However, for the differential formulation, the creep function needs to be properly represented. For the case of hardening viscoelastic materials, Bazant and Jirásek (2018) suggest approximating the function by the following Dirichlet series:

$J\left(t,t^{\mathrm{\text{'}}}\right)=\left[\frac{1}{E_{\mathrm{i}}}+\sum _{\mu =1}^m\frac{1-e^{-\frac{(t-t^{\mathrm{\text{'}}})}{{\tau }_{\mu }}}}{D_{\mu }\left(t^{\mathrm{\text{'}}}\right)}\right]H(t-t^{\mathrm{\text{'}}})$ (7)

where $E_i$ is the age-independent modulus, m is the number of Kelvin units in a chain, $\mu$ is the index referring to the Kelvin unit in a chain, $D_{\mu }$ is the parameter of μth hardening Kelvin unit, ${\tau }_{\mu }$ is the μth retardation time and $H$ is the Heaviside step function.

The least squares method can be applied to obtain the values of the $D_{\mu }$ functions. To this end, the values of the retardation time ${\tau }_{\mu }$, the age at loading $t^{\text{'}}$ and the time increment $t-t^{\text{'}}$ must be chosen properly to ensure both that the system is not ill-conditioned, and that the solution has good accuracy (Bazant and Wu, 1973). The value of $m$ is a function of the time interval for which the viscoelastic response is desired. The longer the range, the greater the number of Kelvin chains needed.

In the numerical solution, the problem is treated incrementally, and time is divided into intervals of size ${∆t}_k$. Integration algorithms lose numerical stability or accuracy when such an interval is greater than the shortest retardation time of the Kelvin chain. Using the exponential algorithm for hardening materials proposed by Bazant (1971), this limitation is overcome, allowing greater flexibility in defining time intervals.

The incremental constitutive relation using the second order exponential algorithm in the time interval ${∆t}_k=t_{k+1}-t_k$ can be written as:

$\mathrm{\Delta }{\sigma }_k={\stackrel{-}{E}}_k(\mathrm{\Delta }{\epsilon }_k-\mathrm{\Delta }{\epsilon }_k^{\mathrm{\text{'}}\mathrm{\text{'}}}-\mathrm{\Delta }{\epsilon }_{sh,k})$ (8)

where the incremental modulus ${\stackrel{-}{E}}_k$, the incremental creep strain at constant stress $\mathrm{\Delta }{\epsilon }_k^{\text{'}\text{'}}$ and the incremental shrinkage strain $\mathrm{\Delta }{\epsilon }_{sh,k}$ are given by:

${\stackrel{-}{E}}_k={\left(\frac{1}{E_{\mathrm{i}}}+\sum _{\mu =1}^m\frac{1-{\lambda }_{\mu ,k}}{D_{\mu ,k+1/2}}\right)}^{-1}$ (9)

$\mathrm{\Delta }{\epsilon }_k^{\mathrm{\text{'}}\mathrm{\text{'}}}=\sum _{\mu =1}^m(1-{\beta }_{\mu k})e_{r\mu ,k}$ (10)

$\mathrm{\Delta }{\epsilon }_{sh,k}={\epsilon }_{sh}\left(t_{k+1},t_s\right)-{\epsilon }_{sh}\left(t_k,t_s\right)$ (11)

where $D_{\mu ,k+1/2}$ is the mid-step stiffness modulus, and the factors ${\beta }_{\mu ,k}$, ${\lambda }_{\mu ,k}$ and $e_{r\mu ,k}$ are given by:

${\beta }_{\mu ,k}=e^{-\mathrm{\Delta }{t}_k/{\tau }_{\mu }}$ (12)

${\lambda }_{\mu ,k}=(1-{\beta }_{\mu ,k})\frac{{\tau }_{\mu }}{\mathrm{\Delta }{t}_k}$ (13)

where $e_{r\mu ,k}$ is an internal variable representing the dimensionless strain rate, whose initial value is null and whose updated value at the end of each step is given by:

$e_{r\mu ,k+1}=\frac{{\lambda }_{\mu ,k}}{D_{\mu ,k+1/2}}\mathrm{\Delta }{\sigma }_k+{\beta }_{\mu ,k}{e}_{r\mu ,k}$ (14)

The functions proposed by Comite Euro-International du Beton (CEB-FIP) in Model Code 1990 (1993) were used to model the increase in the elastic modulus (hardening), creep, and shrinkage.

• (b) Steel

The behaviour of steel was considered as linear elastic for both compression and tension, as it is assumed that for working conditions, the yield stress is not reached. Therefore, the increase in stress in the steel $\mathrm{\Delta }{\sigma }_s$ is given by:

$\mathrm{\Delta }{\sigma }_s=E_s\mathrm{\Delta }{\epsilon }_s$ (15)

where $E_s$ is the elastic modulus of the steel and $\mathrm{\Delta }{\epsilon }_s$ is the steel strain increment.

2.2 Constitutive law of the cross-section

For reinforced concrete elements in working conditions, the tension region can either be cracked or not. So, the section response is divided into non-cracked and cracked. In the first one, it is assumed that concrete withstands up to the tensile strength. In the second one, it is assumed that the tensile strength is null. Furthermore, it is well known that the concrete in tension between cracks has a significant contribution to the stiffness of the element (Beeby et al., 2005). This phenomenon is known as tension-stiffening. To consider such phenomenon, the proposal by Fédération Internationale du Béton (fib) in Fib Model Code for Concrete Structures 2010 (2013) was adopted. In the proposed formulation, the strain of the section at instant k is given by:

${\mathbf{e}}_k=\left(1-\zeta \right){\mathbf{e}}_{1,k}+\zeta {\mathbf{e}}_{2,k}$ (16)

where ${\mathbf{e}}_{1,k}$ is the non-cracked cross-sectional strain vector, ${\mathbf{e}}_{2,k}$ is the cracked cross-sectional strain vector, and $\zeta$ is an interpolation coefficient given by:

$\zeta =\left\{\begin{array}{cc}0& se\mathrm{}M\le M_r\\ 1-\beta {\left(\frac{M_r}{M}\right)}^2& se\mathrm{}M\ge M_r\end{array}\right.$ (17)

where $\beta$ is a coefficient that considers the influence of the loading time on the mean strain, $M_r$ is the cracking moment and $M$ represents the bending moment.

According to Beeby et al. (2005), the value of 0.5 for $\beta$ can be used in the calculation of creep strain, since the loss of tension-stiffening occurs in a shorter period of time when compared to creep. To perform the analysis over time, Eq. (16) must be written in incremental form. Using the Taylor series expansion and removing the higher order terms, it follows:

$\mathrm{\Delta }{\mathbf{e}}_k={\mathbf{e}}_{k+1}-{\mathbf{e}}_k\cong \left(1-\zeta \right)\mathrm{\Delta }{\mathbf{e}}_{1,k}+\zeta \mathrm{\Delta }{\mathbf{e}}_{2,k}+\zeta \mathrm{\text{'}}({\mathbf{e}}_{2,k}-{\mathbf{e}}_{1,k})\mathrm{\Delta }M$ (18)

where

${\zeta }^{\mathrm{\text{'}}}=2\beta M_r^2/M^3$ (19)

Using the Euler-Bernoulli hypothesis, it is possible to define an expression relating the increase in internal forces with the increase in the strain of the section. In the case of viscoelastic materials, such expression is:

$\mathrm{\Delta }{\widehat{\mathbf{s}}}_k+\mathrm{\Delta }{\widehat{\mathbf{s}}}_k^{\mathrm{\text{'}}\mathrm{\text{'}}}={\stackrel{-}{\mathbf{k}}}_{s,k}\mathrm{\Delta }{\mathbf{e}}_k$ (20)

where $\mathrm{\Delta }{\widehat{\mathbf{s}}}_k$ is the increment of cross-sectional internal forces vector, $\mathrm{\Delta }{\widehat{\mathbf{s}}}_k^{\text{'}\text{'}}$ is the increment of constrained cross-sectional internal forces vector due to creep and shrinkage and ${\stackrel{-}{\mathbf{k}}}_{s,k}$ is the stiffness matrix of the transformed cross-section. For the case of a section with n materials, the values of $\mathrm{\Delta }{\widehat{\mathbf{s}}}_k^{\text{'}\text{'}}$ and ${\stackrel{-}{\mathbf{k}}}_{s,k}$ are given by:

$\mathrm{\Delta }{\widehat{\mathbf{s}}}_k^{\mathrm{\text{'}}\mathrm{\text{'}}}=\sum _{m=1}^n\left({\int }_{A_m}^{\mathrm{}}{\mathbf{a}}_s^{\mathrm{T}}{\stackrel{-}{E}}_{k,m}{\mathbf{a}}_{s}dA\right)\mathrm{\Delta }{\mathbf{e}}_{k,m}^{\mathrm{\text{'}}\mathrm{\text{'}}}$ (21)

${\stackrel{-}{\mathbf{k}}}_{s,k}=\sum _{m=1}^n\left({\int }_{A_m}^{\mathrm{}}{\mathbf{a}}_s^{\mathrm{T}}{\stackrel{-}{E}}_{k,m}{\mathbf{a}}_{s}dA\right)=\sum _{m=1}^n{\stackrel{-}{E}}_{k,m}\left[\begin{array}{cc}{A}_{s,m}& B_{s,m}\\ B_{s,m}& I_{s,m}\end{array}\right]$ (22)

where $\mathrm{\Delta }{\mathbf{e}}_{k,m}^{\text{'}\text{'}}$ is the increment of cross-sectional strain due to creep and shrinkage, $A_{s,m}$ is the cross-sectional area, $B_{s,m}$ is the section first moment of area and $I_{s,m}$ is the section second moment of area.

Considering the same stress increment $\mathrm{\Delta }{\widehat{\mathbf{s}}}_k$, the strain of the non-cracked section $\mathrm{\Delta }{\mathbf{e}}_{1,k}$ will be different from the strain of the cracked section $\mathrm{\Delta }{\mathbf{e}}_{2,k}$, since the consideration or not of the concrete tensile strength affects the properties of area $A_{s,m}$, $B_{s,m}$ and $I_{s,m}$ of the cross-section. Thus, the response of the non-cracked and cracked section can be given, respectively, by:

$\mathrm{\Delta }{\widehat{\mathbf{s}}}_k+\mathrm{\Delta }{\widehat{\mathbf{s}}}_{1,k}^{\mathrm{\text{'}}\mathrm{\text{'}}}={\stackrel{-}{\mathbf{k}}}_{1,s,k}\mathrm{\Delta }{\mathbf{e}}_{1,k}$ (23)

$\mathrm{\Delta }{\widehat{\mathbf{s}}}_k+\mathrm{\Delta }{\widehat{\mathbf{s}}}_{2,k}^{\mathrm{\text{'}}\mathrm{\text{'}}}={\stackrel{-}{\mathbf{k}}}_{2,s,k}\mathrm{\Delta }{\mathbf{e}}_{2,k}$ (24)

In the cracked section, it is necessary to define the height of the neutral line to calculate the properties of area. Such analysis can be done using an iterative process (Kawakami and Ghali, 1996). Furthermore, because of both creep and concrete shrinkage, the neutral line position tends to change over time. However, if this height is supposed constant over time, the analysis is simplified and errors resulting from this hypothesis are not significant (Ghali, 1986). In a comparison carried out by Torres (2001), the difference in section curvature was only 2%.

2.3 Beam element formulation

The set of equations that must be solved for a beam element are related to the strain-displacement relationship (Eq. 1), constitutive cross-section relationship (Eqs. 16, 23, 24) and equilibrium equations given by:

$\frac{\partial N}{\partial x}+n_x=0$ (25)

$\frac{{\partial }^{2}M}{\partial z^2}+p_y=0$ (26)

where $n_x$ and $p_y$ are the loadings per unit length of the beam in the x and y directions, respectively.

As can be seen in Eq. (16), due to tension-stiffening, the response of the cross-section is non-linear when the cracking moment is exceeded. To model such behaviour with finite elements based on the displacement field, it is necessary to subdivide the element so that the nonlinear behaviour of the material is satisfactorily approximated. When using elements based on the force field, such refinement is not necessary and, thus, elements with less degrees of freedom can be used, leading to greater computational efficiency.

In the present work, the mixed finite element formulation proposed by Taylor et al. (2003) was used. After discretizing the fields and calculating the integrals along the element, the equations obtained correspond to the balance of forces, to the displacement at the ends of the element, and to the balance of the internal forces along the element. Putting the three equations in matrix incremental form, follows:

$\left[\begin{array}{ccc}0& {\mathbf{H}}^{\mathbf{T}}& 0\\ \mathbf{H}& 0& -{\mathbf{b}}_{\mathbf{\eta }}^{\mathbf{T}}\\ 0& -{\mathbf{b}}_{\mathbf{\eta }}& {\stackrel{-}{\mathbf{k}}}_{\mathit{s},\eta }\end{array}\right]\left\{\begin{array}{c}∆\mathbf{a}\\ ∆\mathbf{q}\\ ∆{\mathbf{e}}^{\mathit{\eta }}\end{array}\right\}=\left\{\begin{array}{c}\mathbf{F}-{\mathbf{H}}^{\mathbf{T}}\mathbf{q}\\ {\mathbf{b}}_{\eta }^{\mathbf{T}}{\mathbf{e}}^{\mathit{\eta }}-\mathbf{H}\mathbf{a}\\ {{\mathbf{b}}_{\mathbf{\eta }}^{\mathbf{T}}\mathbf{q}+\mathbf{s}}_{\mathbf{\eta }}^{\mathbf{p}}-{\widehat{\mathbf{s}}}_{\mathbf{\eta }}\left({\mathbf{e}}^{\mathit{\eta }}\right)\end{array}\right\}$ (27)

where index $\eta$ indicates the value of the variable at a point of the quadrature, $\mathbf{a}$ is the vector of nodal displacements, $\mathbf{q}$ is the vector of nodal internal forces, $\mathbf{H}$ is the transformation matrix for eliminating rigid body movements, $\mathbf{b}$ is the interpolation function of the nodal internal forces, $\mathbf{F}$ is the nodal forces vector, ${\mathbf{s}}^{\mathbf{p}}$ are the internal forces due to loading, and $\widehat{\mathbf{s}}$ are the internal forces due to cross-sectional strain. These variables are given by:

$\mathbf{a}={\left\{\begin{array}{cc}\begin{array}{ccc}{u}^1& v^1& {\theta }^1\end{array}& \begin{array}{ccc}{u}^2& v^2& {\theta }^2\end{array}\end{array}\right\}}^T$ (28)

$\mathbf{q}={\left\{\begin{array}{ccc}N& M^1& M^2\end{array}\right\}}^T$ (29)

$\mathbf{H}=\left\{\begin{array}{cc}\begin{array}{ccc}-1& 0& 0\\ 0& \frac{-1}{l}& -1\\ 0& \frac{1}{l}& 0\end{array}& \begin{array}{ccc}1& 0& 0\\ 0& \frac{1}{l}& 0\\ 0& \frac{-1}{l}& 1\end{array}\end{array}\right\}$ (30)

${\mathbf{b}}_{\mathbf{\eta }}=\left[\begin{array}{ccc}1& 0& 0\\ 0& \frac{1-{\xi }_{\eta }}{2}& \frac{1+{\xi }_{\eta }}{2}\end{array}\right]W_{\eta }$ (31)

$\mathbf{F}={\left\{\begin{array}{cc}\begin{array}{ccc}{F}_x^1& F_y^1& M_z^1\end{array}& \begin{array}{ccc}{F}_x^2& F_y^2& M_z^2\end{array}\end{array}\right\}}^T-{\left\{\begin{array}{cc}\begin{array}{ccc}{-N}_p^1& M_p^{\mathrm{\text{'}}1}& {-M}_p^1\end{array}& \begin{array}{ccc}{-N}_p^2& {-M}_p^{\mathrm{\text{'}}2}& M_p^2\end{array}\end{array}\right\}}^T$ (32)

${\mathbf{s}}_{\mathbf{\eta }}^{\mathbf{p}}=\left\{\begin{array}{c}{N}_p^{\eta }{W}_{\eta }\\ M_p^{\eta }{W}_{\eta }\end{array}\right\}$ (33)

${\widehat{\mathbf{s}}}_{\mathbf{\eta }}=\left\{\begin{array}{c}{\widehat{N}}^{\eta }{W}_{\eta }\\ {\widehat{M}}^{\eta }{W}_{\eta }\end{array}\right\}$ (34)

where $l$ is the element length, $W_{\eta }$ is the quadrature weight and length, and $-1\le \xi \le 1$.

Substituting Equations (18), (23) and (24) in Equation (27), the incremental formulation for a viscoelastic material considering the contribution of the concrete between cracks is:

$\left[{\begin{array}{c}0\\ H\\ 0\\ 0\end{array}\begin{array}{c}{\mathbf{H}}^{\mathbf{T}}\\ {\mathit{\zeta }}^{\text{'}}({\mathbf{e}}_2^{\mathbf{\eta }}-{\mathbf{e}}_1^{\mathbf{\eta }}){\mathbf{b}}_{\mathbf{\eta }}\\ -{\mathbf{b}}_{\mathbf{\eta }}\\ -{\mathbf{b}}_{\mathbf{\eta }}\end{array}\begin{array}{c}0\\ -(1-\zeta ){\mathbf{b}}_{\mathbf{\eta }}^{\mathbf{T}}\\ {\stackrel{-}{\mathbf{k}}}_{1,\mathit{s},\mathit{\eta }}\\ 0\end{array}\begin{array}{c}0\\ -\zeta {\mathbf{b}}_{\mathbf{\eta }}^{\mathbf{T}}\\ 0\\ {\stackrel{-}{\mathbf{k}}}_{2,\mathit{s},\mathit{\eta }}\end{array}}^{}\right]\left\{\begin{array}{c}\begin{array}{c}∆\mathbf{a}\\ ∆\mathbf{q}\end{array}\\ \begin{array}{c}\mathrm{\Delta }{\mathbf{e}}_1^{\mathbf{\eta }}\\ \mathrm{\Delta }{\mathbf{e}}_2^{\mathit{\eta }}\end{array}\end{array}\right\}=\left\{\begin{array}{c}\begin{array}{c}\mathbf{F}-{\mathbf{H}}^{\mathbf{T}}\mathbf{q}\\ {{\left(1-\mathbf{\zeta }\right)\mathbf{b}}_{\mathbf{\eta }}^{\mathbf{T}}{\mathbf{e}}_1^{\mathbf{\eta }}+\mathbf{\zeta }\mathbf{b}}_{\mathbf{\eta }}^{\mathbf{T}}{\mathbf{e}}_2^{\mathbf{\eta }}-\mathbf{H}\mathbf{a}\end{array}\\ \begin{array}{c}{{\mathbf{b}}_{\mathbf{\eta }}^{\mathbf{T}}\mathbf{q}+\mathbf{s}}_{\mathbf{\eta }}^{\mathbf{p}}-{\widehat{\mathbf{s}}}_{\mathbf{\eta }}\left({\mathbf{e}}^{\mathit{\eta }}\right)+{\mathrm{\Delta }\widehat{\mathbf{s}}}_{1,\mathit{k},\mathit{\eta }}^{\text{'}\text{'}}\\ {{\mathbf{b}}_{\mathbf{\eta }}^{\mathbf{T}}\mathbf{q}+\mathbf{s}}_{\mathbf{\eta }}^{\mathbf{p}}-{\widehat{\mathbf{s}}}_{\mathbf{\eta }}\left({\mathbf{e}}^{\mathit{\eta }}\right)+{\mathrm{\Delta }\widehat{\mathbf{s}}}_{2,\mathit{k},\mathit{\eta }}^{\text{'}\text{'}}\end{array}\end{array}\right\}$ (35)

Note that, as the constitutive relationship is linear and the force field along the element is known, for both non-cracked and cracked sections it is possible to calculate the respective strain fields along the element and, therefore, the fields related to every element of the Kelvin chain.

The position of the cracked region changes depending on the load level and, thus, the positions of the integration points will also be modified, making it necessary to determine the values of the internal variables at the new points. Since the strain field is known, it is possible to calculate strains at any point. The solution algorithm used is given in appendix.

2.4 Solution algorithm

The solution of the structural system is obtained by iterations. The implemented algorithm follows the steps:

• 1. It is assumed that the process starts after the convergence of the previous step, which is represented by index k. Increment time according to $t_{k+1}=t_k+{∆t}_k$, nodal forces to ${\mathbf{F}}_{k+1}={\mathbf{F}}_k+{∆\mathbf{F}}_k$ and loads along the element to ${\mathbf{f}}_{k+1}={\mathbf{f}}_k+{∆\mathbf{f}}_k$;

• 2. For each material with viscoelastic behaviour, calculate: $D_{\mu ,k+1/2}$ , ${\beta }_{\mu ,k+1}$ , ${\lambda }_{\mu ,k+1}$ , ${\stackrel{-}{E}}_{k+1}$ , ${∆\mathbf{e}}_{k+1}^{\mathrm{\text{'}}\mathrm{\text{'}}}$;

• 3. For each cross-section, calculate: ${\stackrel{-}{\mathbf{k}}}_{s,k+1}\mathrm{},\mathrm{}∆{\widehat{\mathbf{s}}}_{k+1}^{\mathrm{\text{'}}\mathrm{\text{'}}}\mathrm{},\mathrm{}{M}_{r,k+1}$;

• 4. In order to determine the state of the structure subject to the force vectors ${\mathbf{F}}_{k+1}$ and ${\mathbf{f}}_{k+1}$, an iterative solution strategy is started aiming at reaching the global equilibrium with an iterative counter $i=0$; the structure state in k is used as initial guess;

• 5. For each element, assuming that all nodal displacements increments are null, calculate: $x_{r,i}\mathrm{},{\mathbf{k}}_{el,i}\mathrm{},{∆\mathbf{q}}_i\mathrm{},{\mathbf{K}}_{el,i}\mathrm{},\mathrm{}{\mathbf{F}}_{R,el,i}$ ;

• 6. Create the stiffness matrix and the resisting nodal force vector and calculate the nodal displacements: $∆{\mathbf{D}}_{i+1}={\left[{\mathbf{K}}_i\right]}^{-1}\left({\mathbf{F}}_{k+1}-{\mathbf{F}}_{R,i}\right)$;

• 7. Using the displacement increment in each element $\mathrm{\Delta }{\mathbf{D}}_{el,i+1}$, determine the displacements in local coordinates without rigid body movement: ${\mathbf{w}}_{i+1}={\mathbf{w}}_i+\mathbf{H}{\mathbf{a}}_r\mathrm{\Delta }{\mathbf{D}}_{el,i+1}$;

• 8. Using the displacement increment in each element $∆{\mathbf{w}}_i$, calculate the force vector increment: ${∆\mathbf{q}}_i={\mathbf{k}}_{el,i}\left\{∆{\mathbf{w}}_i-{\mathbf{b}}^T{\left[{\stackrel{-}{\mathbf{k}}}_{s,k+1}\right]}^{-1}\left({\mathbf{b}}_p\mathrm{\Delta }{\mathbf{f}}_{el,k+1}+\mathrm{\Delta }{\mathbf{s}}_{k+1}^{\mathrm{\text{'}}\mathrm{\text{'}}}\right)\right\}$;

• 9. Start iteration counter $j=0$ and use system state at $i$ as initial attempt;

• 10. Update the values of the nodal forces applied to the elements: ${\mathbf{q}}_j={\mathbf{q}}_i+{∆\mathbf{q}}_j$;

• 11. For each force field ${\mathbf{q}}_j$, calculate: $x_{r,j+1}\mathrm{},{\mathbf{w}}_{j+1}\mathrm{},{\mathbf{k}}_{el,j+1}$;

• 12. Calculate the nodal force increment using line search to determine the value of α that minimizes the norm $\left|{\mathbf{w}}_{i+1}-{\mathbf{w}}_{j+1}\right|$: ${∆\mathbf{q}}_{j+1}={\mathbf{k}}_{el,j+1}\left({\mathbf{w}}_{i+1}-{\mathbf{w}}_{j+1}\right)$ e ${∆\mathbf{q}}_{ls,j+1}=\alpha {∆\mathbf{q}}_{j+1}$.

• If the residual strain norm is less than the specified tolerance or the number of iterations is greater than the established maximum limit, use the states in $j+1$ for the states in $i+1$ and proceed to the next step; otherwise, increment $j$ and return to step 10;

• 13. Considering the displacement increment $∆{\mathbf{D}}_{i+1}$, calculate the stiffness matrix of the structure ${\mathbf{K}}_{i+1}$ and the resisting nodal force vector ${\mathbf{F}}_{R,i+1}$;

• 14. Calculate the residual of forces ${\mathbf{R}}_{i+1}={\mathbf{F}}_{k+1}-{\mathbf{F}}_{R,i+1}$. if the residual norm is less than the established tolerance or the number of iterations is greater than the maximum limit, update the state variables in $k+1$ with the values of the variables in $i+1$ and proceed to step 15; otherwise, increment in $i$ and return to step 6;

• 15. Calculate the strain increment $\mathrm{\Delta }{\mathbf{e}}_k$ and the total strain ${\mathbf{e}}_{k+1}={\mathbf{e}}_k+\mathrm{\Delta }{\mathbf{e}}_k$;

• 16. For each viscoelastic material of the cross-section of each element, calculate: ${∆\sigma }_k\mathrm{},{∆\gamma }_k\mathrm{},{\epsilon }_{o,\mu ,k+1}\mathrm{},{\psi }_{\mu ,k+1}\mathrm{},\mathrm{}{\epsilon }_{o,k+1}$ e ${\psi }_{k+1}$;

• 17. Increment the value of k and go back to step 1.

3 MODELING FOUNDATION SETTLEMENTS BY CONSOLIDATION

The magnitude of foundation settlements is fundamental for the analysis of the soil-structure interaction. In the proposed modeling, settlement development over time of a foundation derives from soil consolidation and a foundation element will be represented by a truss element (vertical) whose response over time is equivalent to that of the foundation. In this manner, the implementation of such analysis in frame-structure programs, which are used in most design offices, is simplified.

3.1 Long-term settlement calculation

The final settlement of a normally consolidated clay under 1D conditions can be obtained with Terzaghi’s formulation:

$w_c=\frac{C_c}{1+e_0}log\left(\frac{{\sigma }_{v0}^{\mathrm{\text{'}}}+\mathrm{\Delta }{\sigma }_{v0}^{\mathrm{\text{'}}}}{{\sigma }_{v0}^{\mathrm{\text{'}}}}\right)h$ (36)

where $C_c$ is the compression index, $e_0$ is the initial void ratio, ${\sigma }_{v0}^{\text{'}}$ is the initial effective vertical stress, $\mathrm{\Delta }{\sigma }_{v0}^{\text{'}}$ is the stress increase (average in the layer) and $h$ is the layer thickness.

If the stress increment is small in relation to the initial stress, the following approximation can be made so as to linearize the function:

$w_c=\frac{C_c}{1+e_0}\frac{h}{\mathrm{l}\mathrm{n}\left(10\right)}\frac{{\mathrm{\Delta }\sigma }_{v0}^{\mathrm{\text{'}}}}{{\sigma }_{v0}^{\mathrm{\text{'}}}}$ (37)

In order to use a truss element to represent a consolidating clay, it is necessary to define its area $A$, its length l, and the elastic modulus of the material. The length can be adopted as the clay layer thickness. The area can be defined as a function of the average stress increase in the layer, $\mathrm{\Delta }{\sigma }_{v0}^{\text{'}}$, of the mean applied pressure $q$, and of the dimensions $B$ and $L$ of the foundation with:

$A=\frac{qBL}{\mathrm{\Delta }{\sigma }_{v0}^{\mathrm{\text{'}}}}$ (38)

The average stress increase in the consolidating layer can be calculated by the Theory of Elasticity. The average stress increase can be calculated from the profile of stress increases obtained with Boussinesq’s equation.

Based on Eq. (37), the elastic modulus $E_c$ relative to consolidation can be defined as:

$E_c=\frac{(1+e_0)\mathrm{l}\mathrm{n}\left(10\right){\sigma }_{v0}^{\mathrm{\text{'}}}}{C_c}$ (39)

3.2 Settlement development over time

The development of settlement by consolidation can be modelled as a viscoelastic behaviour (Fajman and Sejnoha, 2007). For the 1D case, the consolidation equation is:

$\frac{\partial p}{\partial t}=c_v\frac{{\partial }^{2}p}{\partial z^2}$ (40)

where $p$ is the excess pore-pressure and $c_v$ is the coefficient of consolidation. Assuming that both upper and lower boundaries are draining and 2$h_d$ apart, and that the load application is instantaneous, the settlement $w_c$ can be calculated by:

$w_c\left(t\right)=w_c\left(\mathrm{\infty }\right)\sum _{n=\mathrm{1,3},5,\dots }^{\mathrm{\infty }}\frac{8}{n^2{\pi }^2}\left\{1-e^{-\frac{n^2{\pi }^2}{4h_d^2}{c}_{v}t}\right\}$ (41)

Comparing Eq. (41) with the Dirichlet series (Eq. 7), it can be observed that the term $E_i$ tends to infinity, that is, the material does not present instantaneous responses to loading. The terms $E_{\mu }$ and ${\tau }_{\mu }$ can be written as:

$E_{\mu }=\frac{8}{n^2{\pi }^2}$ (42)

${\tau }_{\mu }=\frac{4h^2}{n^2{\pi }^2{c}_v}$ (43)

The analogy described above allows the rational determination of the number of elements in the Kelvin chain and their respective parameters to accurately represent soil consolidation. It is interesting to note that the delay time value ${\tau }_{\mu }$ decreases with $n$, and that it is the smallest value of ${\tau }_{\mu }$ that indicates the capacity of the model to accurately represent the shortest time interval under analysis. Thus, considering that a finite number $m$ of elements in the Kelvin chain will be used in the simulation, it is necessary to introduce the term $E_i$, which is given by:

$E_{\mathrm{i}}={\left(1-\sum _{n=1}^{2m-1}\frac{8}{{\pi }^2{n}^2}\right)}^{-1}$ (44)

As the viscous model used for the soil is the same as that used for the concrete, it is possible to compare the values of the delay time of each model to define how many Kelvin elements should be used and which minimum time step should be initially adopted in the analysis.

4 VALIDATION TESTS

4.1 Experimental validation of concrete creep

For the experimental validation of creep, Benchmark no. 1 was used, with tests of one span reinforced concrete beams conducted by Subcommittee 3 of the Technical Committee TC 114 of RILEM (Espion, 1993).

The characteristics of the beams and loading can be seen in Figure 1 and Table 1. A total of 7 beams were tested with 5 load levels, 2 of them repeated in order to assess the consistency of the results. The main objective of the tests was to determine the development of the displacement at the middle of the span over time.

Figure 1
Loading and beam cross-section (dimensions in mm).

Table 1
Applied loads and average compressive strength of concrete in each test.

Beam	C11	C12	C22	C13	C14	C24	C15
P (kN)	5.77	12.19	12.19	18.61	25.04	25.04	31.45
fcm (MPA)	28.8	29.4	32.9	30.9	29.4	32.0	29.3

The beams were cured under constant moisture and temperature, thus, allowing to disregard thermal variation. The loading was applied instantly at 28 days and kept constant until the test was concluded.

Some parameters required for modeling had to be estimated, based on standards and recommendations, as creep tests of the concrete were not available. Thus, a perfect match between the experimental results and the numerical predictions should not expect, but rather the latter should reproduce the response of the structure with acceptable accuracy.

The experimental and numerical results can be compared in Figure 2. For the non-cracked case (C11) or the case with higher loading (C15), there is a better agreement with test results than for the other cases, in which the size of the crack formation region is larger. Overall, it can be concluded that the model provides satisfactory results.

Figure 2
Comparison of experimental results with numerical predictions.

4.2 Validation of soil consolidation modeling

In order to verify the accuracy of the method presented in item 3.2, the analytical solution was compared with results using 1, 2, 5 and 7 elements in the Kelvin chain. Figure 3 presents the percentage of settlement as a function of the time factor, T, for the analytical solution (U) and for different numbers of elements in the Kelvin chain (U1, U2, U5 and U7). The parameters were adjusted using Eq. (44) as to obtain 100% of settlement at the end of the analysis. Note that for $T>0.1$, there is no difference between the models used, whereas for smaller values of T, a significant difference between U1 and U can be observed. The difference between U5 and U7 is very small and both models were able to produce the expected response with good accuracy.

Figure 3
Percentage of settlement as a function of time factor for exact solution (U) and for solutions using Kelvin chains.

5 STUDY OF A CONTINUOUS BEAM ON SETTLING FOOTINGS DUE TO SOIL CONSOLIDATION

To assess time-dependent effects on soil-structure interaction, a continuous beam with 3 supports on footings was analyzed. This case was chosen as the beams of a structure are the main elements responsible for load redistribution, and for keeping the problem simple in order to evaluate how each of the intervening phenomena influences the problem.

5.1 The problem in study

The structure in study is shown in Figure 4a. Circular footings had dimensions defined by for an allowable stress of 50 kN/m², with the loads obtained with rigid supports: 112.5 kN and 375 kN. The ground profile is shown in Figure 4b. The water table is at the base of the sand, and the unit weights of the soils are ${\gamma }_{sand}=19\mathrm{k}\mathrm{N}/\mathrm{m}³$ and ${\gamma }_{clay}=14.7\mathrm{k}\mathrm{N}/\mathrm{m}³$. The initial vertical effective stress in the middle of the clay layer is ${\sigma }_{vo}^{\text{'}}=113\mathrm{k}\mathrm{N}/\mathrm{}\mathrm{m}²$. The stress increase profile in the soft clay layer – obtained with Boussinesq’s Equation – is shown in Figure 4b.

Figure 4
(a) Analyzed structure and (b) ground profile with the vertical stress increases.

To predict settlements due to the consolidation of the clay, the vertical stress increases along its thickness were considered. Performing numerical integration, the average stress increases are 1.91 kN/m² and 5.76 kN/m² under the side and center footings, respectively. Considering a normally consolidated clay and with $C_c=1.0$ and $e_0=2.5$, the final settlements are 1.24 cm and 3.70 cm, respectively. The contribution of the dense sand to the settlements of the footings can be neglected and the clay consolidation settlements will be the (long term) footing settlements.

The structural model is shown in Figure 5, where the horizontal element represents the beam and the vertical elements represent the footings (more specifically, the soil behaviour under the footing). The characteristics of the material and cross-section of the beam are shown in Table 2.

Figure 5
Structural model.

Table 2
Material and cross-section characteristics of the beam.

E (GPa)	b (mm)	h (mm)	A (mm2)	I (mm4)
25	300	600	180,000	5.4x109

The cross-sectional area of the side vertical element is 65.1 m² and of the center element is 58.9 m², and correspond to the average stress increase in the soft clay layer, i.e., the area of the elements was adopted so that the stresses remain equal to those indicated above (1.91 kN/m² and 5.76 kN/m²).

5.2 Influence of the viscoelastic behaviour of materials

To assess the effect of the viscoelastic behaviour of materials − namely concrete creep and soil consolidation − 6 cases were analyzed, as indicated in Table 3. E₀ represents the elastic modulus of the material at the initial moment, E_∞ represents the modulus at infinite time, and the viscoelastic material presents a transition from E₀ to E_∞ over time. In relation to concrete, the instantaneous elastic behaviour is represented by cases 1, 2 and 4. Creep is considered in cases 3 and 5; and the elastic behaviour after material creep is considered in case 6. For the soil, the elastic behaviour after consolidation is considered in cases 2, 3 and 6; and consolidation is considered in cases 4 and 5.

Table 3
Behaviour of the material for the different cases analyzed.

Case	1	2	3	4	5	6
Concrete	$E_0$	$E_0$	$E_0\to E_{\infty }$	$E_0$	$E_0\to E_{\infty }$	$E_{\infty }$
Soil	$Rigid$	$E_{\infty }$	$E_{\infty }$	$E_0\to E_{\infty }$	$E_0\to E_{\infty }$	$E_{\infty }$

Figure 6 shows (a) the differential settlement, (b) the reaction force at the central support (normal force in the element representing the foundation), and (c) the bending moment in the beam at the central support.

Figure 6
Result for the central support, over time, for the 6 cases analyzed.

At the beginning of the analysis, case 3 is equal to 2, and cases 4 and 5 are both equal to 1. It is also noted that there is a smooth transition over time. Hence, at the end of the analysis, case 4 is equal to case 2, and cases 3 and 5 are both equal to case 6. Cases 4 and 2 present the same result at the end of the analysis because in both cases the material of the beam is elastic, and that the only difference is the development of the supports’ settlement over time. The equality of cases 3, 5 and 6 indicates that the viscoelastic material presents a smooth transition from instantaneous behaviour to the behaviour at infinity time. Therefore, both the behaviour at the beginning and at the end of the analysis can be obtained using the respective elastic modulus, rendering unnecessary an analysis over time.

Knowing that the differential settlement calculated with rigid supports (case 1) is 24.6mm, it can be noted that the consideration of flexible supports (case 2) causes three different reductions: in differential settlement (36%), in the reaction of support (16%) and in the negative (tension in upper fiber) bending moment (81%). Those results were expected since the soil-structure interaction tends to cause more uniform settlements and to transfer loads from central to the external columns. However, when considering the elastic modulus of the concrete in infinite time together with the flexible support (case 6), it can be seen that the reductions in differential settlement (11%), in support reaction (6%), and in bending moment (30%) are much less when compared to case 2. Thus, it can be stated that creep tends to reduce the effects of the soil-structure interaction, since it reduces the elastic modulus of the structure, making it more flexible.

It is interesting to note that the redistribution of internal forces is quite significant and that these forces are used for the design of the beams. Therefore, not considering the interaction can lead to results on the unsafe side.

5.3 Influence of the coefficient of consolidation

To evaluate the influence of settlement rate on internal forces, analyzes were performed by varying the value of the soil coefficient of consolidation, $c_v$, and comparing them to the cases of elastic and rigid supports. Evaluating the results for the central support (Figure 7), it can be noted that at $t=0$, all analyzes indicate the same response as that of rigid supports, as expected. Over time, there is a smooth transition so that by the end of the analysis all results coincide with that of the elastic supports. It can be concluded that the $c_v$ of the soil does not influence neither the initial nor final response, but rather the way in which the transition between these two states occurs. Furthermore, it seems that the cases of rigid and elastic supports act as an envelope for all other cases.

Figure 7
Result for the central support, over time, for different values of $c_v$.

5.4 Influence of concrete hardening

Analyzing the cases considering concrete hardening on cure (Figure 8), it can be observed that when $c_v$ is high ($c_v=\mathrm{1,7}\bullet {10}^7$ and $c_v=\mathrm{1,7}\bullet {10}^6{m}^2/day$), consolidation occurs in such a short time interval that hardening does not make a difference, the result remaining equal to those of the elastic supports. For the other cases, it can be observed that the smaller the coefficient of consolidation $c_v$, the longer the differential settlement take to occur. Hence, the structure will be more rigid, causing both a reduction in the final differential settlement, and an increase in the internal forces redistribution. Therefore, it can be concluded that the rate of consolidation affects the differential settlements and the internal forces. Furthermore, it is interesting to note that the cases of elastic and rigid supports cease to function as an envelope for possible results.

Figure 8
Result for the central support, over time, for different values of cv, considering concrete hardening.

5.5 Influence of cracking with tension-stiffening

The reinforcement adopted for the analysis of the cracked structure was obtained by designing the cross-section with the usual safety factors (e.g., from Brazilian standard NBR 6118:2014). As the beam is subjected to positive and negative bending moments, it was divided into three sections with two different cross-sections as indicated in Table 4.

Table 4
Reinforcement of the beam (Section 1 along the span and Section 2 at the center).

Bending Moment	Reinforcement (cm2)
	Section 1	Section 2
Negative	0.00	14.73
Positive	8.59	8.59

The 6 cases shown in Figure 6 are also presented in Figure 9, this time considering cracking and tension-stiffening. Comparing both figures, it is notable that there is a reduction of approximately 50% in the difference in support reaction and bending moment for cases 1 and 2, leaving the two analyzes with close results. This reduction in the redistribution of forces is due to the increased flexibility of the elements when considering concrete cracking.

Figure 9
Results for the central support, over time, for the 5 cases considering tension-stiffening

Analyzing cases 2 and 6, it is observed that both the difference in support reaction and in bending moment were reduced by approximately 70%. As such cases represent the initial and the final behaviour of the structure, and the creep of the concrete defines the transition between them, it can be concluded that in the cracked structure creep has a smaller effect on the redistribution of internal forces. This stems from the fact that, due to cracking, the area of concrete that is mechanically contributing is smaller.

The results of cases 2 and 5 are compared in Figure 10 for the central support of the non-cracked (nc) and cracked (c) beams. Note that, without cracking, the difference in the support reaction varies from 15.8% at the initial time to 10.8% at the final time. Considering cracking, the variation goes from 7.6% to 2.8% over time. The same behaviour is observed for the bending moment. Such reductions indicate that there is less force redistribution in the structure. Furthermore, it is observed that both creep and cracking have a tendency to reduce the redistribution of forces.

Figure 10
Comparison of results for the central support of cracked and non-cracked structure.

Observing the influence of the $c_v$ value on the results for the central support (Figure 11), it can be noticed that the behaviour of the cracked beam is similar to that of the non-cracked one, as the values tend to remain the same as the results of the elastic support when time goes to infinity, with the exception of rigid support. In this case, such a difference is due to the redistribution of forces between the concrete and the reinforcement caused by creep, which ends up inducing a change in the internal forces. However, the magnitude of the redistribution of forces in the cracked beam is much lower than in the non-cracked one.

Figure 11
Result for the central support of the cracked structure, over time, for different values of $c_v$.

Including the concrete hardening effect (Figure 12), it is possible to observe that this phenomenon influences the response of the structure for cases where the $c_v$ is smaller, similar to the observation made from Figure 8. However, it is noted that the contribution of such an effect is smaller in the cracked structure, where the variation of the support reaction and of the bending moment were 0.2% and 5.9%, respectively. As for the non-cracked structure, the values were 4.1% and 25.4%, in the same order.

Figure 12
Result for the central support, over time, for different values of $c_v$ considering both concrete hardening and cracking.

6 CONCLUSIONS

The incorporation of hardening on cure, shrinkage, creep, cracking and tension-stiffening of the concrete, and the consolidation of the soil in the soil-structure interaction allowed an assessment of these phenomena not found in previous publications. The modelling proposed to consider all those phenomena in a coupled analysis also contains some novelties.

The results obtained in the present study enabled the following conclusions:

• When considering the creep of concrete, the structure tends to become more flexible, reducing the uniformity of settlements and the redistribution of internal forces. The results observed for bending moment and support reactions are found between the results obtained in the analyses with elastic support and with rigid support.

• The value of the coefficient of consolidation influences the rate of transition between the behaviour with rigid supports and the behaviour with elastic supports. For high values of the coefficient of consolidation, settlement occurs so quickly that the effects of concrete hardening are not noticeable. However, for smaller values, there is a gain in the structure's stiffness and, consequently, a change in differential settlements and internal forces.

• Cracking has a relevant contribution in reducing the stiffness of the structure, significantly reducing the redistribution of forces. Creep effects tend to be attenuated by cracking, as the concrete area that is mechanically stressed is reduced.

• Considering both cracking and creep, the reduction in the stiffness of the structure can be considerable, so that settlements and internal forces can be considerably different from those obtained in the more usual elastic soil-structure analysis.

From the above results, it can be said that the consideration of both creep and cracking of the concrete significantly reduces the stiffness of the structure, attenuating the effects of soil-structure interaction. Thus, it is recommended that such phenomena be considered in analyses, for more realistic results.

Acknowledgments

The senior author thanks the Brazilian Air Force for the leave for his doctoral studies. The studies presented herein were developed in the Civil Engineering Department of COPPE (Graduate School of Engineering), UFRJ (Federal University of Rio de Janeiro). CAPES - Coordenação de Aperfeiçoamento de Pessoal de Nivel Superior, Brazil, partly financed this research under Finance Code 001.

References

• ABNT NBR 6118 (2014) Design of concrete structures - procedure. Associação Brasileira de Normas Técnicas, Rio de Janeiro, RJ.
• Arapakou, A.E., and Papadopoulos, V.P., (2020). Effects of soil simulation on the interaction analyses of framed structures under 2-D and 3-D conditions. Geotechnical and Geological Engineering 38: 5389-5407.
• Bazant, Z.P., (1971). Numerically stable algorithm with increasing time steps for integral-type aging creep. In Proceedings of the 1st International Conference for Structural Mechanics in Reactor Technology. IASMiRT, Berlin, DE, pp. 119-128.
• Bazant, Z.P., and Jirásek, M., (2018). Creep and hygrothermal effects in concrete, Springer (Netherlands).
• Bazant, Z.P., and Wu, S.T., (1973). Dirichlet series creep function for aging concrete. Journal of the Engineering Mechanics Division ASCE 99: 367-387.
• Beeby, A.W., Scott, R.H., Jones, A.E.K., (2005). Revised code provisions for long-term deflection calculations. Proceedings of the Institution of Civil Engineers – Structures and Buildings 158: 71-75.
• Brown, P.T., and Yu, S.K.R., (1986). Load sequence and structure-foundation interaction. Journal of Structural Engineering ASCE 112: 481-488.
• Comite Euro-International du Beton (1993). CEB-FIP Model Code 1990: Design Code.
• Cosenza, E., and Greco, C., (1990). Comparison and optimization of different methods of evaluation of displacements in cracked reinforced concrete beams. Materials and Structures 23: 196-203.
• Dalili S, M., Alkarni, A., Noorzaei, and Paknahad, M., (2011). Numerical simulation of soil-structure interaction in framed and shear-wall structures. Interaction and Multiscale Mechanics 4: 17-34.
• Dalili S., M., Huat, B.B.K, Jaafar, M.S., and Alkami, A., (2014). Soil-Framed Structure Interaction Analysis – A New Interface Element. Latin American Journal of Solids and Structures 12:226-249
• Espion, B., (1993). Benchmark Examples for Creep and Shrinkage Analysis Computer Programs. In Proceedings of th 5th Intertational RILEM Symposium on Creep and Shrinkage of Concrete. RILEM. Barcelona, ES, pp. 877-888.
• Espion, B., and Halleux, P., (1990). Long-term deflection of reinforced concrete beams: Reconsideration of their variability. ACI Structural Journal 87: 232-236.
• Fajman, P., and Sejnoha, J., (2007). A simplified approach to time-dependent subsoil-structure interaction. Computers and Structures 85: 1514-1523.
• Fédération Internationale du Béton (2013). Fib Model Code for Concrete Structures 2010.
• Ghali, A., (1986). A unified approach for serviceability design of prestressed and non prestressed reinforced concrete structures. PCI Journal 31: 118-137.
• Gonçalves, J.C., Santa Maria, P.E.L., and Danziger, F.A.B., (2016). Study of beams and frames subject to creep and shrinkage to support analyses of soil-structure interaction in buildings (in Portuguese). In Proceedings of the XVIII Congresso Brasileiro de Mecânica dos Solos e Engenharia Geotécnica. ABMS, Belo Horizonte, Brazil.
• Goschy, B., (1978). Soil-Foundation-Structure interaction. Journal of the Structural Division ASCE 104: 749-761.
• Gusmão, A.D., (1990). Study of the soil-structure interaction and its influence on building settlements. MSc dissertation (in Portuguese), Universidade Federal do Rio de Janeiro, Brazil.
• Hora, M.S., (2014). Nonlinear interaction analysis of infilled frame-foundation beam-homogeneous soil system. Coupled System Mechanics 3: 267-289.
• Jardine, R.J., Potts, D.M., Fourie, A.B., and Burland, J.B., (1986). Studies of the influence of non-linear stress-strain characteristics in soil-structure interaction. Géotechnique 36: 377-396.
• Kawakami, M., and Ghali, A., (1996). Cracking, ultimate strength and deformations of prestressed concrete sections of general shape. PCI Journal 41: 114-122.
• Lee, I.K., and Brown, P.T., (1972). Structure-Foundation interaction analysis. Journal of Structural Division ASCE 98: 2413-2431.
• Lopes, F.R., and Gusmão, A.D., (1991), On the influence of soil-structure interaction in the distribution loads settlements. In Proceedings of the 10th European Conference on Soil Mechanics and Foundation Engineering. Firenze, Italy, pp. 475-478.
• Majid, K.I., and Cunnell, M.D., (1976). A theoretical and experimental investigation into soil-structure interaction. Géotechnique 26: 331-350.
• Masih, R., (1993). Structural stiffness influence on soil consolidation. Journal of Geotechnical Engineering ASCE 119: 168-172.
• Noorzaei, J., Viladkar, M.N. and Godbole, P.N., (1995). Elasto-plastic analysis for soil-structures interaction in framed structures. Computers & Structures 55: 797-807.
• Rosa, L.M.P., Danziger, B.R. and Carvalho, E.M.L., (2018). Soil-structure interaction analysis considering concrete creep and shrinkage. IBRACON Structures and Materials Journal 11: 564-585.
• Santa Maria, P.E.L., Santa Maria, F.C.M., and Santos, A.B., (1999), Analysis of continuous beams resting on viscoelastic supports and its applications to problems of soil-structure interaction. Solos e Rochas 22: 179-194.
• Taylor, R.L., Filippou, F.C., Saritas, A., and Auricchio, F., (2003). A mixed finite element method for beam and frame problems. Computational Mechanics 31: 192-203.
• Torres, L.L., (2001). Modelo numérico y verificación experimental del comportamiento en servicio de estructuras de hormigón. DSc thesis. Universitat Politécnica de Catalunya, Spain.
• Viladkar, M.N., Ranjan, G., and Sharma, R.P., (1993). Soil-structure interaction in the time domain. Computers and Structures 46: 429-442.

Edited by

Editor: Pablo Andrés Muñoz Rojas

Publication Dates

• Publication in this collection
  22 July 2022
• Date of issue
  2022

History

• Received
  12 Apr 2022
• Reviewed
  18 Apr 2022
• Accepted
  30 June 2022