《1. Introduction》

1. Introduction

The alkali–silica reaction (ASR), which takes place in concrete material, is a mechanism that can be explained as a chemical reaction between the alkali ions contained in cement paste and the silica existing in aggregate pieces in the presence of water. The ASR results in an expansive gel, generally called "ASR gel,” which causes damage and cracking in concrete over time [1]. As a result of this long-term material deterioration, both the stiffness and the strength of concrete decay. The key factor in ASR occurrence is the existence of water, which is why the ASR mostly affects concrete structures that are located in environments with high humidity (usually above 60%). Temperature is another crucial factor controlling the driving force of the ASR [2].

Stanton [3] was the first to study the ASR effect: He investigated the chemistry of the reaction and the effect of the ASR at the material and structural levels, developed experimental techniques to assess the effect of concrete mix and aggregate distribution on ASR, and devised new solutions to mitigate ASR damage. Since then, numerous research efforts, including experimental, analytical, and numerical studies, have scrutinized ASR and considered various aspects of its effect on structural deterioration.

Laboratory ASR experiments are predominantly performed on small-scale specimens in accelerated ASR situations. The accelerated mortar bar test (AMBT) [4] is a popular experimental method that is performed over a period of approximately 16 d. In AMBT experiments, specimens cast from crushed aggregate pieces are immersed in alkali-rich solutions, and the temperature is increased to 80 °C. Another ASR experimental protocol is the concrete prism test (CPT) [5], which provides more comprehensive data on the influence of the ASR on concrete structures. The CPT is generally conducted on concrete samples cast from aggregates used in industrial applications without any adjustment. In addition, samples are kept in near-saturation humidity conditions (97% or above) in totally sealed containers to provide maximum moisture content and minimize possible leaching. CPTs fall into two categories with respect to experiment time duration: accelerated tests that occur over a period of 6 months are often performed at 60 °C, whereas longer tests that occur over a period of 1–2 years are performed at 38 °C.

Many researchers have explored numerical modeling of the ASR effect at various length scales and with different accuracy levels. Bazˇant et al. [6] first explored ASR damage and deterioration through a fracture mechanics approach and predicted the pessimum aggregate size.

The literature contains several macroscale continuum models that were developed to investigate the ASR mechanism and its relation to the degradation of concrete mechanical properties; of these, one of the earliest is a phenomenological model developed by Charlwood et al. [7] and Thompson et al. [8]. More detailed models, such as that portrayed in Ref. [9], along with others that take creep effect into account [10], were developed later on, and can be used to successfully predict stress and deformation history throughout the structure. In addition, researchers formulated computational models based on the kinetics of ASR and implemented them into various finite-element software programs using the crack band approach [11,12]. These models successfully replicated some of the experimental data available on ASR expansion [13]. The stress-state effect on ASR deterioration mechanism was incorporated into the models developed by Saouma and Perotti [14] and Multon et al. [15]. Comi et al. [16,17] developed a damage model that took into account certain mechanical and chemical aspects of the ASR process in a consistent thermodynamic fashion. Temperature and humidity factors were integrated into the kinetics of the ASR by Poyet et al. [18]. Finally, Bazˇant et al. [19] and Rahimi-Aghdam et al. [20] formulated a microplane model for ASR degradation in which the main degradation mechanism is the pressure induced by ASR gel flow. While all of the aforementioned models were deterministic, Capra and Sellier [21] proposed a probabilistic approach to study ASR, built upon the key parameters of ASR and on the mechanical behavior of concrete.

The mechanical models listed above share a common drawback that is the inability to accurately reproduce the ASR-induced crack pattern and damage distributions. This restricts the realistic prediction of the ASR degradation of concrete mechanical properties and necessitates the use of phenomenological formulations relating concrete mechanical properties to ASR gel expansion. These models also require phenomenological relationships to couple the applied stress state to the ASR gel expansion in order to reproduce the ASR effect in concrete under confinement. These restrictions are mainly due to the assumption that considers concrete as a homogeneous and isotropic material [22]. The first mechanical model able to successfully overcome these limitations is ASRLDPM [23], which is employed in the present work. In this model, ASR effect is combined with the lattice discrete-particle model (LDPM) [24,25]. LDPM is a discrete mechanical model that simulates concrete as an ensemble of three-dimensional (3D) rigid polyhedral cells, each of which represents a spherical aggregate piece and the surrounding mortar. These cells are linked through lattice struts, whose mechanical behavior is formulated in a vectorial form. In previous studies, it was shown that ASR-LDPM could successfully replicate the crack distribution and pattern in ASR-free expansion experiments, as well as multi-uniaxial compression experiments [23]. Furthermore, ASR-LDPM successfully reproduced the degradation of concrete mechanical properties in relation to the effect of temperature and the level of alkali content. It must be noted that in the ASR-LDPM formulation, no phenomenological relationship exists between the ASR expansion and the degradation of concrete mechanical properties. Recently, creep and shrinkage effects have also been incorporated into the model, as these are necessary to improve the model prediction under variable temperature and humidity [26,27].

Despite the significant capabilities of the ASR-LDPM framework, this model carries a major limitation: the immense computational cost of the simulation of large concrete structural members. Thus, in order to solve practical problems, it is imperative to develop a multiscale framework to decrease the computational cost of the analysis, while preserving the accuracy of the results. Wu et al. [28] established a homogenization framework based on the averaging theorem to analyze the ASR effect in concrete structures. However, no two-scale homogenization analysis was performed in their work, and the numerical accuracy and computational efficiency of the multiscale model was not studied. In the present work, a recently developed multiscale homogenization framework [29–32] is employed to investigate the ASR effect on concrete structures. The homogenization of an ASR-LDPM representative volume element (RVE) into a macroscopic tensorial constitutive law allows the simulation of macroscopic problems by means of the classical finite-element method, and leads to a significant reduction of the computational cost required to solve these problems.

《2. Review of the ASR-LDPM model》

2. Review of the ASR-LDPM model

The ASR-LDPM model [23] is formulated by introducing an ASRinduced mesoscale strain into the concrete mesoscale mechanical interaction modeled by LDPM. The evolution of the ASR strain over space and time is related to the chemo-physics of the problem, and the degradation induced is left solely to the LDPM constitutive behavior. The following subsections review the formulations of both the LDPM and ASR models, and then explain how they are coupled.

《2.1. The lattice discrete-particle model》

2.1. The lattice discrete-particle model

By using the mix design of concrete as input (i.e., cement content c, water-to-cement ratio w/c, and aggregate-to-cement ratio a/c), the mesoscale geometry of a generic LDPM system is constructed through the following steps:

(1) Coarse aggregate pieces, which are assumed to have a spherical shape, are inserted into the concrete sample volume through a trial-and-error procedure. The largest aggregates are placed in the sample first, and the aggregate sizes are calculated by sampling a Fuller curve F , in which d is the aggregate nominal diameter, da is the maximum aggregate nominal diameter, and n is the Fuller curve exponent. A certain minimum aggregate nominal diameter d0 is chosen as a lower threshold defining the model resolution. Following this procedure, a dogbone specimen filled with spherical aggregate particles can be generated, as shown in Fig. 1(a).

《Fig. 1》

 

Fig. 1. (a) Dogbone specimen including distributed spherical aggregates; (b) generic LDPM tetrahedron with four aggregate particles and its related tessellation; (c) tetrahedron segment related to node  ; (d) two adjacent LDPM polyhedral cells encompassing the associated spherical aggregates; (e) dogbone specimen generated from polyhedral cells.

(2) Computational nodes with zero radius are placed on the outer surfaces of the sample; these are the entities used to apply the boundary condition.

(3) A Delaunay tetrahedralization is performed based on the centers of the generated aggregates and the external nodes in order to generate a system of tetrahedra. A subsequent 3D tessellation of this tetrahedral system builds a network of triangular facets inside each tetrahedron, as depicted in Fig. 1(b). As a result of the tessellation, each tetrahedral element is divided into four segments, and each segment is connected to one node. The combination of the tessellation subdomains connected to a generic node yields a polyhedral cell that encloses that node and the corresponding spherical aggregate. Fig. 1(c) shows the subdomain associated to the node PI. The triangular facets are the locations where a vectorial constitutive equation is defined in terms strain to stress vectors. These facets are the potential planes on which cracking and damage can form in concrete, and are the locations through which two adjacent polyhedral cells interact, as shown in Fig. 1(d). Fig. 1(e) reports the polyhedral cell representation of the dogbone specimen depicted in Fig. 1(a). It should be noted that the spherical aggregates are only created to construct a geometry resembling the real mesoscale material structure of concrete, and they do not play a role in the numerical solution procedure. On the other hand, the centroid of the spherical particles and the associated polyhedral cells, which are respectively referred to as the “node” and “cells” hereinafter, are the geometrical entities that are used in the numerical analysis. In this paper, the cells and nodes are symbolized as and , respectively, where I = 1, ..., N, and N is equal to the total number of particles.

2.1.1. Facet strain vector definition

Considering a generic facet shared between two neighboring polyhedral cells with nodes and (Fig. 1(b)), the displacement jump vector can be computed at the centroid of the facet on the basis of the rigid-body kinematics describing the deformation of the system of cells. is then used to define the facet strain vector as follows:

where = N, M, L, as shown in Fig. 1(c). The facet strain vector formulation in Eq. (1) is written by assuming small displacements, rotations, and strains. In the rest of this paper, as in Eq. (1), the superscript IJ is written only if the sign of that quantity changes by interchanging I and J, such as for facet normal or tangential vectors . As explained in Ref. [24], the local system of reference defining the normal and tangential orientations is relevant to the triangular facet projected orthogonally to the lattice line connecting the two cells’ nodes, as opposed to the original facet orientation. It is worth mentioning that the strain vector definition in Eq. (1) is consistent with the projection of the continuum micropolar strain tensor onto a generic orientation in space [33].

In Eq. (1), given that  is the normal component of the facet strain, and  are the tangential components of the facet strain. As shown in Fig. 1(b), is the distance between the centers of two adjacent particles and , and it can be calculated as  of the facet. are the unit vectors describing the facet local Cartesian reference system:  is the unit vector orthogonal to the projected facet, such that ,and and  are the unit vectors in the tangential direction to the facet. The vectors that connect the two nodes and to the centroid of the facet are denoted by  , respectively. Figs. 1(b) and (c) can be referred to for all these geometrical entities. are the displacement vectors of the nodes  and , and are the rotation vectors of the nodes and .

2.1.2. Facet constitutive equations

A constitutive equation, formulated at each facet centroid, describes the material behavior and relates the facet strain vector to the facet stress vector. The LDPM constitutive equation is briefly summarized below but reader is referred to previous work [24,34] for additional details. To simulate the elastic behavior, the facet normal and shear stress components are directly proportional to the corresponding facet normal and shear strain components.

To simulate fracture under tension and shear-tension (), an effective strain  and an effective stress  are calculated to define the normal and shear stresses on the facet as , and . The effective stress t is incrementally elastic  and it is bounded by the inequality . The stress boundary is , in which , and =, where  and   is the parameter that governs the coupling degree between the normal and the shear loading.  is the maximum effective strain experienced by the facet, which is  when there is no unloading. The softening modulus for the post-peak behavior ,  in which nt is the softening exponent and Ht is the softening modulus in pure tension (when = π /2), formulated as Ht = is the tensile characteristic length and Gt is the mesoscale fracture energy. The transition between pure shear  to pure tension occurs through a smooth function over which the strength varies in a parabolic form given by , where  is the shear-to-tensile strength ratio.

On the contrary, when the facet is in compression , the normal stress is calculated considering the inequality  where  is the boundary stress function that depends on the element volumetric strain, ,  and the facet deviatoric strain, . The volumetric strain is calculated by the volume change of the tetrahedral element as , which is considered to be equal for all facets belonging to that tetrahedral element. The function , describing pore collapse and rehardening under high confinement, is defined as follows:  for ;  and  for   and  for . In the previous equations,  and   for ; , and  are material parameters.

For facets in compression, the shear stresses are calculated incrementally as  and , where  and   is the plastic multiplier for which the loading and unloading conditions are  and . The yielding surface and the plastic potential are written as  with .  and , respectively.  is the transitional normal stress, is the initial internal coefficient of friction, and is the final internal coefficient of friction.

2.1.3. Cell equilibrium equations

Finally, the translational and rotational equilibrium equations of a generic LDPM cell can be written as:

where  is the series of facets that cell  is composed of; A is the area of the projected facet;  is the cell volume;  is the body force vector, which is considered to be uniform over the cell; and  is the resultant stress vector acting on the triangular facet.

LDPM is a part of a computational package named MARS [35] and has been widely used for the simulation of concrete mechanical behavior [25]. In addition, LDPM has shown great success in the simulation of concrete response in a dynamic regime [36,37], fracture modeling of concrete reinforced by fiber-reinforced polymer (FRP) and steel [34], simulation of fiber-reinforced concrete [38,39], and concrete perforation response [40].

《2.2. ASR governing equations》

2.2. ASR governing equations

Over the past few decades, researchers have broadly investigated the chemical and physical processes that occur during ASR progression. It has been shown that those processes significantly depend on aggregate mineralogy and chemistry, cement composition, cement replacement products, and other additives. Furthermore, it has been reported that silica is not uniformly distributed over each aggregate piece; rather, it appears in the form of discrete inclusions, pockets, and veins [41,42]. Various types of alkali ions, such as Na+ , K+ , and Ca2+, which exist in the cement paste [43,44], react with the silica inside the aggregate pieces in the presence of hydroxide ions (OH-) and water (H2O) [45-47], mainly later in the concrete lifetime [48,49]. The ASR results in an amorphous expansive gel that varies widely in chemical composition due to different aggregate mineralogies and available alkali ions [50,51]. In addition, researchers advise that the ASR gel flowability and expansive characteristics significantly depend on the calcium and water content [52–54], and that further chemical reactions take place at later stages [55]. The produced ASR gel shows a considerable expansive characteristic, which results in significant cracking and damage in the aggregates and cement paste [56,57]. Experimental inspections have clearly demonstrated the presence of ASR gel, particularly on the outer surface of aggregates, inside aggregate pieces, and inside the cracks in the case of very reactive aggregates. It should be noted that this does not necessarily indicate gel flow, since it may be related to the fact that cracks facilitate the transport of water and ions, which would later produce more gel inside the crack [52,53,58].

Nevertheless, during the ASR, the alkali ions available in the cement paste must permeate into the aggregate in order to react and form ASR gel; furthermore, additional water must be transported to the gel in order for it to imbibe water and subsequently expand. The aforementioned statements clearly indicate that a mesoscale diffusion process must occur. Therefore, in previous research, a mesoscale concrete model was selected that can capture the major ASR phenomena and the averaged subscale ones [23]. In this approach, the overall average expansion rate of a single aggregate piece, which depends on basic gel formation and water imbibition, is converted into an inhomogeneous gel strain and applied as an eigenstrain on the concrete meso-structure within the LDPM framework. ASR-LDPM was originally formulated for saturated conditions and uniform alkali concentration [23], and it has recently been generalized for the case of variable environmental conditions and variable alkali ion concentration [26].

As mentioned earlier, in order to generate the ASR gel, alkali ions and water must reach the silica available in the aggregates by means of a diffusion process through the volume of aggregates. Therefore, the mass of the gel Mg produced around a generic aggregate piece of diameter d is calculated by solving a steady-state mass balance of the radial diffusion process into the aggregate body. A schematic of the diffusion process that takes place during the ASR is plotted in Figs. 2(a) and (b), in which z is the position of the diffusion front calculated from the aggregate center, which is equal to the radius of the unreacted segment of the aggregate. Thus, the produced ASR gel mass Mg can be formulated as follows:

《Fig. 2》

Fig. 2. (a) A generic aggregate piece and the surrounding cement paste; (b) ASR gel created around a generic aggregate piece; (c) the triangular element generated by three aggregate pieces in a 2D model under the ASR effect.

where is a parameter that accounts for the level of alkali content in the cement paste. is the lower bound of the alkali content, at which little or no ASR expansion is recorded, and is the saturation alkali content sufficient for a complete ASR. is a free parameter used for the silica-to-gel conversion. The speed of the ASR diffusion front can be written as follows:

where T0 and T are the reference and current temperatures, respectively; Eag is the diffusion process activation energy; and R is the universal gas constant. ww is the water content density in the cement paste in the vicinity of the aggregate pieces, which can be estimated as  at saturation [59]; here, = is the asymptotic hydration degree.  is the permeability of alkali-rich water into the aggregate at room temperature (296 K).

Once created, ASR gel tends to absorb water and expand as a result, thus applying an expansive pressure to the surrounding material that leads to damage and cracks. The process of water imbibition is formulated by the following evolution law, as presented in Ref. [26]:

where is the mass of imbibed water, is the ratio between the maximum amount of imbibed water and the mass of gel, is the initial micro-diffusivity of the water in the vicinity of the aggregate surface, and is the effective distance for water to travel from the concrete surrounding an aggregate into the ASR gel. The exponential term reflects the fact that the diffusivity of the ASR gel decreases as the water imbibition process proceeds. Finally, the temperaturedependence of the imbibition process is formulated in the form of an Arrhenius-type equation expressed in terms of the imbibition process activation energy .

The calculation of the mass of imbibed water permits then the calculation of an average expansion of the corresponding aggregate piece, whose radius becomes

where the brackets <> extract the positive value of the enclosed expression. In Eq. (6),  is the imbibed water associated to a gel expansion that can be accommodated in the interfacial transition zone (ITZ) pores without generating any pressure, r = D/2 is the initial radius of the aggregate piece,  characterizes the available pore space in the ITZ, and is the mass density of water.

Considering Eq. (6), the time rate of the radius of a generic aggregate piece can be written as follows:

Therefore, the normal eigenstrain time rate generated by the ASR gel expansion  on a generic facet in the LDPM mesh is determined as follows:

where  are the time rate of the radius of aggregate pieces I and J, which share the facet under consideration, as illustrated in Fig. 2(c).

It should be noted that the shear eigenstrains generated by the ASR gel expansion are assumed to be zero, , This is only an approximation of the actual situation due to the irregular shape of real aggregate pieces.

《3. The multiscale homogenization method》

3. The multiscale homogenization method

This research uses the homogenization technique formulated in Ref. [30] to calculate the average mechanical response of an LDPM RVE affected by ASR. This average response is then used to perform standard finite-element calculations that account for the effect of concrete heterogeneity at the mesoscale.

With reference to Fig. 3, two coordinate systems, X and x, are considered in the macroscale; these respectively represent a global macroscopic coordinate system and a local coordinate system defining the position within the RVE, as shown in Fig. 1(e). In addition, a mesoscale coordinate system y is also considered in which the material heterogeneity is fully resolved. In X and x, the material domain is assumed to be homogeneous, with no material heterogeneity being visible at this scale. According to the scale separation assumption [60,61], the relation holds when g is a small positive scalar, . It should be noted that in Fig. 1(e), the separation of scales assumption is visually violated in order to provide a clear representation of the two-scale problem.

《Fig. 3》

Fig. 3. The two-scale homogenization scheme. (a) Macroscopic material domain; (b) mesoscale domain including the material heterogeneities; (c) two adjacent polyhedral cells.

Two independent field variables, u(x, y) and (x, y), which respectively represent the displacement and the rotation vectors of a generic computational node , are considered. These variables are functions of the macroscale and mesoscale coordinate systems and they correspond to the mesoscale degrees-of-freedom (DOFs) when evaluated for . These variables can be expanded in regard to the two separate scales as follows:

where are respectively the macroscale and fine-scale displacement fields, which are continuous in x and discontinuous in y. Furthermore, since the rotation vector can be related to the displacement vector through a curl operator, it is straightforward to show that  and  are fine-scale rotation vectors, while and are macroscale rotation fields [30].

Using the asymptotic expansion of the displacement and rotation vectors in Eq. (1) and considering the rigid body motion of the fine-scale RVE [30], the asymptotic expression of the LDPM facet strains is obtained in the following form:

By neglecting the macroscopic curvatures effect, the following is obtained:

The facet strains of order zero, ,  consist of three different terms:  are the fine-scale facet strains, which correspond to Eq. (1), while the second and third terms,  ,account for the projections of the symmetric macroscopic strain tensor , which is defined as , and the antisymmetric part of the macroscopic strain tensor , which is defined as . The tensor is the Levi–Civita permutation symbol. and  are the projection operators mathematically derived as  and  . In essence, the second and third terms in Eq. (12) represent how the macroscopic strain tensors at each finite-element integration point should be projected onto the corresponding LDPM RVE facets in order to solve the RVE problem. A detailed derivation of the theory can be found in Ref. [30]. 

The asymptotic expansion of the LDPM facet stresses can be derived by using the asymptotic form of the facet strains; together with Eq. (2), this results in the derivation of translational and rotational equilibrium equations for two separate scales: order zero, which is the RVE problem, and order one, which is the macroscopic or coarse-scale problem.

The RVE translational and rotational equilibrium equations are as follows:

where the zero-order facet stresses  are calculated using the constitutive equations presented in Section 2.1 and according to the zero-order facet strains . Eq. (12) can be rewritten as , which implies that the projection of the macroscopic strain tensor is applied in the form of facet eigenstrains that drive the RVE problem. By solving the RVE equilibrium equations using the periodic boundary conditions, the zero-order traction can be obtained; this is then used to calculate the macroscopic symmetric and antisymmetric stress tensors as follows:

where V0 is the volume of the RVE. The homogenization theory also provides the macroscopic couple stress tensor, which, however, is not reported here because it has been shown to be negligible [62].

Finally, the macroscopic equilibrium equations can be obtained by averaging the equilibrium equations of order one [30], as follows:

Since the antisymmetric part of the stress tensor and, as mentioned earlier, the couple stress tensor are zero, the antisymmetric part of the strain tensor must be also zero: . This permits the numerical implementation of the LDPM-based homogenization framework with no macroscopic rotational DOFs and with standard displacement-based finite elements.

Fig. 4 reports an overview of the adopted homogenization framework, which includes the mesoscale ASR effect according to the following steps:

《Fig. 4》

Fig. 4. General framework of the multiscale homogenization method with the ASR effect. FE: finite-element.

(1) The macroscopic material domain is discretized by finite elements, and an LDPM RVE is assigned to each macroscopic finite-element Gauss point.

(2) For each macroscopic loading step, an increment of the macroscopic strain tensor is calculated at every Gauss point of all the macroscopic finite elements.

(3) At each finite-element Gauss point, the calculated macroscopic strain tensor increment is projected onto all the corresponding LDPM RVE facets using each facet orientation using (Eq. (12)).

(4) The increments of the ASR-generated normal strain are calculated for all facets inside each LDPM RVE.

(5) The RVE problem is solved by applying to all RVE facets under the periodic boundary conditions applied on the RVE. The ASR-generated normal strains, along with the projected macroscopic strains, are applied as eigenstrains on all of the RVE facets, leading to the calculation of the fine-scale solution governed by the LDPM constitutive equations:, in which  is 1 if = N and 0 if = M or L. The mesoscale displacement and rotation vectors of all RVE particles and the traction vectors t 0 a of all RVE facets are then calculated.

(6) The macroscopic stress tensor is calculated using Eq. (14), based on the RVE problem solution.

(7) The RVE homogenized stress tensor is then transferred back to the corresponding Gauss point in order to update the finiteelement nodal displacements and forces.

《4. Numerical results》

4. Numerical results

In this section, the performance of the developed multiscale framework to assess the ASR effect on concrete in terms of volumetric expansion, along with mechanical properties degradation, is investigated. First, the ASR-LDPM parameters are calibrated according to the reported experimental data by Shehata and Thomas [63]. They studied the free expansion of sealed concrete prisms with the dimensions 75 mm × 75 mm × 265 mm for 320 d. In their work, five different levels of alkali content ( = 2.89, 3.15, 3.90, 4.20, and 5.25 kg·m-3 ) were considered. The reported aggregate size distribution of the concrete used in the experiments is as follows: one third of the aggregate particles are between 4.75 and 9.5 mm, one third are between 9.5 and 12.5 mm, and the last third are between 12.5 and 19 mm. For the generation of LDPM specimens, the following parameters are used: minimum aggregate size = 4.75 mm, maximum aggregate size = 19 mm, and Fuller curve exponent n = 0.55; and mix composition parameters with cement content c = 420 kg·m-3 , water-tocement ratio w/c = 0.45, and aggregate-to-cement ratio a/c = 4.25. To calibrate the ASR-LDPM parameters, concrete prisms are fully simulated with LDPM and their volumetric expansion is studied over time. A cell representation of an LDPM prism is shown in Fig. 5(a). For each level of alkali content, a full LDPM analysis is performed with three different LDPM particle realizations in order to investigate the RVE mesh realization effect, and the average behavior is computed. The experimental volumetric strain versus time and the corresponding full fine-scale LDPM analysis results are plotted in Fig. 5(b). It can be seen that the numerical simulation results match well with the experimental data. It is evident that the volumetric strain rate decreases and plateaus as the ASR continues. It should be noted that the ASR-LDPM parameters are first calibrated for an alkali content level of 2.89 kg·m-3 , and are then used to validate the response of the model for other alkali content values. The calibrated ASR parameters are: =2.7 kg·m-3 , = 4.37 kg·m-3 ,  = 3.45 × 10 ×13 m2·s-1 ,  = 2.8 ×10-10 m2·s-1 ,  = 439540 kg-1 , kg = 689 kg·m-3 , Eag = Ead = 500 J·mol-1 , and =1. The reference temperature e T0 and the current temperature T are considered to be equal, which means that the temperature effect is neglected. It should be noted that the specimens are sealed, and that they maintain a high relative humidity close to 100%. The LDPM parameters used in the simulation are: EN = 60 GPa, = 4.75 MPa,  = 150 MPa, = 0.25, = 0.2, = 75 mm, = 2.6, = 0.4,  = 0.4,  = 0, = 1, = 5, and = 600 MPa. These parameters are those that have been calibrated for standard concrete mechanical behavior [23], since no mechanical properties were reported in Ref. [63].

《Fig. 5》

Fig. 5. (a) Polyhedral cell representation of a concrete prism; (b) volumetric expansion of the concrete prism obtained from the full fine-scale LDPM analysis and the experiment (the values in image (b) are the alkali content with unit of kg·m-3 ).

Since the fine-scale RVE is the principal element of the homogenization framework, the LDPM RVE homogenized response under tension and compression and affected by the ASR is examined first in the next section. The concrete prism studied experimentally in Ref. [63] is then modeled using the homogenization method in order to determine the efficiency of the multiscale framework in comparison with the full fine-scale simulations and experimental data. However, it must be noted that the aim of this paper is to scrutinize the efficiency of the multiscale framework in order to accurately reproduce the results of the full LDPM analysis with a lower computational cost, rather than validation with respect to the experimental data. The latter task was already presented in Refs. [23,27]. In addition, other phenomena such as shrinkage and creep that occur during the concrete lifetime are neglected in this study, although their incorporation into the multiscale model lies within the scope of our future work.

《4.1. Analysis of RVE homogenized mechanical properties under the ASR effect》

4.1. Analysis of RVE homogenized mechanical properties under the ASR effect

In this section, LDPM RVEs are generated with the mix composition parameters presented above. To study the effect of RVE size, four different RVE sizes of 35, 50, 75, and 100 mm are considered. The cell and particle distributions of the four RVE sizes are presented in Figs. 6 and 7. For each RVE size, five different particle distributions are considered in order to examine the effect of the RVE mesh realizations. In Fig. 8, five different particle realizations are presented for a 75 mm LDPM RVE.

《Fig. 6》

Fig. 6. Polyhedral cell representations of RVEs of different sizes. (a) 35 mm; (b) 50 mm; (c) 75 mm; (d) 100 mm.

《Fig. 7》

Fig. 7. Spherical particle distribution inside RVEs of different sizes. (a) 35 mm; (b) 50 mm; (c) 75 mm; (d) 100 mm.

《Fig. 8》

Fig. 8. Five different particle realizations inside 75 mm RVE.

The generated RVEs are tested under uniaxial tension at three different stages of the ASR evolution. First, they are tested with no ASR effect. Next, the RVEs are tested after 60 and 120 d of ASR free expansion. The ASR expansion of the RVEs is simulated using the parameters presented in Section 4, and the alkali content considered is 2.89 kg·m-3 . The ASR expansion of the samples at 60 and 120 d is approximately 0.015% and 0.025%, respectively. It should be noted that these values are about the same for all RVEs of four different sizes. The effect of RVE size on the free-expansion curves and its comparison with the experimental data are presented in the next section. The crack-opening contours of the RVEs due to free expansion after 120 d of the ASR are presented in Fig. 9. It can be seen that the ASR damage is diffused over the specimen and the resulting cracks are randomly oriented, which corresponds to the damage patterns observed in free-expansion experiments.

《Fig. 9》

Fig. 9. Crack-opening contours of the RVEs due to free expansion after 120 d. (a) 35 mm; (b) 50 mm; (c) 75 mm; (d) 100 mm.

Damaged RVEs at 60 and 120 d are then subjected to uniaxial tension, and the obtained stress–strain curves are compared with that obtained from the specimen with no ASR damage. The tensile stress–strain curves are averaged for different particle realizations of each RVE size, as presented in Fig. 10. From each curve, the Young’s modulus Et and the tensile peak stress are calculated and plotted versus the specimen ASR-affected time, as shown in Fig. 11. All curves are slightly shifted along the time axis to provide a clearer view of the error bars. At each time, it can be seen that the calculated average Et and  for different RVE sizes match very well. Therefore, it demonstrates that the degradation of the Young’s modulus and tensile strength of concrete due to ASR are perfectly captured by the homogenization scheme with a minimum effect of the RVE size. For each RVE size, the error bar shows the standard deviation of the quantity of interest from five different particle realizations. As expected, the standard deviation is higher for a smaller RVE size, since the ratio of the RVE size to the maximum aggregate size is smaller. The error bar value is negligible for the 100 mm RVE size, which implies that particle distribution inside the RVE has no effect on the homogenized quantity. These results clearly show that smaller RVE sizes can be employed in multiscale homogenization analysis to successfully capture the degradation trend of the tensile mechanical properties.

《Fig. 10》

Fig. 10. Average curves of RVE nonlinear behavior under tension (a) with no ASR evolution, (b) at 60 d of the ASR evolution, and (c) at 120 d of the ASR evolution for different RVE sizes.

《Fig. 11》

Fig. 11. Variation of (a) Young’s modulus in tension and (b) tensile peak stress over time due to ASR.

The different post-peak responses observed in Fig. 10 for different RVE sizes are known as the ‘‘size effect” in the fracture mechanics of quasi-brittle materials such as concrete. This effect is due to the fact that damage localization occurs in quasi-brittle materials under tension, while the elastic energy of the undamaged parts of the domain is released. The damage localization patterns for different RVE sizes are plotted in Fig. 12. This issue is a well-known problem in the multiscale homogenization of softening materials that requires special treatment in order for the multiscale framework to correctly capture the dissipated energy [64].

《Fig. 12》

Fig. 12. Crack-opening contours in a direct tension test under an axial strain of 1 × 10-3 . (a) 35 mm; (b) 50 mm; (c) 75 mm; (d) 100 mm.

Next, the evolution of the compressive strength and the Young’s modulus of concrete due to ASR is calculated by means of uniaxial compression tests performed on RVEs after 182 d and 2 years of the ASR evolution. The compressive response of the RVEs is compared with that obtained from the RVEs with no ASR effect. All the LDPM material parameters, as well as the ASR parameters, are those that were used in the previous section. Since degradation of the concrete compressive strength due to ASR is less pronounced than the degradation of the tensile strength [56], the RVEs are subjected to longer periods of ASR evolution compared with those described in the previous section.

The compressive stress–strain curves are averaged for different particle realizations of each RVE size, and are presented in Fig. 13 for no ASR, 182 d of ASR free expansion, and 2 years of ASR free expansion. The evolutions of the Young’s modulus Ec and compressive peak stress due to ASR are plotted in Fig. 14, in which all curves are slightly shifted along the time axis to provide clear visibility of the error bars. It can be seen that the homogenized values of Ec and obtained from the analysis of different RVE sizes match very well at the three different times. Therefore, smaller RVE sizes can be used in multiscale homogenization analysis while maintaining the accuracy of the results. Similar to the previous section, the standard deviation of the homogenized Ec and decreases as the ratio of the RVE size to the maximum aggregate size is increased. However, the standard deviation of these quantities is still negligible for even the smallest RVE size. It should be noted that, as opposed to the tensile behavior, the post-peak behavior of all RVE sizes under compression match very well, as shown in Fig. 13. This is due to diffused crack distribution and a lack of damage localization in the RVEs under uniaxial compression, as depicted in Fig. 15.

《Fig. 13》

Fig. 13. Average curves of RVE nonlinear behavior under compression (a) with no ASR evolution, (b) at 182 d of the ASR evolution, and (c) at 2 years of the ASR evolution for different RVE sizes.

《Fig. 14》

Fig. 14. Variation of (a) the Young’s modulus in compression and (b) the compressive peak stress over time due to ASR.

《Fig. 15》

Fig. 15. Crack-opening contours in an unconfined compression test under an axial strain of 8  103 . (a) 35 mm; (b) 50 mm; (c) 75 mm; (d) 100 mm.

《4.2. Two-scale homogenization analysis of an ASR-affected concrete prism》

4.2. Two-scale homogenization analysis of an ASR-affected concrete prism

In this section, the concrete prisms that were studied experimentally in Ref. [63] and presented in Section 4 are simulated by means of both the full fine-scale LDPM and the two-scale homogenization method in order to investigate the ability of the multiscale model to reproduce the full fine-scale analysis results. The ASR-LDPM parameters that were calibrated with respect to the experimental record of the prisms’ free expansion, as shown in Fig. 5(b), are employed here for both the fine-scale and multiscale homogenization analyses. In the homogenization analyses, the concrete prism is discretized with three single-integration-point hexahedral finite elements with an edge length of 75 mm, as shown in Fig. 16. To investigate the effect of RVE size, the homogenization analysis is repeated with three different RVE sizes of 35, 50, and 75 mm, which are assigned to the integration point of each finite element. In addition, for each RVE size, a two-scale analysis is carried out for six different particle realizations inside the RVE in order to examine the effect of the particle distribution inside the RVE.

《Fig. 16》

Fig. 16. Homogenization model of the concrete prism with three solid finite elements and three different RVE sizes.

The volumetric expansion of the concrete prisms obtained from the two-scale homogenization analysis is plotted versus time for the three different assigned RVE sizes in Figs. 17(a)–(c). A freeexpansion analysis is performed for the five different alkali contents, as reported in the experiments [63]. For each RVE size and alkali content value, six different curves are plotted corresponding to six different particle distributions inside the RVE. It is clear that the curves are more scattered for smaller RVE sizes, as was observed in the previous section. Therefore, for a higher ratio of RVE size to maximum particle size, the standard deviation of the results is smaller. Moreover, it can be observed that the expansion curves are more scattered for greater level of alkali content, since greater alkali content leads to a higher level of ASR damage and to subsequently enhanced inhomogeneity of damage distribution inside the RVE.

These curves are averaged for each RVE size and different alkali contents, and are compared with the full fine-scale LDPM analysis results in Figs. 17(d)–(f). It can be seen that the homogenization and full fine-scale analysis results match well for different RVE sizes and different alkali content values. The homogenization results agree more closely with the full fine-scale simulations as the RVE size is increased. This result is expected, since larger RVEs are better representatives of the homogenized lower-scale material structure.

《Fig. 17》

Fig. 17. Volumetric expansion of prisms with differing alkali contents obtained from the homogenization framework using RVE sizes of (a) 35 mm, (b) 50 mm, and (c) 75 mm; comparison of the homogenization results using RVE sizes of (d) 35 mm, (e) 50 mm, and (f) 75 mm with the full fine-scale results (the values in images (d)–(f) are the alkali content with unit of kg·m-3 ).

Next, the same prisms are loaded under uniaxial compression after different periods of ASR evolution in order to investigate the ability of the multiscale homogenization method to capture the degradation of the concrete compressive strength in comparison with the full fine-scale LDPM analysis. First, the specimen is subjected to the ASR effect over a certain time interval, while being constrained in the axial direction. Next, the prism is loaded under uniaxial compression until it reaches failure. To determine the trend of compressive strength degradation, specimens are subjected to the ASR effect for 0, 122, 243, 365, 488, 610, and 730 d while being axially constrained. All tests are carried out for ca equal to 2.89 and 3.15 kg·m-3 in order to examine the effect of alkali content. Figs. 18(a) and (b) compare the evolution of the concrete prism compressive strength obtained from the full LDPM simulations with the homogenization analyses for three different RVE sizes. It can be seen that the results of the homogenization analysis for different RVE sizes correspond well to the results obtained from the full fine-scale LDPM analysis. In addition, it is clear that the degradation of the prism compressive strength is more significant for greater alkali content, and that this finding is successfully captured by the homogenization framework, in agreement with the full LDPM results. Furthermore, the effect of particle realization decreases for larger RVE sizes, as observed in previous results, which is shown by the decrease of the scatter on the curves.

《Fig. 18》

Fig. 18. Comparison of the compressive strength obtained by the full fine-scale LDPM analysis and the homogenization framework for (a) ca = 2.89 kg·m-3 and (b) ca = 3.15 kg·m-3 ; (c) comparison of the computational cost of the full fine-scale LDPM analysis and the homogenization analysis for different RVE sizes.

The computational cost study is one of the most important aspects of the multiscale methods, without which the development of multiscale models has no purpose. In the literature, the computational cost efficiency of developed multiscale models has not been investigated as much as their numerical accuracy. The significant advantage of the current developed homogenization framework is its enormous saving in computational cost, while retaining the numerical accuracy of the results. Fig. 18(c) compares the computational cost of the full fine-scale LDPM analysis with that of the homogenization method for different RVE sizes. The computational time of the homogenization analysis of the concrete prism using the 75 mm RVE is approximately 30% greater than that of the full fine-scale LDPM analysis. This is because the total number of the computational nodes and LDPM tetrahedral elements in the homogenization analysis for the 75 mm RVE is greater than in the full fine-scale LDPM analysis. In this case, the volume of the RVE is equal to the macroscopic finite-element size, while the RVE surface nodes are duplicated for the finite-element internal surface. On the other hand, the computational time of the multiscale analysis for the 50 and 35 mm RVEs is 67% and 23%, respectively, of the full fine-scale LDPM analysis computational time. It should be noted that only three macroscopic finite elements are used in the homogenization analysis. Therefore, the computational cost saving in the analysis of larger structures discretized by a large number of finite elements will be more substantial. This significant computational cost saving, which is one of the major purposes of developing a multiscale framework, is successfully accomplished by means of the developed homogenization scheme.

Finally, the failure pattern of the prisms obtained from the full LDPM and homogenization analyses is depicted in Fig. 19 under both free expansion and uniaxial compression. Since the strain field is uniform in both the free expansion and uniaxial compression cases, the damage contour of just one RVE of each size is illustrated for both cases. Fig. 19(a) shows the crack-opening contour of a specimen after 730 d of ASR free expansion (without axial constraint) for ca = 3.9 kg·m-3 . It can be seen that the damage is distributed throughout the specimen in a random pattern, and a similar damage pattern is captured by the homogenization analysis. The crack-opening contour of a prism loaded under uniaxial compression for 243 d of ASR expansion is illustrated in Fig. 19(b). It can be observed that both the fine-scale and homogenization analyses capture the vertical splitting cracks due to the specimen lateral expansion that is generated by uniaxial deformation.

《Fig. 19》

Fig. 19. Crack-opening contours of a concrete prism and RVEs of different sizes under (a) free expansion due to the ASR effect with ca = 3.9 kg·m-3 and (b) uniaxial compression after ASR expansion.

《4.3. The four-point bending test (4pbt)》

4.3. The four-point bending test (4pbt)

In this numerical example, a concrete beam with a span of 2.5 m and a height of 0.5 m is simulated by both the full LDPM and the two-scale homogenization method, as shown in Fig. 20(a). For the homogenization analysis, the concrete specimen is discretized by hexahedral finite elements with an edge length of 50 mm and a single Gauss point, and an LDPM RVE of 35 mm is assigned to each macroscopic integration point. Both the LDPM RVE and the full concrete beam LDPM specimens are generated using the parameters presented in Section 4. The concrete beams are tested under the four-point bending loading condition for three different ASRaffected times: no ASR, 100 d of ASR free expansion, and 150 d of ASR free expansion. The ASR-LDPM parameters used in both the full fine-scale and the homogenization analyses are the same as those presented in Section 4. An alkali content of 2.89 kg·m-3 is used in all simulations. The force–displacement (F–u) curves obtained from the full LDPM and homogenization analyses are compared in Fig. 20(b). It can be seen that the F–u curves for concrete beams with no ASR effect match very well for the full LDPM and the two-scale homogenization analyses in terms of the force peak value, elastic branch, and area under the F–u curve. The difference in the post-peak slope is due to the fact that the RVE size is smaller than the macroscopic finite-element size in the homogenization analysis. Gitman et al. [64] has shown that the postpeak branch of the softening response is controlled by the RVE size, and that it accounts for the correct dissipated energy if the RVE size is equal to the finite-element size. 

In Fig. 20(b), the F–u curves obtained from the four-point bending test (4pbt) on the specimens after 100 and 150 d of ASR are also presented. It can be seen that the peak force drop due to the ASR calculated by the full LDPM analysis is higher than that obtained from homogenization analysis for both 100 and 150 d of ASR. An important factor with a significant impact on the mechanical response is the crack propagation path through the specimen during failure. The crack-opening contours in LDPM specimens after 4pbt on the beams with different ASR periods are shown in Figs. 21(a)–(c). In Figs. 21(b) and (c), distributed microcracks can be seen throughout the beam volume generated during the preloading ASR free-expansion phase. Fig. 21(d) presents the contour of the strain tensor component in the direction of the beam length obtained from the homogenization analysis at failure. This contour looks similar in all three cases (no ASR, 100 d of ASR, and 150 d of ASR). It is clear from this figure that in the full LDPM analysis, due to the fine spatial discretization of the specimen, macrocracks follow the path that minimizes the potential energy. In the homogenization analysis, however, the macrocracks are forced to grow straight in the vertical direction following the macroscopic finite element discretization. This difference is more important for the specimens with ASR free expansion before mechanical loading, because of the existence of distributed microcracks in the specimen. Therefore, homogenization constrains the macroscopic crack path by the macroscopic finite-element discretization, which is the likely reason for the difference between the full LDPM and multiscale F–u curves plotted in Fig. 20(b) for the cases of 100 and 150 d of ASR. Regarding the computational cost, the simulation time of the homogenization analysis in all three cases is approximately 40% of the full LDPM analysis run time, which confirms the significant improvement in computational efficiency provided by the multiscale analysis.

《Fig. 20》

Fig. 20. (a) Concrete beam modeled by full LDPM and finite element for homogenization analysis; (b) force–displacement curves obtained from the full LDPM and the homogenization analyses after no ASR, 100 d of ASR free expansion, and 150 d of ASR free expansion.

《Fig. 21》

Fig. 21. Crack-opening contours of the 4pbt on the concrete beam after (a) no ASR, (b) 100 d of ASR, and (c) 150 d of ASR; (d) longitudinal strain contour of the homogenization simulation.

《5. Conclusions》

5. Conclusions

This paper presents the theory and numerical implementation of a multiscale homogenization analysis of the ASR in concrete. In the homogenization framework, a mesoscale particle model, which has shown significant success in modeling the ASR effect in concrete, is linked to the continuum finite-element model. Damage and cracks in concrete due to mechanical loading, along with ASR damage, are all simulated in the fine-scale RVE, and no predefined nonlinear behavior is considered at the macroscale. First, it is shown that RVEs of different sizes are able to reproduce the degradation of the mechanical properties of concrete, both in tension and compression. Therefore, the size of the RVE can be chosen to be smaller than the corresponding finite-element size, which leads to a considerable saving in computational cost. Next, free expansion of a concrete prism due to the ASR is analyzed using the multiscale homogenization framework.

It is shown that the multiscale framework can accurately reproduce the expansion curves obtained from experiments and full fine-scale LDPM analysis. In addition, degradation of the compressive strength of a concrete prism is evaluated after several ASRaffected time periods, and it is demonstrated that the results of the multiscale homogenization and full fine-scale LDPM analyses agree well. The superior advantage of this multiscale analysis is its considerable reduction in computational cost. It is shown that the computational time of the prism analysis can be significantly reduced with respect to the full fine-scale LDPM analysis by decreasing the RVE size, while the numerical accuracy of the results is greatly preserved. Finally, the numerical accuracy and computational efficiency of the multiscale framework in comparison with the full fine-scale LDPM analysis are examined and confirmed by a simulation of the 4pbt on large concrete beams.

《Acknowledgement》

Acknowledgement

This material is based upon work supported by the Nuclear Regulatory Commission (NRC-HQ-60-14-G-0003).

《Compliance with ethics guidelines》

Compliance with ethics guidelines

Roozbeh Rezakhani, Mohammed Alnaggar, and Gianluca Cusatis declare that they have no conflict of interest or financial conflicts to disclose.