Next Article in Journal
Advances in the Fabrication of Antimicrobial Hydrogels for Biomedical Applications
Next Article in Special Issue
Upscaling Cement Paste Microstructure to Obtain the Fracture, Shear, and Elastic Concrete Mechanical LDPM Parameters
Previous Article in Journal
Textural, Structural and Biological Evaluation of Hydroxyapatite Doped with Zinc at Low Concentrations
Previous Article in Special Issue
Mesoscale Characterization of Fracture Properties of Steel Fiber-Reinforced Concrete Using a Lattice–Particle Model
 
 
Font Type:
Arial Georgia Verdana
Font Size:
Aa Aa Aa
Line Spacing:
Column Width:
Background:
Article

Lattice Modeling of Early-Age Behavior of Structural Concrete

1
Department of Civil and Environmental Engineering, University of California, Davis, CA 95616, USA
2
School of Civil Engineering, University of Castilla-La Mancha, 13071 Ciudad Real, Spain
*
Author to whom correspondence should be addressed.
Materials 2017, 10(3), 231; https://doi.org/10.3390/ma10030231
Submission received: 27 December 2016 / Revised: 12 February 2017 / Accepted: 18 February 2017 / Published: 25 February 2017
(This article belongs to the Special Issue Numerical Analysis of Concrete using Discrete Elements)

Abstract

:
The susceptibility of structural concrete to early-age cracking depends on material composition, methods of processing, structural boundary conditions, and a variety of environmental factors. Computational modeling offers a means for identifying primary factors and strategies for reducing cracking potential. Herein, lattice models are shown to be adept at simulating the thermal-hygral-mechanical phenomena that influence early-age cracking. In particular, this paper presents a lattice-based approach that utilizes a model of cementitious materials hydration to control the development of concrete properties, including stiffness, strength, and creep resistance. The approach is validated and used to simulate early-age cracking in concrete bridge decks. Structural configuration plays a key role in determining the magnitude and distribution of stresses caused by volume instabilities of the concrete material. Under restrained conditions, both thermal and hygral effects are found to be primary contributors to cracking potential.

1. Introduction

Early-age cracking is often a root cause of premature loss of safety and serviceability of concrete structures, including reinforced concrete bridge decks and slabs. Such cracking may occur due to plastic settlement and shrinkage, prior to concrete setting, or as a consequence of volumetric changes associated with thermal, hygral, or chemical phenomena that occur after setting [1,2]. The factors that affect cracking potential can be categorized in different ways, such as the grouping presented in Table 1. These factors can be viewed as design parameters, since they are controlled (at least to some degree) during the process of design, construction, and curing of the structural components.
Many of the parameters in Table 1 influence autogenous and drying shrinkages, which are prime contributors to early-age cracking of concrete bridge decks and overlays. Control of these shrinkages, possibly through the use of shrinkage reducing admixtures (SRA), has proven to be effective in reducing cracking in such structures [3,4,5,6]. Nonetheless, other factors, acting alone or in combination, may strongly influence the probability of cracking. In particular, thermal effects due to heat of hydration and heat exchange with the environment may be important [7,8,9,10,11].
Knowledge of early-age cracking has been gained primarily through laboratory testing, often using small-scale specimens, and through the monitoring of field performance. Laboratory testing allows for the examination of multiple specimens under controlled conditions, but it only approximates the circumstances of concrete bridge decks. In particular, laboratory testing might not adequately account for the restraint mechanisms and volumetric instabilities associated with actual structural and environmental boundary conditions. On the other hand, data from the monitoring of field performance are scarce and typically provide an incomplete picture of the reasons for good or poor performance. Furthermore, the consequences of variations in material composition, construction processes, and environmental conditions are difficult to anticipate without appropriate models.
This research involves: (1) the development of a methodology for simulating the early-age behavior of structural concrete; and (2) application of the methodology to assessing early-age cracking of concrete bridge decks. The methodology relies on a simple form of lattice model [12], which has been extended to represent the salient thermal, hygral, and mechanical processes that affect the early-age behavior. A primary goal of most lattice analyses of concrete has been the simulation of fracture at the meso-scale, at which the matrix, a coarse fraction of the aggregates, and the matrix-aggregate interface are explicitly modeled [13,14,15,16,17,18,19]. Such lattice models provide a quantitive linking of meso-scale processes and structural behavior. Lattice models are also effective in performing multi-field analyses with emphases on the durability mechanics of concrete materials and structures [20,21,22,23]. The lattice-based analyses presented herein differ in that concrete is regarded as a homogeneous continuum. There is no inherent disadvantage relative to approaches based on finite element technology, which may serve the research needs equally well [24,25]. Property development of the aging concrete is based on a modeling of cementitious materials hydration. The model resolves the spatial fields of temperature, relative humidity, and displacement. From these bases, realistic connections can be made between the design parameters (including materials composition, construction practice, and environmental conditions) and properties that strongly affect susceptibility to cracking, such as strength, stiffness, and creep potential. The simple, effective representation of tensile fracture is one attribute of this modeling approach.
Validation exercises are conducted using measurements taken from small-scale laboratory tests and an instrumented bridge deck. The validated model is used to assess the various contributors to the early-age cracking potential of concrete bridge decks. Along with autogenous deformations and hygral gradients caused by the drying of exposed surfaces, thermal effects may be a prime contributor to deck cracking. The importance of restraint conditions on cracking potential is examined.

2. Modeling Framework

2.1. Program Structure

The analyses of concrete early-age behavior are based on lattice models, which are composed of randomly positioned nodal points that are interconnected by lattice elements (Figure 1). The three relevant fields (i.e., temperature, relative humidity, and displacement) are represented at the nodal points of the lattice. With respect to temperature and relative humidity analyses, the lattice elements may be viewed as conduits that transport heat or moisture between each i j pair of nodes [21]. With respect to the structural analyses, the lattice elements are based on the rigid-body-spring concept [12,26] and are akin to conventional frame elements.
The program is driven by a modeling of cementitious materials hydration, along with the influence of environmental boundary conditions that may vary with time. As described later, the hydration model also serves to determine mechanical properties, including tensile strength, elastic modulus, and creep behavior. After each time step of the thermal-hygral-structural analyses, the incremental changes in temperature and humidity are calculated. Corresponding thermal and hygral strain increments are determined and then introduced into the structural analysis.

2.2. Domain Discretization

The lattice topology is defined by the Delaunay tessellation of a set of randomly placed nodes (Figure 1); the Voronoi tessellation of the same set of nodes defines element properties [21]. The discretization process is semi-automated and involves the sequential definition of Voronoi vertices, edges, and faces, followed by random filling of the domain interior with nodal points [27]. Since the vertices, edges, and faces of the body are defined by sets of four, three, and two nodes, respectively, the discretization of complex geometries is cumbersome relative to use of the finite element method.
The use of irregularly spaced nodal points, connected by 1D elements, typically introduces spurious heterogeneity into the modeling effort. Such irregular lattices can be rendered homogenous through special techniques [28] or, as done herein, by scaling the element properties according to the Voronoi tessellation of the domain [12]. As demonstrated in previous studies, this form of lattice model is also adept at representing heterogeneity, either explicitly or by using probabilistic methods [29,30]. For modeling structural concrete, however, large computational demands are incurred by the necessity of finer, three dimensional discretizations. Alternatively, lattice simulations of meso-scale heterogeneity can be used to support macro-scale analyses of homogeneous continua [8].

2.3. Cementitious Materials Hydration

Desirable features of the hydration model include the ability to simulate the hydration of blended cements, accounting for the chemical composition of the blend constituents. Herein, we use a versatile model of cementitious materials hydration [31,32]. Model development and calibration, done by Riding et al. [31], was based on a dataset of 204 concrete mixtures, which were tested using semi-adiabatic calorimetry. The cementitious materials considered by the model include portland cement, slag, fly ash and/or silica fume. The model has been validated through comparisons with a separate set of 57 semi-adiabatic calorimetry tests.
The influence of temperature on the rate of hydration is represented by a maturity relation
t e = 0 t exp E a R ( 1 T c 1 T r ) Δ t
in which t e (h) is the equivalent age for a material hydrating at reference temperature T r (K); R is the universal gas constant (8.314 J/mol/K); T c is the temperature of the concrete (K); and E a is the apparent activation energy (J/mol), which depends on composition and proportioning of the cementitious materials [33].
The degree of cementitious materials hydration is assumed to be
α ( t ) = H ( t ) H u
where H ( t ) is the cumulative amount of heat produced by the hydration reaction (J/g); and H u is the total heat available for reaction (J/g), which depends on the composition and proportioning of the cementitious materials.
H u = H cem p cem + 461 p slag + 1800 p FA CaO p FA + 330 p SF
in which each term represents the heat of hydration of component i multiplied by its weight fraction, p i . The heat of hydration of the cement (J/g) is
H c e m = 500 p C 3 S + 260 p C 2 S + 866 p C 3 A + 420 p C 4 AF + 624 p SO 3 + 1186 p free Ca + 850 p MgO
The relationship between degree of hydration and equivalent age is
α ( t e ) = α u exp τ t e β
in which β and τ are hydration parameters; and α u is the ultimate degree of cementitious materials reaction. The rate of heat release with respect to time is calculated by combining the information presented in Equations (1), (2) and (5).
Q ( t ) = H u c τ t e β β t e α u exp τ t e β exp E a R 1 T r + 1 T c
where c is the cementitious materials content (g/m 3 ). The rate of heat release depends on β, τ, and α u , which themselves depend on the composition and proportioning of the cementitious materials. These dependencies have been established through multi-variate regression analyses based on concrete mixtures with a wide range of compositions [31].

2.4. Primary Analysis Modules

2.4.1. Thermal Analysis

The governing differential equation for heat conduction is
· ( λ T ) + ρ c p Q c p = ρ c p T t
in which T is temperature (K), Q is the rate of heat production (J/(kg · s)), λ is thermal conductivity (W/(m · K)), ρ is mass density (kg/m 3 ), and c p is specific heat capacity (J/(kg · K)). The boundary conditions for Equation (7) involve either prescribed temperatures or heat flux across the boundary. For bridge deck systems, heat exchange with the environment is of primary importance. Herein, the following types of heat exchange are modeled:
  • Convection—Convective heat exchange across exposed surfaces depends on the difference between the solid surface temperature T s and that of the surrounding ambient medium T a
    q c o n v = Λ T ( T s T a )
    where Λ T is the coefficient of convective heat transfer, which depends on wind speed and other factors [34].
  • Solar radiation—The amount of solar radiation reaching the concrete surface depends on several factors including the structure’s location, surface orientation, altitude, atmospheric conditions, time of the day, and day of the year. Incoming heat due to solar radiation is
    q s u n = γ a b s q i n c
    in which γ a b s is the solar absorptivity of the concrete, which is influenced by the color and texture of the concrete surface, and q i n c is the incident solar radiation acting on a horizontal surface (W/m 2 ). The latter of the two values can be obtained from recorded weather data for the time period and region of interest.
  • Thermal radiation - Heat loss to the surroundings due to grey-body radiation is calculated using
    q s k y = σ ϵ ( T s K 4 T s k y 4 )
    where σ is the Stefan-Boltzmann constant ( 5.669 × 10 2 W/(m 2 · K 4 )), ϵ is the emissivity of the concrete (=0.9, in this study), and T s K is the temperature of the concrete surface (K). T s k y is the sky temperature, which depends on the sky emissivity, the dew point temperature, and the cloud conditions [34].
Heat exchange also occurs due to evaporation and condensation, but those mechanisms are considered to be of secondary importance for ordinary concrete structures. Some other relevant aspects of the thermal analyses are as follows.
  • The heat capacity of the cement paste is estimated using an approach given by Bentz [35], in which heat capacity is a function of degree of reaction of the cement. The heat capacity of the concrete is then determined from the heat capacities of the cement paste and aggregates, according to the mass fractions of each using an ordinary rule of mixtures.
  • Thermal conductivity of the concrete is estimated by taking the average of the Hashin-Shtrikman bounds for a two-phase composite formed of paste and aggregates [35].
For the lattice modeling of heat conduction, the material continuum is represented by a collection of 1D elements connected on a semi-random set of nodal points (Figure 1). Each edge of the Delaunay tessellation, connecting adjacent nodes i and j, acts as a lineal conduit element. A semi-discrete form of Equation (7) serves as the basis for this lattice modeling of heat conduction.
M T ˙ + K T = f
where the capacity, M , and conductivity, K , matrices are assembled from their respective elemental contributions:
M e = h i j A i j 6 d 2 1 1 2 ,
K e = D A i j h i j 1 1 1 1
in which A i j is the area of the Voronoi facet associated with nodes i and j, h i j is the distance between the nodes, and D is thermal diffusivity (= λ / ρ c p ). In the element capacity matrix, d = 1, 2, and 3 for 1D, 2D, and 3D networks, respectively. For each time step, Equation (11) is solved using the Crank-Nicholson method [36] in conjunction with a fixed-point algorithm [37] to achieve convergence. That is, Equation (11) is solved and corrected (for the dependencies of Q on T and c p on α) until the L 2 norm of the difference of nodal temperature values between iterations is satisfactorily small.
Thermal strains in the elastic material arise from changes in temperature at the lattice nodes:
Δ ε T = β T Δ T
where β T is the coefficient of thermal expansion (CTE) and Δ T is the difference in nodal temperature over the time step. Equation (14) uses the average of Δ T calculated at the nodes of a given element; Δ ε T is then introduced in the axial direction of the corresponding element.
The value of β T depends on several factors, including moisture content and the type of aggregate. The dependence on moisture content is particularly strong, such that β T is significantly larger during the setting process and for several hours thereafter [38,39,40]. To account for this behavior, the coefficient of thermal expansion is modeled as indicated in Table 2, in which α T h is the CTE of hardened concrete and χ represents the magnification of CTE prior to final setting. The relationship between β T , degree of reaction α, and the setting parameters ( α 0 i , α 0 , and φ , as described by Krauß and Rostásy [41]) is shown in Figure 2.

2.4.2. Hygral Analysis

The moisture field within the concrete varies with time due to both autogenous and external drying effects. Using relative humidity h as the field variable [42,43], the governing differential equation can be expressed as
h t = div ( D h ( h ) grad h ) + h s t
where D h is hygral diffusivity and h s is relative humidity associated with self-desiccation. Herein, h s is assumed to be a function of degree of hydration [44]
h s = ( h s u 1 ) α α u s + 1
where h s u is the asymptotic value of relative humidity as α / α u approaches unity under sealed conditions. Parameter s controls the rate of humidity drop and is thought to depend on w / c [44].
Various expressions have been used to represent the dependence of hygral diffusivity on moisture content of the material. Herein, the expression of di Luzio and Cusatis is employed [45].
D h ( h , T ) = ψ ( T ) D 1 1 + D 1 D 0 1 ( 1 h ) n 1
where
ψ ( T ) = exp E a R 1 T 0 1 T ,
which represents the influence of temperature on diffusivity. D 1 (m 2 /h) is the diffusivity in the saturated state (h = 1), D 0 (m 2 /h) is the diffusivity of the fully dried material (h = 0), and T 0 is the reference temperature (K). Parameter n controls the position of the diffusivity drop as the concrete dries from a saturated state (Figure 3). The temperature dependence of D h is deemed to be important, since concrete bridge decks may experience significant solar heating during the day and cooling at night. For a portland cement concrete, the dependence of diffusivity on water-to-cement ratio, w / c , can be represented by [45]
D 0 = D 0 ˜ w c 3 and D 1 = D 1 ˜ w c 2.5
where D 0 ˜ (m 2 /h) and D 1 ˜ (m 2 /h) are material parameters.
Convective transport of moisture across exposed surfaces is governed by
q h = Λ h ( h h a )
in which Λ h is the hygral convection coefficient and h a is the ambient relative humidity.
The lattice representation of the hygral field, and its solution procedure, are analogous to that of the thermal field (i.e., the lattice elements can be viewed as conduits that transfer moisture between nodes of differing relative humidity). This approach was developed by Sadouki and van Mier [46] and later adapted to this form of random lattice [21,23]. Alternatively, the hygral field can be represented by a dual lattice network defined by the Voronoi edges of the tessellated volume [47,48,49]. The main advantage of this alternative approach comes in the post-cracking state: elements are aligned with the crack surfaces, which more naturally accounts for crack-assisted transport and diffusion from the crack faces into the bulk material.
At this stage of model development, it is assumed that the diffusion process is uncoupled from the mechanical behavior of the material and any damage incurred during mechanical or hygral loading. Shrinkage strains in the elastic material arise from changes in relative humidity at the lattice nodes. A constant hygral shrinkage coefficient β h is used in this one-way coupling of the hygral and mechanical analyses
Δ ε h = β h Δ h
where Δ h is the difference in nodal relative humidity over the time step. Equation (21) uses the average of Δ h calculated at the nodes of a given element and Δ ε h is then introduced in the axial direction of the corresponding element.
Coupling of the thermal and hygral fields is not considered herein, but has been included in other studies [50,51]. Such coupling is particularly relevant when simulating moderate to high-temperature loading of concrete [52,53].

2.4.3. Structural Analysis

A structural element is defined by two neighboring nodes, i and j, and their common Voronoi facet (Figure 1). The element stiffness relations are based on a zero-size spring set, located at the area centroid (point C) of the Voronoi facet, and connected to the element nodes via rigid-arm constraints. The spring set consists of three axial springs, oriented normal and tangential to the facet, and three rotational springs about the same local (n-s-t) axes. The stiffness coefficients of the axial springs are:
k n = ξ k s = ξ k t = E A i j h i j ,
where A i j is the Voronoi facet area; h i j is the distance between nodes i and j; factor ξ relates the normal and shear spring stiffnesses, which can be adjusted to simulate macroscopic Poisson ratio of the material [18]. In most cases of random discretization, there will be boundary layer effects on the mechanical properties [54]. For ξ = 1, which is used herein, the lattice is elastically homogeneous under uniform modes of straining [55,56], although Poisson ratio ν = 0. Correct representation of the Poisson effect, both macroscopically and in an element-local sense, has recently been accomplished using an iterative procedure [57,58]. The stiffness coefficients of the rotational springs are:
k ϕ n = E J p h i j , k ϕ s = E I 11 h i j , k ϕ t = E I 22 h i j ,
where J p is the polar second moment of the facet area, and I 11 and I 22 are the two principal second moments of the facet area. Directions s and t are aligned with the facet principal axes. The spring constants appear on the diagonal of the material matrix, D , given by
D = ( 1 ω ) diag k n , k s , k t , k ϕ n , k ϕ s , k ϕ t ,
where ω is a scalar damage parameter used to model material fracture. Prior to fracture initiation, ω = 0. The element stiffness matrix (with respect to element local coordinates) is
K e = B T D B .
in which B relates the generalized spring displacements and nodal displacements [55]. After transforming K e to global coordinates, the direct stiffness approach is used to assemble element stiffness matrices and internal force contributions into the structural equation set
K Δ u = ( R + F )
in which Δ u is the increment in generalized nodal displacements; K is the system stiffness matrix assembled from the elemental components; R and F are the external and internal nodal force vectors, respectively.
After each equilibrium iteration within a load increment, the spring set forces of each element i j are updated based on Δ u and the spring stiffnesses (Equations (22) and (23)), which are related through the rigid-body constraints. The tensile stress resultant acting in a given element is
σ R = F n 2 + F s 2 + F t 2 A i j P
where F n , F s , and F t are the forces in the normal and two tangential springs, respectively; and A i j P is the projection of the Voronoi facet area of element i j on a plane perpendicular to the force resultant. A planar representation of this concept is shown in Figure 4a, where θ R indicates the inclination of the force resultant relative to the element axis; and s i j is the length of the Voronoi segment common to nodes i and j, such that A i j P = b s i j cos θ R where b is the element thickness. The dimension of h i j cos θ R represents the crack band width [59].
Tensile fracture is simulated using an event-by-event approach [19,60], in which the element with the largest stress ratio κ = max σ R / f ( α , w ) 1 undergoes fracture. At degree of hydration α, residual tensile strength f ( α , w ) is defined by a softening relation in terms of crack opening w (Figure 4b). This approach provides energy conserving, grid-insensitive simulations of tensile fracture, as confirmed by previous studies [48,55,56].
The spring set properties evolve with degree of cementitious materials hydration. Creep and stiffness development are based on solidification and microprestress theory, as described in the next section. The development of tensile strength is represented by
f ( α ) = f ( α u ) α α 0 α u α 0 ζ for α α 0
where f ( α u ) is the tensile strength at the ultimate degree of hydration; α 0 is the degree of hydration associated with setting; and coefficient ζ =1 when tensile strength is being modeled [41,61]. The softening relation evolves as a function of degree of hydration, as shown in Figure 4b where r is the fraction of ultimate tensile strength acting at hydration degree, α. By holding w c constant, fracture energy (i.e., the area under the softening curve) grows with degree of hydration, whereas the characteristic length decreases [62,63,64]. Kinetics of the reaction are also influential and can be modeled [25,65], but are neglected herein for the sake of simplicity and in light of the uncertainties in the model inputs. In Equation (28), strength development depends phenomenologically on the degree of hydration. Alternatively, property development can be based on explicit representations of the cementitious materials and their hydration products [66,67,68]. The increased physical bases of such microstructural models are attractive, especially when simulating early-age behavior.
If f is based on a short duration tensile test, cracking occurs at σ R / f < 1 under sustained loading due the occurrence of tertiary creep. From uniaxial tension tests, Domone [69] found that creep rupture occurs when the sustained stress exceeds 0.75 f. Altoubat and Lange [70] measured the fraction of tensile strength acting at the time of cracking to be about 0.60 to 0.64 for split cylinder tests, and 0.75 to 0.80 for uniaxial tension tests. van Breugel found that, under longer duration loading caused by self-induced stress, early-age cracking occurs at about 0.75 f [71]. From beam tests, Wittmann et al. [72] determined the stress at failure to be less than the tensile strength (measured by a short duration test) by a factor of 0.6 to 0.8. Emborg [7,62] reports that for low loading rates, such as those associated with the cooling of mass concrete, the reduction factor applied to tensile strength is about 0.7. A reduction factor of 0.7 is assumed for the parametric analyses presented in Section 4.

2.5. Stiffness and Creep Representation

The stiffness and creep of concrete are modeled using solidification and microprestress theory [25,73,74]. In effect, the axial (n-component) springs of the rigid-body-spring network have been replaced with the series construct shown in Figure 5. The construct consists of: (1) a spring representing instantaneous elastic deformation, ε i ; (2) a chain of m solidifying Kelvin units representing viscoelastic creep, ε v ; (3) an aging dashpot with viscosity that is dependent on microprestress, representing viscous creep, ε f ; and 4) a unit representing the combination of thermal strain, ε T , and hygral strain, ε h . The shear (s- and t-component) springs of the rigid-body-spring network are replaced by the same construct, but without the thermal and hygral unit. Superposition of the individual strain units, and differentiating their values with respect to time, produces the total strain rate:
ε ˙ = ε ˙ i + ε ˙ v + ε ˙ f + ε ˙ h + ε ˙ T
The relationship between instantaneous strain rate and uniaxial stress rate, σ ˙ , can be expressed by
ε ˙ i = q 1 σ ˙
where q 1 (MPa 1 ) is an age independent parameter.
According to solidification theory, viscoelastic strains arise from creep of the solidified hydration products (which are mainly calcium-silicate-hydrate gels) whose properties do not change with age. Viscoelastic microstrain has the following formulation:
γ = 0 t Φ ( t r ( t ) t r ( τ ) ) σ ˙ d τ
in which t r ( t ) is a reduced time function that accounts for the effects of temperature and humidity on fine-scale creep processes; Φ ( t t 0 ) represents the non-aging micro-compliance function of cement gel, which can be defined as
Φ ( t t 0 ) = q 2 ln 1 + ( t t 0 ) 0.1
where t t 0 is the load duration, and q 2 (MPa 1 ) is an adjustable model parameter [73].
Aging of the cement paste matrix is attributed to accumulation of the solidified hydration products. The viscoelastic micro- and macro-strain rates are related as follows:
ε ˙ v ( t ) = 1 v ( α ) γ ˙
where v ( α ) is an aging function, representing the volume fraction of solidified material at degree of hydration, α. Following di Luzio and Cusatis [25], the aging function is expressed as:
v ( α ) = α α u n α
where n α is an adjustable model parameter.
The determination of viscoelastic microstrains, γ, using the history integral formulation of Equation (31) is computationally expensive, even for modest numbers of degrees of freedom. The problem is therefore recast as a rate formulation based on a Dirichlet series expansion of Φ ( t t 0 ) . Details regarding this approximation, the determination of viscoelastic compliance, and other aspects of the numerical implementation of the viscoelastic creep model are given elsewhere [25,73].
According to the microprestress theory [73], hindered adsorbed water produces a disjoining pressure within the micropores, which induces a microprestress in the solidified material. For modeling purposes, microprestress S is defined as the average stress acting on bonds that bridge the micropores. Viscous creep results from shear slip caused by the breakage of overstressed bonds. Under varying temperature and humidity conditions, the viscous creep flow is
ε ˙ f ( t ) = ψ ( t ) σ ( t ) η S
where the viscosity η is a decreasing function of the microprestress, which is given by
1 η S = q 4 κ 0 S
where q 4 (MPa 1 ) and κ 0 are adjustable parameters [25]. The evolution of microprestress is governed by
S ˙ + ψ S κ 0 S 2 = κ 1 T ˙ ln h + T h ˙ h
where κ 1 is a model parameter. In the above formulation, ψ and ψ S are reduced time coefficients representing the effects of temperature and relative humidity on the creep processes and microprestress evolution, respectively. The coefficient applied to the creep processes is
ψ = ( 0.1 + 0.9 h 2 ) exp E v R 1 T 0 1 T t
where E v is the activation energy ( E v / R 5000 K). The coefficient applied to microprestress evolution, ψ S , takes the same form except that activation energy E S replaces E v ( E S / R 3000 K). Microprestress can be determined from a central difference approximation of Equation (37). Details of this determination, and other aspects of the viscous creep modeling, are provided elsewhere [25,73].

3. Validation Exercises

3.1. Stiffness and Basic Creep Development

The implementation of solidification and microprestress theory [25,74] is verified for the case of basic creep. Although numerous creep test data are available for comparison, cases of early age loading are of primary interest herein. The results of Laplante [75] have been commonly used as a benchmark for creep behavior beginning with loading at 24 h [25,65].
Cylindrical specimens were loaded at 1, 3, 7, and 28 days after casting [75]. The specimens were sealed to prevent moisture exchange between the concrete and environment. The specimens were loaded at 30% of the concrete strength at each respective age.
For the creep simulations, the cylindrical specimens are discretized as shown in Figure 6. The cylinder dimensions of 100 mm by 200 mm differ from those of the experiment [75], but that is irrelevant for the simulation of basic creep of a material with uniform properties. A compressive load is applied uniformly over the top and bottom faces of the cylinder without lateral restraint. Each element of the lattice is a spring-dashpot construct that also accounts for hygral and thermal effects, as shown in Figure 5. At t = 0, the freshly cast concrete is assumed to be fully saturated (h = 1). The radial surface of the concrete exchanges heat with the environment (at T = 21 C, as recorded in the experiment) via a convective boundary condition. The cement composition is required as model input, as indicated in Section 2.3. For this simulation, a typical ASTM Type I cement was assumed [76]: C 3 S = 55%; C 2 S = 22%; C 3 A = 10%; C 4 AF = 8%; and Blaine fineness = 365 m 2 /kg. A single set of creep parameter values has been used for the series of simulations: q 1 = 11.5 × 10 6 /MPa; q 2 = 28 × 10 6 /MPa; q 4 = 0.18 × 10 6 /MPa; and n α = 1.6.
Comparisons of the model and experimental results are presented in Figure 7. The simulation results are reasonably accurate except for the case of 3-day loading, for which the simulated creep values are significantly less than the measured values. Creep associated with microprestress is activated by changes in temperature and internal humidity of the concrete [25]. Although moisture exchange with the environment is prevented by sealing, free water is consumed by the hydration process and such drying activates creep due to microprestress (i.e., the model simulates the effects of self-desiccation on creep).
Due to the lack of friction on the planes of load application, all material points are in a state of uniaxial compression. Neglecting the effects of non-uniform temperature within the specimen due to heat exchange with the environment, the cylinder model of Figure 6 should provide precisely the same results as a single element aligned with the loading direction. Such a comparison is given in Figure 8. The small difference between the 1-day loading curves is attributable to the accounting of heat of hydration in the cylinder model, whereas the single element was held at the ambient temperature of 21 C. Heat of hydration increases the degree of reaction, albeit only slightly for the case of a small cylindrical test specimen. At later ages of loading, the differences in degree of hydration are insignificant, such that the creep deformation curves are indistinguishable. These results highlight the ability of the irregular lattice to simulate uniaxial creep, using solidification and microprestress theory, without mesh bias.

3.2. Strength Development

For compressive strength development, Equation (28) becomes
f c ( α ) = f c ( α u ) α α 0 α u α 0 ζ for α α 0
where f c ( α u ) is the asymptotic limit of compressive strength with degree of hydration. Herein, exponent ζ is assigned the typical value of 1.5 [41]. Figure 9 compares estimated strength values with those of the Laplante study, which was considered in the previous section. Strengthening commences for α α 0 . A typical range of α 0 provides reasonably good estimates of strength development, which supports the validation of the hydration model.

3.3. Autogenous and Drying Shrinkage Tests

The series of tests of Yang et al. [77] are simulated. The test program involved the autogenous and drying shrinkage testing of blended cement concrete at three water-to-cementitious materials ratios: 0.25, 0.35, and 0.45. The binder consisted of ordinary portland cement and blast furnace slag at a 50:50 proportion. The blast furnace slag was finely ground, having a Blaine fineness of 600 m 2 /kg, which likely accentuated autogenous shrinkage. For each w / c value, both autogenous shrinkage and drying shrinkage tests were conducted using 100 × 100 × 400 mm 3 prisms. When drying was considered, it was initiated at 1, 3, or 7 days. The 3-day and 7-day tests are simulated hereafter. Results from the 1-day drying experiments exhibit a qualitatively different trend that could not be captured using the same set of model parameters.
For the lattice modeling of the tests, symmetry conditions are exploited such that 1/8 of the prism volume is discretized (Figure 10). Heat and moisture exchanges occur across the exterior surfaces, but are prevented across the planes of symmetry. In addition, nodes along the planes of symmetry do not move normal to those planes of symmetry. Nodes were placed at ends of the prism, along the longitudinal axis, to serve as gage points for calculating longitudinal strain. Nodal density has been reduced within the core of the specimen where temperature and hygral gradients are relatively small.
Chemical composition of the portland cement fraction was not reported, so that values typical for Type I OPC have been used for the hydration model calculations. Humidity change due to self-desiccation and concrete diffusivity are represented by Equations (16) through (19) with the parameter values given in Table 3. The large h s u values are due, in part, to the high volumes of finely ground slag in the concrete mixtures.
The simulated shrinkage strains are compared with the experimental values in Figure 11. The overall trends are captured well particularly for the autogenous shrinkage strain curves, which are relatively easy to fit. For the cases of drying, however, the experimental results exhibit sustained high rates of shrinkage in the first few days. That behavior can be simulated more precisely by increasing the diffusivity and/or hygral coefficient of the concrete, but then the ultimate shrinkage values are too large. The model is limited to the service range of loading (i.e., linear behavior, including creep) and does not account for the nonlinear effects of microcracking and creep damage under sustained high loading. Ultimately, those effects need to be incorporated into the analyses. The fitting exercise presented here should be understood in that context.

3.4. Analysis of Concrete Bridge Decks

As part of a Caltrans sponsored project [78], the early-age behavior of concrete bridge decks was studied through their instrumentation and data acquisition for a period of time after concrete casting. Sensors positioned within the freshly cast concrete monitored temperature, internal relative humidity, and strain. Along with these sensor readings, ambient temperature, relative humidity, and wind speed were recorded on-site by a weather station. Observations of cracking behavior were also made. Herein, simulated temperatures are compared with such field measurements for the Markham Ravine Bridge on Route CA-55 near Lincoln, California. In addition, concrete cylinders and prisms were cast on-site for measuring strength and drying shrinkage, respectively. Simulations of these specimens, along with the deck temperatures, serve to calibrate the lattice model for the parametric study of deck cracking potential in the next chapter.

3.4.1. Model Definition

Discretization of the bridge deck/girder system is shown in Figure 12. Nodes are placed at the thermocouple locations, as indicated in the figure. Concrete forming the soffit and girder stems is assumed to be mature. Simulation of cementitious materials hydration is limited to the freshly cast deck. Inputs to the hydration model are indicated in Figure 13. Whereas the computational framework is three-dimensional, the simulations presented herein are for planar characterizations of the structure. The effects of reinforcing steel on the thermal and mechanical components of the analyses have not been included. Whereas the restraint provided by reinforcing steel reduces the free shrinkage of concrete, the modeling of this effect in a planar analysis framework is not straightforward. Reinforcement elements can be readily incorporated into this form of lattice model [79] and we plan to do so after extending the analysis to three dimensions. Boundary conditions for the thermal, hygral, and mechanical analyses are described in the following subsections.
Thermal analysis inputs: The initial temperature of the simulated, cast concrete was set to its measured value at the time of casting (19.5 C). The initial temperature of the simulated, mature concrete was set to the value measured at mid-height within one of the girder stems at the time of casting (14.9 C).
Figure 14 shows typical measured solar radiation intensities for clear and overcast days for the geographic proximity and time of the year when casting occurred [80]. Along with an additional profile, representing broken cloud cover, these intensity plots are used as templates for the daily solar radiation input to the model.
One point of interest is the use of thermal/curing sheets on the cast surface. A combination of burlap and plastic membrane (Transgard 4000) was used to cover the bridge deck for the time interval 0.25 t 9 days, where t = 0 is the time of concrete casting. Placement of this curing medium modifies heat exchange with the environment.
q c o n v = η 1 q c o n v q s u n = η 2 q s u n for 0.25 t 9 days q s k y = η 3 q s k y
in which η 1 , η 2 , and η 3 are reduction factors associated with each respective form of heat flux. According to experimental measurements [81], placement of a fabric/plastic membrane can reduce convective heat transfer by roughly a factor of 5 for the case of wind velocity v = 0 m/s, which results in η 1 = 0.2. The product sheet for Transgard 4000 indicates a light reflectance of 0.85. This is roughly twice the reflectance value of ordinary portland cement concrete, which ranges from about 0.34 to 0.48 [82], such that η 2 = 0.23. For lack of information, the same reduction factor was assumed for the sky radiation boundary condition (i.e., η 3 = 0.23).
The coefficient of convective heat transfer, Λ T , was a function of measured wind speed [34]. Based on a series model [24], a modified value for Λ T was utilized to account for the presence of 1/2 (12.7 mm) plywood formwork along the lower face of the fresh concrete. The thermal conductivity of plywood was assumed to be 0.12 W/(m · K).
Measurements of the CTE of the hardened deck concrete provided β T h = 8.6 × 10 6 / C. Variation in CTE at early age was modeled using the expressions given in Table 2, in which χ = 3.5 based on measurements of Kada et al. [38] on various high-performance concrete mixtures. Similar large variations were observed for cement pastes and mortars in other studies [83].
Hygral analysis inputs: The ultimate value of relative humidity, associated with self-desiccation, was assumed to be h s u = 0.9 with s = 3 governing the rate of self-desiccation according to Equation (16). The parameter values expressing hygral diffusivity, according to Equations (17) and (19), are: D 0 ˜ = 0.017 mm 2 /h; D 1 ˜ = 9 mm 2 /h, and n = 5.
The concrete surfaces bounded by plywood formwork are assumed to be sealed (i.e., no moisture exchange occurs between the concrete and plywood formwork.) While the curing sheeting is in place (0.25 d t 9 d), along with periodic wetting, no moisture is exchanged between the upper surface of the concrete deck and the environment. Otherwise, the exposed surface exchanges moisture with the environment according to Equation (20), using Λ h = 0.25 mm/h. Soon after finishing of the cast concrete, a light-colored membrane curing compound was sprayed on the deck surface. Based on laboratory measurements [84], this membrane is assumed to reduce convective transport across the surface by a factor of 2 (i.e., q h = 0.5 q h ). The hygral shrinkage coefficient β h = 0.0030 was calibrated with experimental measurements given in Section 3.4.2. A zero-flux condition was enforced across the plane of symmetry (z = 0).
Mechanical analysis inputs: The parameters for the solidification and microprestress modeling of concrete behavior are presented in Table 4. The creep parameter values were set according to the B4 model [85], based on the actual mixture composition. We have made no effort to distinguish between tensile creep, of interest herein, and compressive creep (even though tensile creep may be much larger than compressive creep for the same levels of applied stress [86]).
All stresses originate from thermal and hygral strains. Externally applied loads (e.g., due to construction activities) were not considered at this stage of the analysis. Elements spanning the cast deck and the mature concrete substrate were given properties of the substrate concrete. This is a reasonable approximation given the small element size and our relative lack of interest in behavior at that bi-material interface. A zero-displacement condition was enforced across the plane of symmetry (z = 0), along with a single support along that plane to prevent rigid-body motion in the vertical direction.
The degree of hydration at concrete setting, α 0 , is required for the strength development model expressed by Equation (28). Within the test program associated with the Markham Ravine Bridge deck concrete, penetration resistance was measured according to ASTM C403 [87] to determine initial and final setting times, t i s and t f s , respectively. As per standard, mortar was extracted from the deck concrete and subjected to penetration testing. The results are shown in Figure 15, in which σ p is the measured penetration resistance and σ f s is the penetration resistance associated with final set. The prescribed resistance levels associated with initial and final set are 3.5 MPa (500 psi) and 27.6 MPa (4000 psi), respectively. The initial and final set times ( t i s = 4.9 h and t f s = 6.33 h) are determined by the intersections of these respective stress levels with the resistance curve.
To determine the degrees of hydration at initial and final set, a volume of mortar was simulated using the measured initial temperature of the mortar (18.89 C) and ambient temperature (21.1 C). For the thermal calculations, the coarse aggregates were removed within the hydration routine and in the calculation of specific heat capacity. The simulated degree of hydration is plotted as a function of time in Figure 15. The degrees of hydration associated with t i s and t f s are α 0 i = 0.085 and α 0 = 0.13, respectively. Some studies suggest the development of the mechanical threshold (e.g., tensile resistance) of hardened concrete begins rather close to the time of initial setting [88]. Other sources use the point of final set to define the mechanical threshold [89,90]. The α 0 value is used for strength development in the bridge deck simulations that follow. The influences of w / c and aggregate content on the mechanical threshold have been investigated via computational modeling [68].

3.4.2. Simulation Results

Deck temperature histories: The simulated temperature histories are compared with the field measurements in Figure 16, for the case of thermocouples TC2 and TC5 located within the mid-deck region and above the girder stem, respectively. The recorded ambient temperature history is also plotted in the figures. Comparisons with the entire set of thermocouple readings are presented in Figure 17. Several comments can be made.
  • The influence of environmental factors is evident from the oscillatory behavior of the temperature history recorded by each thermocouple. After the first day, locations closer to the surface exhibit larger temperature swings, whereas deeper locations are less affected by environmental changes. This meets expectations.
  • Peak temperatures occur at about 10 h after concrete casting (Figure 18). The lower temperatures over the supporting girders are due to conduction of heat toward the cooler substrate concrete. Conversely, the insulative properties of the plywood formwork give rise to higher temperatures between the supporting girders. These temperatures are significantly higher than the ambient temperature. The differences between the ambient and measured deck temperatures largely diminish over the first two days after casting, in contrast to mass concrete applications in which large temperature differences can exist for several weeks. Despite the discrete, irregular discretization of the domain, the iso-contours of temperature do not exhibit artifacts associated with mesh bias.
  • After removal of the curing sheets, not only do diurnal variations in deck temperature increase, but also the temperature gradient through the deck thickness tends to increase. The implications of the larger thermal gradients are discussed later.
Strength development: Concrete cylinder specimens were cast on site and kept local to the bridge deck for 7 days, prior to laboratory storage and measurement of splitting tensile strength. Strength development was simulated using Equation (28), with α 0 = 0.20, and the same lattice model adopted for the creep test simulations (Figure 6). Ambient temperature for the modeling exercise was constant and set equal to the average recorded temperature over the curing period. Figure 19a compares simulated strength development with the measured values, where the 28-day splitting tensile strength is used as a normalizing factor.
Shrinkage development: For measuring the drying shrinkage properties of the concrete mixture, prisms were cast on-site and transferred to the laboratory within 24 h, after which they were kept in a moist condition until exposure to a drying environment at t = 7 days. Testing was done according to ASTM C157, except for the 7 days of moist pre-conditioning according to Caltrans recommendations. The measured and simulated shrinkage strains are compared in Figure 19b. Initial expansion of the prisms was simulated by prescribing h = 0.97 as an initial condition prior to moist curing at 1 day. Whereas the sources of expansion may be due to other factors [91], this simulation of swelling demonstrates the expected workings of the model. As seen in the previous examples of autogenous and drying shrinkage, the shape of the simulated shrinkage curve does not conform to that of the experimental curve: the simulated rate of shrinkage over the first several days is lower. This could be remedied by increasing hygral diffusivity and/or the shrinkage coefficient, but the ultimate shrinkage strain would then be overestimated. Here, too, it is apparent that effects of microcracking need to be incorporated into the analyses.

4. Parametric Study

Through the preceding examples, basic workings of the thermal, hygral, and mechanical components of the model have been validated. The model is now used to assess the early-age cracking potential of the Markham Ravine Bridge. In particular, the relative importances of thermal and hygral contributions to cracking potential are evaluated for two design factors: (1) the use of curing media; and (2) restraint to deck movement. Cracking of the concrete is not simulated, except for the analysis results presented near the end of Section 4.2. The scope of this exercise is limited to planar analyses of the bridge deck using the lattice model shown in Figure 12. The implications of this simplification of the bridge deck configuration and boundary conditions are discussed later.

4.1. Curing Protocol

As noted above, curing compound and polymer/fabric sheets were applied on the concrete deck as part of the curing protocol. Concrete stresses at the mid-deck location, normalized by the 28-day tensile strength, are plotted in Figure 20 for the case when these curing media are absent. To compare the thermal and hygral stress contributions with their sum, the z-component of the stress tensor has been evaluated at the TC1 and TC3 thermocouple locations. These σ Z values are roughly equal to the principle tensile stress values at the same locations. The hygral and thermal contributions to total stress are isolated by setting β T = 0 and Λ h = 0, respectively.
As expected, exposure to the drying environment produces tension near the deck surface, due to the hygral gradient in the y-direction (Figure 20a). Location TC3 is in compression, as is much of the depth of the deck, in response to the surface tension. The oscillatory nature of the stress curves is due to diurnal variations in ambient relative humidity: relative humidity climbs at night due to decreasing temperature such that moisture intake and swelling occurs near the upper surface of the deck. The large drops in stress magnitudes beginning at about t = 6 and 13 days are due to rain events. Rainfall was not modeled directly, but rather captured through the relative humidity boundary condition. Sustained, large relative humidity values were measured on-site during those rain events. Based on simulated field testing of concrete, Asamoto et al. [92] have found that shrinkage of concrete is significantly reduced by rainfall and that continuous rainy days have such influence, even if the amount of precipitation is small.
Figure 20b plots the stresses associated with thermal effects. During daytime heating, the upper surface of the deck expands relative to the lower surface, producing compression at the TC1 location and tension at the TC3 location. This tendency is reversed during nighttime cooling. For either case, the tensile stresses are significantly smaller than those produced by hygral loading. For both the hygral and thermal loading cases, stress variations over the girder stem (at TC4 and TC6) are qualitatively similar, but slightly smaller in magnitude.
Tensile strength of the concrete at the TC1 location, as determined from the degree of cementitious material reaction (Equation (28)), is indicated in each plot within Figure 20. The evolution of tensile strength, scaled by a reduction factor of 0.7 (as discussed in Section 2.4.3), is also indicated in Figure 20c. The tensile stresses in Days 2, 3, and 4 cross over this curve, which suggests the occurrence of cracking.
While in place, the curing sheets act to significantly reduce both the hygral and thermal stresses, as shown in Figure 21. Upon removal of the sheets, however, sudden drying of the deck surface from the nearly saturated state causes high tensile stresses. For this case, as well, the tensile stress values cross over the scaled tensile strength curve indicating the likelihood of cracking.
Knowing the stress state and tensile strength at each nodal location, we can plot maps of cracking potential (i.e., σ 1 / f , where σ 1 is principle tensile stress and f is the current tensile strength value at that location) over the computational domain. Figure 22 provides such plots for the thermal stress contributor and total stress at t = 11.42 days, which is near the time of maximum σ 1 / f for the TC1 location. To highlight thermal stresses, different scales are used for each plot. The highly stressed region extends about 40 mm (or 1/5 of the deck thickness) from the exposed surface. Similarly steep stress gradients have been found near the surface of concrete pavements [93] and bridge decks [9]. Cracking would tend to relieve stress along the drying surface and, due to the sharp gradient in stress, it appears that cracks would not propagate through the deck thickness. That hypothesis is evaluated using concrete fracture mechanics, as described later.
The strong potential for cracking indicated by Figure 20 through Figure 22 is accentuated by at least three factors: (1) the high 28-day shrinkage measured in the laboratory and used for the analyses (Figure 19b); (2) the low tensile strength of the test prisms used for this comparison; and (3) lack of consideration of membrane-type curing compound applied in the field study. These simulations have been rerun after: (1) scaling the drying shrinkage curve in Figure 19 to meet the 28-day limit of 450 μm/m prescribed by several transportation agencies [94]; and (2) calculating the 28-day tensile strength according to f t = 0.3 ( f c ) 2 / 3 , where f c is the measured 28-day compressive strength [89]. Based on these assumptions, the likelihood for cracking is greatly reduced, as shown in Figure 23a,b for each respective simulation case.
Furthermore, if a curing compound has been applied, the cracking potential is further reduced. For the simulation results presented in Figure 23c, application of the curing compound is assumed to reduce the hygral convective coefficient by 50%, as described in Section 3.4.1. Analyses of the effects of structural configuration, which follow, are based on the parameter combination corresponding to Figure 23b (in which the prism 28-day shrinkage strain is 450 μ m/m and the tensile strength is based on the measured compressive strength).

4.2. Structural Configuration

The structural configuration used for the previous analyses is representative of in-plane conditions away from supports or diaphragms placed along the span length. As a rough approximation of the conditions local to such stiffening elements, the cell volumes have been filled with concrete (Figure 24). The concrete within the cells is given the same properties as those of the soffit and girder stems. Otherwise, the following analyses retain the same set of parameter values used to produce the results in Figure 23b.
Stress values due to the thermal component are plotted in Figure 25a. The amplitudes of the diurnal variations are larger than for the unrestrained cases. Furthermore, the entire deck cross-section experiences tension, in contrast to the previous results where temperature gradients through the thickness produced mainly tension-compression couples. The increase in deck tension over the first several days corresponds to the overall cooling of the deck, as shown in Figure 17.
The deck stresses due to hygral effects are plotted in Figure 25b,c. The stresses due to autogenous shrinkage are significant, in contrast to the unrestrained cases where such stresses were negligible and therefore not plotted in the previous examples. In Figure 23b, for example, there is no upward drift into tension over the first 9 days when the curing sheets are in place. For the case of moisture loss to the environment (Figure 25c), the stress amplitudes are slightly higher than for the case without restraint. Otherwise, the results for the two cases are qualitatively similar.
Total stresses at the mid-span location, which include the thermal, autogenous, and drying shrinkage components, are plotted in Figure 25d. The total stresses are not a simple summation of the component values, since the autogenous and drying shrinkages are coupled through the dependence of hygral diffusivity on internal relative humidity (Equation (17)). Comparing with the case where restraint is absent (Figure 23b), the stress amplitudes plotted in Figure 25d are greater and the entire cross-section experiences tension.
By allowing fracture to occur, according to the procedure outlined in Section 2.4.3, consequences of the restraint conditions become evident (Figure 26). When restraint is absent, the sharp humidity gradient (with secondary contributions from the thermal gradient) produces cracking near the drying surface. This has been observed in other studies, as well [95]. When the restraint is present, both autogenous shrinkage and global temperature change of the deck introduce tension over the deck cross-section. Crack density is higher and some cracks propagate through the deck cross-section, as shown in Figure 26b. Similar cracking behavior has been noted in field studies where concrete abutments provided fixity to the deck ends [9,96]. Transitions from diffuse distributed cracking to localized fracture are simulated simply and effectively, which is desirable for evaluating most concrete durability problems.
The surface cracking due to hygral and thermal gradients is more extensive than shown in Figure 26, which only depicts cracks larger than a prescribed threshold. Furthermore, the openings of such surface cracks are much smaller than those of the through cracks that form when restraint is present. This is seen in the histogram of cracked-element openings given in Figure 26c. As the drying front progresses with time, many of the surface cracks have closed and are therefore not counted in the histogram. The time of day also influences the opening counts. For the 12:00 pm counting shown in Figure 26c, solar heating of the top surface has reduced or closed some of the surface cracks. The scatter in openings for the through cracks, about an average opening of about 0.7 mm, is mainly due to flexure of the deck under thermal/hygral loading. Crack opening is proportional to depth within the deck cross-section.
Whereas these results indicate deck cracking would occur, the restraint condition associated with the internal diaphragm is an extreme case. Three-dimensional analysis of the deck is necessary to more accurately assess the effects of structural configuration on deck stresses. Furthermore, the effects of restraint and thermal straining in the longitudinal direction of the bridge are not captured by these planar analyses, even though such effects might contribute significantly to cracking potential.

5. Conclusions

This research pertains to modeling the early-age behavior of concrete materials with realistic structural and environmental boundary conditions. The discrete modeling approach, presented herein, is a novel assemblage of works of the authors and many other researchers. The relevant field quantities (temperature, relative humidity, and displacement) are represented by lattice models, which are defined by Delaunay/Voronoi tessellations of the computational domain. The simulations begin when the concrete materials are mixed. Based on a modeling of cementitious materials hydration, the simulations provide a spatial description of property (e.g., stiffness, strength) and stress development. This information can be used to assess the potential for early-age cracking. The following conclusions were deduced from this study:
  • This form of discrete model is capable of simulating the multi-field quantities associated with early-age concrete behavior, despite its discontinuous representation of the problem domain. There are no major disadvantages of this discrete approach with respect to continuum approaches, such as the finite element method. Advantages of this form of discrete model include its simplicity and adeptness at simulating the transition from diffuse damage to localized cracking.
  • Stress-to-strength ratio is lacking as a practical measure of cracking potential. Sharp hygral or thermal gradients near exposed surfaces typically cause high stresses, which is indicative of cracking. However, the cracking may only be superficial. Knowledge of the stress conditions through the structural cross-section is also necessary for evaluating the severity and potential consequences of cracking. In this sense, numerical modeling nicely complements knowledge gained by laboratory testing and field observations.
  • Structural configuration plays a key role in determining the magnitude and distribution of stresses caused by volume instabilities of the concrete material. Under largely restrained conditions, both thermal and hygral effects were found to be primary contributors to cracking potential, leading to crack propagation through the depth of the deck. Three-dimensional simulations are needed to assess the influence of longitudinal restraints, thermal flexing of the mature concrete girders, and the effects of reinforcing bars.
  • Realistic simulations of the early-age behavior of structural concrete require a wealth of information regarding the material constituents, production/curing processes, and structural and environmental boundary conditions. As shown by the parametric studies conducted herein, cracking potential is sensitive to input quantities that are typically not well defined, especially in field applications. Ultimately, the assessment of cracking potential needs to be cast in a probabilistic framework that accounts for uncertainties in the various inputs to the modeling effort.

Acknowledgments

The authors gratefully acknowledge support provided by Caltrans through technical agreement #65A0532 with the University of California, Davis.

Author Contributions

The objectives and calculation methods form bases for the ongoing PhD dissertation research of Yaming Pan and Armando Prado As a visiting scholar at UC Davis, Rocio Porras Soriano actively contributed to experiments that helped guide the simulation work. Omar Hafez made key contributions to the routines for simulating cementitious materials hydration and the thermal-hygral-mechanical couplings. The research was coordinated by John Bolander All authors participated in writing the manuscript.

Conflicts of Interest

The authors declare no conflicts of interest.

Abbreviations

The following symbols are used in this manuscript:
SymbolDefinition
h i j Euclidean distance between nodes i and j [m]
ttime [s, h, d]
A i j area of Voronoi facet common to nodes i and j [m 2 ]
Thermal analysis
a / c aggregate to cementitious materials ratio [kg/kg]
ccementitious materials content [kg/m 3 ]
c p specific heat capacity [J/(kg · K)]
q c o n v , q c o n v heat flux due to convection and its modified value [W/m 2 ]
q s u n , q s u n heat flux due to solar radiation and its modified value [W/m 2 ]
q s k y , q s k y heat flux due to radiation and its modified value [W/m 2 ]
t e concrete equivalent age [h]
w / c water-to-cementitious materials ratio [kg/kg]
Dthermal diffusivity [m 2 /s]
E a apparent activation energy [J/mol]
Hcumulative amount of heat produced by hydration [J/kg]
H c e m heat of hydration of cement [J/kg]
H u total heat available for reaction [J/kg]
Qrate of heat production by cementitious materials hydration [J/(kg · s)]
Runiversal gas constant [J/(mol · K)]
Ttemperature [ C]
T a , T c , T r , T s ambient, concrete, reference, and surface temperatures [ C]
T s k y sky temperature [ C]
α, α u degree of cementitious materials hydration and its ultimate value
β, τhydration model parameters
β T , β T h coefficient of thermal expansion and its long-term value [1/ C]
γ a b s solar absorptivity of concrete
ϵemissivity of concrete
η 1 , η 2 , η 3 heat flux reduction factors
λthermal conductivity [W/(m · K)]
ρmass density of concrete [kg/m 3 ]
σStefan-Boltzmann constant [W/(m 2 · K 4 )]
χmagnifying factor for coefficient of thermal expansion at early-ages
Λ T coefficient of convective heat transfer [W/(m 2 · K)]
M e , M elemental and system capacity matrices
K e , K elemental and system conductivity matrices
Hygral analysis
hrelative humidity [-]
h a ambient relative humidity [-]
h s , h s u relative humidity associated with self-desiccation and its ultimate value [-]
nparameter relating hygral diffusivity and relative humidity
q h , q h moisture flux due to convection and its modified value [m/s]
sparameter relating humidity decrease and self-desiccation
D 0 , D 0 ˜ diffusivity values associated with fully dried material [m 2 /s]
D 1 , D 1 ˜ diffusivity values associated with saturated material [m 2 /s]
D h hygral diffusivity [m 2 /s]
β h hygral shrinkage coefficient
Λ h hygral convection coefficient [m/s]
Mechanical analysis
belement thickness [m]
f c , f c u compressive strength and its asymptotic limit [MPa]
ftensile strength [MPa]
k n , k s , k t stiffness coefficients of the uniaxial springs [N/m]
k ϕ x , k ϕ y , k ϕ z stiffness coefficients of the rotational springs [N · m]
n α parameter relating rate of solidification to degree of hydration
q 1 , q 2 , q 4 instantaneous elastic, viscoelastic, and viscous strain parameters [1/MPa]
t 0 age of loading [h]
t e equivalent age [h]
t r reduced time representing the thermal-hygral effects on creep [h]
t i s , t f s times of initial and final concrete setting [h]
wcrack opening [m]
w c traction-free crack opening [m]
A i j P projected area of element facet [m 2 ]
Eelastic modulus of concrete [MPa]
E v , E S activation energies for creep and microprestress processes [J/mol]
F n , F s , F t spring-set forces in the normal and two tangential directions [N]
I 11 , I 22 , J p principal planar and polar second moments of facet area [m 4 ]
D material constitutive matrix
K e , K elemental and global stiffness matrices
α 0 i , α 0 degrees of hydration at initial and final setting
φ coefficient governing concrete setting behavior
γviscoelastic microstrain
ε h , ε T hygral and thermal strains
ε i , ε v , ε f instantaneous elastic, viscoelastic, and viscous flow strains
ζcoefficient governing material property development
ηeffective microscopic viscosity [MPa/s]
κstress-to-strength ratio
κ 0 microprestress model parameter [1/(MPa · d)]
κ 1 microprestress model parameter [MPa/K]
νPoisson’s ratio
θ R inclination of tensile force resultant relative to the element axis
ϑ 1 , ϑ 2 parameters defining the break-point in the tension softening relation
σaxial stress [MPa]
σ R stress resultant [MPa]
σ p , σ f s penetration resistance and its value at final setting [MPa]
ψ , ψ S reduced time coefficients for creep and microprestress processes
ξfactor relating normal and shear spring stiffnesses
ωdegree of damage
Φ non-aging micro-compliance function [1/MPa]

References

  1. Bentz, D.P. A review of early-age properties of cement-based materials. Cem. Concr. Res. 2008, 38, 196–204. [Google Scholar] [CrossRef]
  2. Comite Euro-International du Beton (CEB). Durable concrete structures—Design Guide; Thomas Telford Ltd.: London, UK, 1992. [Google Scholar]
  3. Banthia, N.; Gupta, R. Plastic shrinkage cracking in cementitious repairs and overlays. Mater. Struct. 2009, 42, 567–579. [Google Scholar] [CrossRef]
  4. Bentz, D.P.; Geiker, M.R.; Hansen, K.K. Shrinkage-reducing admixtures and early age desiccation in cement pastes and mortars. Cem. Concr. Res. 2001, 31, 1075–1085. [Google Scholar] [CrossRef]
  5. Maggenti, R.; Knapp, C.; Fereira, S. Controlling shrinkage cracking. Concr. Int. 2013, 35, 36–41. [Google Scholar]
  6. Weiss, W.J.; Yang, W.; Shah, S.P. Shrinkage cracking of restrained concrete slabs. J. Eng. Mech. 1998, 124, 765–774. [Google Scholar] [CrossRef]
  7. Emborg, M.; Bernander, S. Assessment of risk of thermal cracking in hardening concrete. J. Struct. Eng. 1994, 120, 2893–2912. [Google Scholar] [CrossRef]
  8. Šavija, B.; Schlangen, E. Use of phase change materials (PCMs) to mitigate early age thermal cracking in concrete: Theoretical considerations. Constr. Build. Mater. 2016, 126, 332–344. [Google Scholar] [CrossRef]
  9. Subramaniam, K.V.L. Identification of early-age cracking in concrete bridge decks. J. Perform. Constr. Facil. 2016, 30. [Google Scholar] [CrossRef]
  10. Subramaniam, K.V.L.; Kunin, J.; Curtis, R.; Streeter, D. Influence of early temperature rise and movements and stress development in concrete bridge decks. J. Bridge Eng. 2010, 15, 108–116. [Google Scholar] [CrossRef]
  11. Zanotti, C.; Meda, A.; Plizzari, G.; Cangiano, S. Evaluation of the risk of cracking in thin concrete walls due to hydration heat. In Concrete under Severe Conditions: Environment and Loading (CONSEC’10); CRC Press: Boca Raton, FL, USA, 2010; pp. 257–264. [Google Scholar]
  12. Bolander, J.E.; Saito, S. Fracture analyses using spring networks with random geometry. Eng. Fract. Mech. 1998, 61, 569–591. [Google Scholar] [CrossRef]
  13. Benkemoun, N.; Poullain, P.H.; Al Khazraji, H.; Choinska, M.; Khelidj, A. Meso-scale investigation of failure in the tensile splitting test: Size effect and fracture energy analysis. Eng. Fract. Mech. 2016, 168, 242–259. [Google Scholar] [CrossRef]
  14. Cusatis, G.; Pelessone, D.; Mencarelli, A. Lattice discrete particle model (LDPM) for failure behavior of concrete. I: Theory. Cem. Concr. Compos. 2011, 33, 881–890. [Google Scholar] [CrossRef]
  15. Eliáš, J.; Vořechovský, M.; Skoček, J.; Bažant, Z.P. Stochastic discrete meso-scale simulations of concrete fracture: Comparison to experimental data. Eng. Fracture Mech. 2015, 135, 1–16. [Google Scholar] [CrossRef]
  16. Grassl, P.; Jirásek, M. Meso-scale approach to modelling the fracture process zone of concrete subjected to uniaxial tension. Int. J. Solids Struct. 2010, 47, 957–968. [Google Scholar] [CrossRef]
  17. Kozicki, J.; Tejchman, J. Modelling of fracture process in concrete using a novel lattice model. Granul. Matter 2008, 10, 377–388. [Google Scholar] [CrossRef]
  18. Nagai, K.; Sato, Y.; Ueda, T. Mesoscopic simulation of failure of mortar and concrete by 3D RBSM. J. Adv. Concr. Technol. 2005, 3, 385–402. [Google Scholar] [CrossRef]
  19. Schlangen, E.; van Mier, J.G.M. Experimental and numerical analyses of micromechanisms of fracture of cement-based composites. Cem. Concr. Compos. 1992, 14, 105–118. [Google Scholar] [CrossRef]
  20. Alnaggar, M.; Cusatis, G.; di Luzio, G. Lattice discrete particle modeling (LDPM) of alkali silica reaction (ASR) deterioration of concrete structures. Cem. Concr. Compos. 2013, 41, 45–59. [Google Scholar] [CrossRef]
  21. Bolander, J.E.; Berton, S. Simulation of shrinkage induced cracking in cement composite overlays. Cem. Concr. Compos. 2004, 26, 861–871. [Google Scholar] [CrossRef]
  22. Liu, L.; Ye, G.; Schlangen, E.; Chen, H.; Qian, Z.; Sun, W.; van Breugel, K. Modeling of the internal damage of saturated cement paste due to ice crystallization pressure during freezing. Cem. Concr. Compos. 2011, 33, 562–571. [Google Scholar] [CrossRef]
  23. Šavija, B.; Pacheco, J.; Schlangen, E. Lattice modeling of chloride diffusion in sound and cracked concrete. Cem. Concr. Compos. 2013, 42, 30–40. [Google Scholar] [CrossRef]
  24. Faria, R.; Anzenha, M.; Figueiras, J.A. Modeling of concrete at early ages: Application to an externally restrained slab. Cem. Concr. Compos. 2006, 28, 572–585. [Google Scholar] [CrossRef]
  25. Di Luzio, G.; Cusatis, G. Solidification-microprestress-microplane (SMM) theory for concrete at early age: Theory, validation and application. Int. J. Solids Struct. 2013, 50, 957–975. [Google Scholar] [CrossRef]
  26. Kawai, T. New discrete models and their application to seismic response analysis of structures. Nucl. Eng. Des. 1978, 48, 207–229. [Google Scholar] [CrossRef]
  27. Asahina, D.; Bolander, J.E. Voronoi-based discretizations for fracture analysis of particulate materials. Powder Technol. 2011, 213, 92–99. [Google Scholar] [CrossRef]
  28. Schlangen, E.; Garboczi, E.J. New method for simulating fracture using an elastically uniform random geometry lattice. Int. J. Eng. Sci. 1996, 34, 1131–1144. [Google Scholar] [CrossRef]
  29. Asahina, D.; Landis, E.N.; Bolander, J.E. Modeling of phase interfaces during pre-critical crack growth in concrete. Cem. Concr. Compos. 2011, 33, 966–977. [Google Scholar] [CrossRef]
  30. Kang, J.; Kim, K.; Lim, Y.M.; Bolander, J.E. Modeling of fiber-reinforced cement composites: Discrete representation of fiber pullout. Int. J. Solids Struct. 2014, 51, 1970–1979. [Google Scholar] [CrossRef]
  31. Riding, K.A.; Folliard, K.J.; Juenger, M.C.G.; Schindler, A.K. Modeling hydration of cementitious systems. ACI Mater. J. 2012, 109, 225–234. [Google Scholar]
  32. Schindler, A.K.; Folliard, K.J. Heat of hydration models for cementitious materials. ACI Mater. J. 2005, 102, 24–33. [Google Scholar]
  33. Poole, J.L. Modeling Temperature Sensitivity and Heat Evolution of Concrete. Ph.D. Thesis, University of Texas at Austin, Austin, TX, USA, 2007. [Google Scholar]
  34. Bentz, D.P. A Computer Model to Predict the Surface Temperature and Time of Wetness of Concrete Pavements and Bridge Decks; Technical Report NISTIR 6551; National Institute of Standards and Technology: Gaithersburg, MD, USA, 2000. [Google Scholar]
  35. Bentz, D.P. Transient plane source measurements of the thermal properties of hydrating cement pastes. Mater. Struct. 2007, 40, 1073–1080. [Google Scholar] [CrossRef]
  36. Lewis, R.W.; Morgan, K.; Thomas, H.R.; Seetharamu, K.N. The Finite Element Method in Heat Transfer Analaysis; John Wiley & Sons, Inc.: Hoboken, NJ, USA, 1996. [Google Scholar]
  37. Burden, R.L.; Faires, J.D. Numerical Analysis; PWS Publishers: Boston, MA, USA, 1985; 676p. [Google Scholar]
  38. Kada, H.; Lachemi, M.; Petrov, N.; Bonneau, O.; Aïtcin, P. Determination of the coefficient of thermal expansion of high performance concrete from initial setting. Mater. Struct. 2002, 35, 35–41. [Google Scholar] [CrossRef]
  39. Maruyama, I.; Teramoto, A.; Igarashi, G. Strain and thermal expansion coefficients of various cement pastes during hydration at early ages. Mater. Struct. 2014, 47, 27–37. [Google Scholar] [CrossRef]
  40. Sellevold, E.J.; Bjøntegaard, ∅. Coefficient of thermal expansion of cement paste and concrete: Mechanisms of moisture interaction. Mater. Struct. 2006, 39, 809–815. [Google Scholar] [CrossRef]
  41. Krauß, M.; Rostásy, F.S. Determination of initial degree of hydration by means of ultrasonic measurements. In Control of Cracking in Early Age Concrete, Proceedings of the International Workshop on Control of Cracking in Early Age Concrete; Mihashi, H., Wittmann, F.H., Eds.; AA Balkema Publishers: Sendai, Japan, 2002; pp. 19–28. [Google Scholar]
  42. Bažant, Z.P.; Najjar, L.J. Drying of concrete as a non-linear diffusion problem. Cem. Concr. Res. 1971, 1, 461–473. [Google Scholar] [CrossRef]
  43. Jonasson, J.-E. Moisture Fixation and Moisture Transfer in Concrete; Concrete Structures, SMiRT 8; IASMiR: Brussels, Belgium, 1985; pp. 235–242. [Google Scholar]
  44. Oh, B.H.; Cha, S.W. Nonlinear analysis of temperature and moisture distributions in early-age concrete structures based on degree of hydration. ACI Mater. J. 2003, 100, 361–370. [Google Scholar]
  45. Di Luzio, G.; Cusatis, G. Hygro-thermo-chemical modeling of high performance concrete, I: Theory. Cem. Concr. Compos. 2009, 31, 301–308. [Google Scholar] [CrossRef]
  46. Sadouki, H.; van Mier, J.G.M. Simulation of hygral crack growth in concrete repair system. Mater. Struct. 1997, 203, 518–526. [Google Scholar] [CrossRef]
  47. Grassl, P. A lattice approach to model flow in cracked concrete. Cem. Concr. Compos. 2009, 31, 454–460. [Google Scholar] [CrossRef]
  48. Grassl, P.; Bolander, J.E. Three-dimensional network model for coupling of fracture and mass transport in quasi-brittle geomaterials. Materials 2016, 9, 782. [Google Scholar] [CrossRef]
  49. Nakamura, H.; Srisoros, W.; Yahiro, R.; Kunieda, M. Time-dependent structural analysis considering mass transfer to evaluate deterioration process of RC structures. J. Adv. Concr. Technol. 2006, 4, 147–158. [Google Scholar] [CrossRef]
  50. Černy, R.; Rovnaniková, P. Transport Processes in Concrete; CRC Press: Boca Raton, FL, USA, 2002; 560p. [Google Scholar]
  51. Klemczak, B.; Knoppik-Wróbel, A. Reinforced concrete tank walls and bridge abutments: Early-age behaviour, analytic approaches and numerical models. Eng. Struct. 2015, 84, 233–251. [Google Scholar] [CrossRef]
  52. Bary, B.; Durand, S.; Ranc, G.; Carpentier, O. A coupled thermo-hydro-mechanical model for concrete subjected to moderate temperatures. Int. J. Heat Mass Transfer 2008, 51, 2847–2862. [Google Scholar] [CrossRef]
  53. Gawin, D.; Pesavento, F.; Schrefler, B.A. Modeling of hygro-thermal behavior of concrete at high temperature with thermo-chemical and mechanical material degradation. Comput. Methods Appl. Mech. Eng. 2003, 192, 1731–1771. [Google Scholar] [CrossRef]
  54. Eliáš, J. Boundary layer effect on behavior of discrete models. Materials 2017, 10, 157. [Google Scholar] [CrossRef]
  55. Berton, S.; Bolander, J.E. Crack band model for fracture in irregular lattices. Comput. Methods Appl. Mech. Eng. 2006, 195, 7172–7181. [Google Scholar] [CrossRef]
  56. Bolander, J.E.; Sukumar, N. Irregular lattice model for quasistatic crack propagation. Phys. Rev. B 2005, 71, 094106. [Google Scholar] [CrossRef]
  57. Asahina, D.; Aoyagi, K.; Kim, K.; Birkholzer, J.; Bolander, J.E. Elastically-homogeneous lattice models of damage in geomaterials. Comput. Geotech. 2017, 81, 195–206. [Google Scholar] [CrossRef]
  58. Asahina, D.; Ito, K.; Houseworth, J.; Birkholzer, J.; Bolander, J.E. Simulating the Poisson effect in lattice models of elastic continua. Comput. Geotech. 2016, 70, 60–67. [Google Scholar] [CrossRef]
  59. Bažant, Z.P.; Oh, B.H. Crack band theory for fracture of concrete. Mater. Struct. 1983, 16, 155–176. [Google Scholar] [CrossRef]
  60. Rots, J.G.; Bellitti, B.; Invernizzi, S. Robust modeling of RC structures with an event-by-event strategy. Eng. Fract. Mech. 2008, 75, 590–614. [Google Scholar] [CrossRef]
  61. De Schutter, G.; Taerwe, L. Degree of hydration-based description of mechanical properties of early age concrete. Mater. Struct. 1996, 29, 335–344. [Google Scholar] [CrossRef]
  62. Emborg, M. Development of mechanical behaviour at early ages. In Prevention of Thermal Cracking in Concrete at Early Ages; Springenschmid, R., Ed.; RILEM Report 15; E & FN Spon: London, UK, 1998; pp. 76–148. [Google Scholar]
  63. ∅stergaard, L.; Lange, D.; Stang, H. Early-age stress-crack opening relationships for high performance concrete. Cem. Concr. Compos. 2004, 26, 563–572. [Google Scholar] [CrossRef]
  64. De Schutter, G.; Taerwe, L. Fracture energy of concrete at early ages. Mater. Struct. 1997, 30, 67–71. [Google Scholar] [CrossRef]
  65. Cervera, M.; Oliver, J.; Prato, T. Thermo-chemo-mechanical model for concrete. II: Damage and creep. J. Eng. Mech. 1999, 125, 1028–1039. [Google Scholar] [CrossRef]
  66. Bentz, D.P. Three-dimensional computer simulation of cement hydration and microstructure development. J. Am. Ceramic Soc. 1997, 80, 3–21. [Google Scholar] [CrossRef]
  67. Van Breugel, K. Numerical simulation of hydration and microstructural development in hardening cement-based materials: (I) Theory. Cem. Concr. Res. 1995, 25, 319–331. [Google Scholar] [CrossRef]
  68. Torrenti, J.M.; Benboudjema, F. Mechanical threshold of cementitious materials at early age. Mater. Struct. 2005, 38, 299–304. [Google Scholar] [CrossRef]
  69. Domone, P.L. Uniaxial tensile creep and failure of concrete. Mag. Concr. Res. 1974, 26, 144–152. [Google Scholar] [CrossRef]
  70. Altoubat, S.A.; Lange, D.A. Creep, shrinkage, and cracking of restrained concrete at early age. ACI Mater. J. 2001, 4, 323–331. [Google Scholar]
  71. Van Breugel, K.; Lokhorst, S.J. Stress-based crack criterion as a basis for the prevention of through cracks in concrete structures at early-ages. In International RILEM Conference on Early Age Cracking in Cementitious Systems; Kovler, K., Bentur, A., Eds.; RILEM Publications SARL: Paris, France, 2003; pp. 229–236. [Google Scholar]
  72. Wittmann, F.H.; Roelfstra, P.E.; Mihashi, H.; Huang, Y.-J.; Zhang, X.-H.; Nomura, N. Influence of age at loading, water-cement ratio and rate of loading on fracture energy of concrete. Mater. Struct. 1987, 20, 103–110. [Google Scholar] [CrossRef]
  73. Bažant, Z.P.; Cusatis, G.; Cedolin, L. Temperature effect on concrete creep modeled by microprestress-solidification theory. J. Eng. Mech. 2004, 130, 691–699. [Google Scholar] [CrossRef]
  74. Bažant, Z.P.; Prasannan, S. Solidification theory for concrete creep. I: Formulation. J. Eng. Mech. 1989, 115, 1691–1703. [Google Scholar] [CrossRef]
  75. Laplante, P. Mechanical Properties of Hardening Concrete: A Comparative Analysis of Classical and High Strength Concretes. Ph.D. Thesis, Ecole Nationale des Ponts et Chaussées, Paris, French, 1993. [Google Scholar]
  76. Mindess, S.; Young, J.F.; Darwin, D. Concrete; Prentice Hall: Upper Saddle River, NJ, USA, 2003; 644p. [Google Scholar]
  77. Yang, Y.; Sato, R.; Kawai, K. Evaluation of autogenous shrinkage and drying shrinkage based on bound water content of cementitious materials. Proc. JSCE 2001, 53, 193207. [Google Scholar] [CrossRef]
  78. Wiss, Janney, Elstner Associates, Inc. On-Call Structural Concrete Bridge Deck Cracking Investigation Services; Technical Report WJE No. 2009.2643; Wiss, Janney, Elstner Associates, Inc.: Northbrook, IL, USA, 2011. [Google Scholar]
  79. Bolander, J.E.; Le, B. Modeling crack development in reinforced concrete structures under service loading. Constr. Build. Mater. 1999, 13, 23–31. [Google Scholar] [CrossRef]
  80. National Renewable Energy Laboratory, National Solar Radiation Database, U.S. Department of Energy. Available online: http://rredc.nrel.gov/solar/old_data/nsrdb/ (accessed on 21 February 2017).
  81. Lee, Y.; Choi, M.-S.; Yi, S.-T.; Kim, J.-K. Experimental study of the convective heat transfer coefficient of early-age concrete. Cem. Concr. Compos. 2009, 31, 60–71. [Google Scholar] [CrossRef]
  82. Marceau, M.L.; VanGeem, M.G. Solar reflectance values for concrete. Concr. Int. 2008, 30, 52–58. [Google Scholar]
  83. Wyrzykowski, M.; Lura, P. Moisture dependence of thermal expansion in cement-based materials at early ages. Cem.Concr. Res. 2013, 53, 25–35. [Google Scholar] [CrossRef]
  84. Fattuhi, N.I. Curing compounds for fresh or hardened concrete. Build. Environ. 1986, 21, 119–125. [Google Scholar] [CrossRef]
  85. RILEM Technical Committee TC-242-MDC. RILEM draft recommendation: TC-242-MDC multi-decade creep and shrinkage of concrete: Material model and structural analysis. Mater. Struct. 2015, 48, 753–770. [Google Scholar]
  86. Forth, J.P. Predicting the tensile creep of concrete. Cem. Concr. Compos. 2015, 55, 70–80. [Google Scholar] [CrossRef]
  87. ASTM International. ASTM C403-08 Standard Test Method for Time of Setting of Concrete Mixtures by penetration Resistance; Technical Report; ASTM International: West Conshohocken, PA, USA, 2008. [Google Scholar]
  88. Hammer, T.A.; Fosså, K.T.; Bjøntegaard, ∅. Cracking tendency of HSC: Tensile strength and self generated stress in the period of setting and early hardening. Mater. Struct. 2007, 40, 319–324. [Google Scholar] [CrossRef]
  89. Mehta, P.K.; Monteiro, P.J.M. Concrete: Microstructure, Properties, and Materials; Mc Graw Hill: New York, NY, USA, 2014; 675p. [Google Scholar]
  90. Roziere, E.; Cortas, R.; Loukili, A. Tensile behaviour of early age concrete: New methods of investigation. Cem. Concr. Compos. 2015, 55, 153–161. [Google Scholar] [CrossRef]
  91. Baroghel-Bouny, V.; Mounanga, P.; Khelidj, A.; Loukili, A.; Rafai, N. Autogenous deformations of cement pastes, Part II. W/C effects, micro-macro correlations, and threshold values. Cem. Concr. Res. 2006, 36, 123–136. [Google Scholar] [CrossRef]
  92. Asamoto, S.; Ohtsuka, A.; Kuwahara, Y.; Miura, C. Study on effects of solar radiation and rain on shrinkage, shrinkage cracking and creep of concrete. Cem. Concr. Res. 2011, 41, 590–601. [Google Scholar] [CrossRef]
  93. Grasley, Z.C.; Lange, D.A.; D’Ambrosia, M.D. Internal relative humidity and drying stress gradients in concrete. Mater. Struct. 2006, 39, 901–909. [Google Scholar] [CrossRef]
  94. Fu, T.; Deboodt, T.; Ideker, J.H. Development of shrinkage limit specification for high performance concrete used in bridge decks. Cem. Concr. Compos. 2016, 72, 17–26. [Google Scholar] [CrossRef]
  95. De Sa, C.; Benboudjema, F.; Thiery, M.; Sicard, J. Analysis of microcracking induced by differential drying shrinkage. Cem. Concr. Compos. 2008, 30, 947–956. [Google Scholar] [CrossRef]
  96. Darwin, D.; Browning, J.; Lindquist, W.D. Control of cracking in bridge decks: Observations from the field. Cem. Concr. Aggreg. 2004, 26. [Google Scholar] [CrossRef]
Figure 1. Lattice model: (a) domain discretization based on Delaunay and Voronoi tessellations; (b) lattice element i j .
Figure 1. Lattice model: (a) domain discretization based on Delaunay and Voronoi tessellations; (b) lattice element i j .
Materials 10 00231 g001
Figure 2. Assumed relationship between coefficient of thermal expansion and degree of cementitious materials reaction.
Figure 2. Assumed relationship between coefficient of thermal expansion and degree of cementitious materials reaction.
Materials 10 00231 g002
Figure 3. Dependence of hygral diffusivity on relative humidity according to Equation (17) with ψ ( T ) = 1.
Figure 3. Dependence of hygral diffusivity on relative humidity according to Equation (17) with ψ ( T ) = 1.
Materials 10 00231 g003
Figure 4. (a) Determination of tensile stress for planar analyses; and (b) tension softening relation. Parameters ϑ 1 and ϑ 2 define the break point of the softening diagram in terms of tensile strength f ( α u ) and traction-free crack opening w c , respectively.
Figure 4. (a) Determination of tensile stress for planar analyses; and (b) tension softening relation. Parameters ϑ 1 and ϑ 2 define the break point of the softening diagram in terms of tensile strength f ( α u ) and traction-free crack opening w c , respectively.
Materials 10 00231 g004
Figure 5. Series construction of total strain components.
Figure 5. Series construction of total strain components.
Materials 10 00231 g005
Figure 6. Discretization for creep test simulations: (a) lattice representation; and (b) volume rendering of cylindrical specimens.
Figure 6. Discretization for creep test simulations: (a) lattice representation; and (b) volume rendering of cylindrical specimens.
Materials 10 00231 g006
Figure 7. Basic creep curves for loading at different ages.
Figure 7. Basic creep curves for loading at different ages.
Materials 10 00231 g007
Figure 8. Basic creep curves produced from a single element and from a fully discretized cylindrical specimen.
Figure 8. Basic creep curves produced from a single element and from a fully discretized cylindrical specimen.
Materials 10 00231 g008
Figure 9. Strength development as a function of degree of hydration: dependence on setting parameter, α 0 .
Figure 9. Strength development as a function of degree of hydration: dependence on setting parameter, α 0 .
Materials 10 00231 g009
Figure 10. Discretization for shrinkage test simulations: (a) dimensions (in mm); (b) lattice representation; and (c) volume rendering of a symmetric portion of the prism specimens.
Figure 10. Discretization for shrinkage test simulations: (a) dimensions (in mm); (b) lattice representation; and (c) volume rendering of a symmetric portion of the prism specimens.
Materials 10 00231 g010
Figure 11. Autogenous and drying shrinkage simulations.
Figure 11. Autogenous and drying shrinkage simulations.
Materials 10 00231 g011
Figure 12. Lattice representation of symmetric portion of Markham Ravine bridge deck. The locations of thermocouples TC1, TC2, TC3 (midspan) and TC4, TC5, TC6 (over girder stem) are indicated. Dimensions are in meters.
Figure 12. Lattice representation of symmetric portion of Markham Ravine bridge deck. The locations of thermocouples TC1, TC2, TC3 (midspan) and TC4, TC5, TC6 (over girder stem) are indicated. Dimensions are in meters.
Materials 10 00231 g012
Figure 13. Input parameters for hydration model.
Figure 13. Input parameters for hydration model.
Materials 10 00231 g013
Figure 14. Solar radiation intensities for typical clear and overcast days.
Figure 14. Solar radiation intensities for typical clear and overcast days.
Materials 10 00231 g014
Figure 15. Determination of degrees of hydration at initial and final sets using recorded penetration resistance data [78].
Figure 15. Determination of degrees of hydration at initial and final sets using recorded penetration resistance data [78].
Materials 10 00231 g015
Figure 16. Temperature variation within the Markham Ravine Bridge deck: (a) mid-deck at location TC2; (b) solar radiation profiles selected for each day; and (c) above girder stem at location TC5.
Figure 16. Temperature variation within the Markham Ravine Bridge deck: (a) mid-deck at location TC2; (b) solar radiation profiles selected for each day; and (c) above girder stem at location TC5.
Materials 10 00231 g016
Figure 17. Measured and simulated temperatures at thermocouple locations.
Figure 17. Measured and simulated temperatures at thermocouple locations.
Materials 10 00231 g017
Figure 18. Simulated iso-contours of temperature (in C) at t = 10 h.
Figure 18. Simulated iso-contours of temperature (in C) at t = 10 h.
Materials 10 00231 g018
Figure 19. Property development in on-site cast and simulated specimens: (a) splitting tensile strength; and (b) shrinkage strain.
Figure 19. Property development in on-site cast and simulated specimens: (a) splitting tensile strength; and (b) shrinkage strain.
Materials 10 00231 g019
Figure 20. Simulated stresses at mid-deck location when the curing media is absent: (a) hygral stress component (* daily rainfall, in mm, recorded at the Carmichael 0.9 S meteorological station); (b) thermal stress component; and (c) total stress.
Figure 20. Simulated stresses at mid-deck location when the curing media is absent: (a) hygral stress component (* daily rainfall, in mm, recorded at the Carmichael 0.9 S meteorological station); (b) thermal stress component; and (c) total stress.
Materials 10 00231 g020
Figure 21. Simulated stresses at mid-deck location when the curing media is present during 0.25 < t < 9 days: (a) hygral stress component; (b) thermal stress component; and (c) total stress.
Figure 21. Simulated stresses at mid-deck location when the curing media is present during 0.25 < t < 9 days: (a) hygral stress component; (b) thermal stress component; and (c) total stress.
Materials 10 00231 g021
Figure 22. Spatial maps of cracking potential at t = 11.42 days (8 pm on day 11) for the case where the curing sheets were removed at t = 9 days: (a) thermal stress component; and (b) total stress.
Figure 22. Spatial maps of cracking potential at t = 11.42 days (8 pm on day 11) for the case where the curing sheets were removed at t = 9 days: (a) thermal stress component; and (b) total stress.
Materials 10 00231 g022
Figure 23. Simulated stresses at mid-deck location based on 28-day strain limit of 450 μ m/m: (a) the curing sheet is absent; (b) curing sheet is present; and (c) curing sheet and curing compound used.
Figure 23. Simulated stresses at mid-deck location based on 28-day strain limit of 450 μ m/m: (a) the curing sheet is absent; (b) curing sheet is present; and (c) curing sheet and curing compound used.
Materials 10 00231 g023
Figure 24. Discretization of box girder cross-section including cell volumes.
Figure 24. Discretization of box girder cross-section including cell volumes.
Materials 10 00231 g024
Figure 25. Simulated stresses at mid-deck location in the presence of constraint: (a) thermal stress component; (b) hygral stress component (due to self-desiccation); (c) hygral stress component (due to external drying); (d) total stress.
Figure 25. Simulated stresses at mid-deck location in the presence of constraint: (a) thermal stress component; (b) hygral stress component (due to self-desiccation); (c) hygral stress component (due to external drying); (d) total stress.
Materials 10 00231 g025
Figure 26. Simulated deck cracking: (a) without restraint; (b) with restraint associated with an internal diaphragm; and (c) crack opening histogram for the restrained case (at 12 pm on Day 16).
Figure 26. Simulated deck cracking: (a) without restraint; (b) with restraint associated with an internal diaphragm; and (c) crack opening histogram for the restrained case (at 12 pm on Day 16).
Materials 10 00231 g026
Table 1. Design parameters affecting early-age cracking of bridge decks.
Table 1. Design parameters affecting early-age cracking of bridge decks.
Primary CategorySubcategory
Materials composition/proportioningcementitious materials blend
admixtures
water content
aggregate
fiber reinforcement
Environmental boundary conditionsheat exchange: conduction, convection, radiation
moisture exchange
contaminant exposure
Processing and curingconcrete temperature at placement
time of placement
methods of consolidation
curing
Structural boundary conditionsgirder spacing and restraint conditions
deck depth
reinforcing steel
anticipated loading
Table 2. Modeling of coefficient of thermal expansion.
Table 2. Modeling of coefficient of thermal expansion.
β T / β Th Hydration Degree Interval
χ α α 0
χ ( χ 1 ) ( α / α 0 1 ) / ( φ 1 ) α 0 < α < φ α 0
1 α φ α 0
Table 3. Parameter settings for simulations of autogenous and drying shrinkage.
Table 3. Parameter settings for simulations of autogenous and drying shrinkage.
Parameter w / c = 0.25 w / c = 0.45
D 0 ˜ 0.017 mm 2 /h0.017 mm 2 /h
D 1 ˜ 9 mm 2 /h9 mm 2 /h
n65
h s u 0.760.88
s1.43.0
Λ h 0.25 mm/h0.25 mm/h
β h 0.00250.0021
Table 4. Parameter values associated with creep and stiffness development.
Table 4. Parameter values associated with creep and stiffness development.
ParameterValue *
q 1 29.0 × 10 6 /MPa
q 2 69.9 × 10 6 /MPa
q 4 5.50 × 10 6 /MPa
n α 2.2
κ 0 0.01/(MPa · d)
κ 1 5 MPa/K
* q 1 , q 2 , and q 4 values are based on the B4 model [85] for E 28 = 24.1 GPa; w / c = 0.42; a / c = 6.13; and 25% replacement of normal cement with fly ash.

Share and Cite

MDPI and ACS Style

Pan, Y.; Prado, A.; Porras, R.; Hafez, O.M.; Bolander, J.E. Lattice Modeling of Early-Age Behavior of Structural Concrete. Materials 2017, 10, 231. https://doi.org/10.3390/ma10030231

AMA Style

Pan Y, Prado A, Porras R, Hafez OM, Bolander JE. Lattice Modeling of Early-Age Behavior of Structural Concrete. Materials. 2017; 10(3):231. https://doi.org/10.3390/ma10030231

Chicago/Turabian Style

Pan, Yaming, Armando Prado, Rocío Porras, Omar M. Hafez, and John E. Bolander. 2017. "Lattice Modeling of Early-Age Behavior of Structural Concrete" Materials 10, no. 3: 231. https://doi.org/10.3390/ma10030231

APA Style

Pan, Y., Prado, A., Porras, R., Hafez, O. M., & Bolander, J. E. (2017). Lattice Modeling of Early-Age Behavior of Structural Concrete. Materials, 10(3), 231. https://doi.org/10.3390/ma10030231

Note that from the first issue of 2016, this journal uses article numbers instead of page numbers. See further details here.

Article Metrics

Back to TopTop