Macroscopic/Mesoscopic Computational Materials Science Modeling and Engineering
View this Special IssueResearch Article  Open Access
Combined NumericalStatistical Analyses of Damage and Failure of 2D and 3D Mesoscale Heterogeneous Concrete
Abstract
Generation and packing algorithms are developed to create models of mesoscale heterogeneous concrete with randomly distributed elliptical/polygonal aggregates and circular/elliptical voids in two dimensions (2D) or ellipsoidal/polyhedral aggregates and spherical/ellipsoidal voids in three dimensions (3D). The generation process is based on the Monte Carlo simulation method wherein the aggregates and voids are generated from prescribed distributions of their size, shape, and volume fraction. A combined numericalstatistical method is proposed to investigate damage and failure of mesoscale heterogeneous concrete: the geometrical models are first generated and meshed automatically, simulated by using cohesive zone model, and then results are statistically analysed. Zerothickness cohesive elements with different tractionseparation laws within the mortar, within the aggregates, and at the interfaces between these phases are preinserted inside solid element meshes to represent potential cracks. The proposed methodology provides an effective and efficient tool for damage and failure analysis of mesoscale heterogeneous concrete, and a comprehensive study was conducted for both 2D and 3D concrete in this paper.
1. Introduction
Concrete is the most widely used construction material in the world due to its good strength and durability relative to its cost. Numerical simulations coupled with theory and experiment are considered to be an important tool for examining the mechanical behaviour through computational materials science. Wittmann [1] proposed three levels of observation: microlevel (10^{−6} m), mesolevel (10^{−3} m), and macrolevel (10^{0} m). The microlevel represents the most basic level, in which the internal structure of cement paste is considered and where physical and chemical forces dominate. In the mesolevel, concrete is usually represented as a material made up of coarse aggregates, mortar with fine aggregates and cement paste embedded inside, and interfaces between coarse aggregates and mortar. At the macroscale, concrete is treated as a homogenous and continuous material in which the internal structure is not considered [2]. As concrete is generally used in largesized structures and its dependence of mechanical behaviours on mesostructures is significant, mesoscale modelling of concrete becomes essential and was investigated in this study.
An extensive literature review recently carried out by Wang et al. [3] shows that numerical image processing technique and parameterization modelling technique are two main approaches to generate mesostructure models of concrete. Although the mesh generated from images is closest to nature, it is presently very costly and time consuming to generate sufficient scanned images and convert them to mesh for meaningful statistical analyses (e.g., [4–6]). In the parameterization modelling technique, both direct (e.g., [3, 7–9]) and indirect approaches (e.g., [10–13]) to characterize material random heterogeneity have been used. The indirect approach is able to generate a large number of samples with ease but the key multiphasic parameters which could significantly influence the mechanical behaviour cannot be considered. The direct approach can take into account most of the key parameters such as shape, size, gradation and spatial distribution of voids and aggregates, phase volume fractions, and aggregatemortar interfaces. So among all the approaches, it seems that the direct approach is particularly suitable for statistical analysis of concrete samples and was used in this study.
Identification and generation of unit cell geometry are a vital step in the mesoscale analysis of concrete. Both shape and size of aggregates have significant influences on the stress distribution, crack initiation, and damage accumulation up to the macroscopic failure within the concrete material [16]. Most of the existing algorithms for mesostructure generation of concrete only consider random aggregates but neglect voids [2, 7, 17–19]. However, the Xray CT images [4, 14, 20] clearly show that voids exist in concrete at this scale and have severely adverse effects on the specimen strength [3, 9]. So automatic generation of morphological details of materials with randomly distributed different shape of aggregates and voids poses new challenges. In this paper, a computationally efficient and versatile inhouse program of heterogeneous material generator (HMG) has been developed considering both 2D and 3D random aggregates and voids based on the prescribed parameters.
Several numerical models for crack propagation at mesoscopic level have been used to study the heterogeneity of concrete. Continuum based finite element models are the main approaches employed in the literature [21–24], primarily based on cohesive zone model [25]. Other alternative approaches have also been developed, such as discrete element model (DEM) [26] and lattice model [11, 12]. The biggest disadvantage of DEM and lattice model is the difficulty in determining the equivalent model parameters, which is relatively straight forward for continuum based finite element models [27]. The numerical model used in this work for crack propagation at mesoscopic level is based on the cohesive zone model [3]. It is becoming more and more popular in modelling crack propagation due to its simple formation and easy implementation in the form of cohesive interface elements (CIEs).
With the preprocessing functionality in ANSYS [29] and explicit solver in ABAQUS [30], nonlinear modelling of multiple crack propagation with a considerable number of concrete samples has been performed by high performance computing in this study. Statistically, analysis of mechanical behaviour of both 2D and 3D specimens has been conducted using the combined numericalstatistical method, leading to valuable statistical results that may help improve designs of concrete materials.
2. 2D and 3D Heterogeneous Concrete
2.1. Aggregate and Void Size Distribution
The size distribution of aggregates in concrete is often described by the Fuller curve [16], which is discretised into a certain number of segments. The aggregate size distributions found in Hirsch [28] and summarised in Table 1 are used in this study. In this study, concrete is treated as a material consisting of coarse aggregates, mortar incorporating fine aggregates dissolved in it, voids, and interfaces between aggregates and mortar. Here only coarse aggregates are explicitly modelled as mesoscale features. For normal strength concrete cast with a mould in the laboratory, coarse aggregates usually represent 30–50% of the concrete volume.

From Xray CT images, it was found that gravel aggregates have circular/elliptical shapes in 2D and spherical/ellipsoidal shapes in 3D, while crushed aggregates have polygonal shapes in 2D and polyhedral shapes in 3D [14]. A number of voids can be clearly observed in concrete at mesoscale from CT images (see Figure 1). Therefore, voids with size ranging from 2 to 4 mm as presented in [31] are included in the mesostructure of concrete in this study. Coarse aggregates are considered to have elliptical/polygonal shapes in 2D and ellipsoidal/polyhedral shapes in 3D, and the shape of voids is assumed to be circle/ellipse in 2D or sphere/ellipsoid in 3D.
2.2. Specimen Generation
The basic idea is to create the aggregates and voids in the concrete in a repeated manner, until the target area/volume is filled with aggregates and voids. The generation starts with the information input process and is followed by taking and placing voids within size range; the aggregates within the grading segments are produced last. The “input” process reads the controlling parameters for generating a heterogeneous concrete; the “taking” process generates an individual void or aggregate in accordance with the random size and shape descriptions. The “placing” process positions the aggregates and voids into the predefined area/volume in a random manner, subject to the prescribed physical constraints. There are three physical conditions to be satisfied simultaneously: (1) the whole inclusion (void or aggregate) must be within the concrete, (2) there is no overlapping of inclusions, and (3) there is no intersection between any two inclusions. A threelevel hierarchical method is proposed to reduce the computational cost and is outlined in Figure 2. The procedure was programmed using MATLAB [32] in this study, and the detailed algorithms and procedure in generating mesostructures can be found in our previous work [3, 33].
2.3. Mesostructure Models
Elliptical aggregates and circular/elliptical voids in 2D or ellipsoidal aggregates and spherical/ellipsoidal voids in 3D were used for concrete with gravel aggregates. A series of concrete samples with dimensions of 50 mm × 50 mm in 2D and 50 mm × 50 mm × 50 mm in 3D are generated. Figure 3 shows the numerically generated models with gravel aggregates (aggregates in green and voids in red). Aggregate gradation in Table 1, aggregate volume fraction , void content , aspect ratios for elliptical/ellipsoidal aggregates and voids (the ratio of the major axis to the minor axis) , minimum space between the edge of an aggregate and the boundary of the concrete specimen mm, and minimum gap width between any two aggregates mm are adopted here [31].
(a) Elliptical aggregates and circular voids
(b) Elliptical aggregates and elliptical voids
(c) Ellipsoidal aggregates and spherical voids
(d) Ellipsoidal aggregates and ellipsoidal voids
Polygonal aggregates and circular/elliptical voids in 2D or polyhedral aggregates and spherical/ellipsoidal voids in 3D are considered for concrete with crushed aggregates. A series of concrete samples with the same dimensions as the previous gravel aggregate samples are generated. Figure 4 shows the numerically generated models with crushed aggregates. The minimum and maximum edges of the polygons/polyhedrons are set to , . All the other parameters are the same as the values used for gravel aggregates.
(a) Polygonal aggregates and circular voids
(b) Polygonal aggregates and elliptical voids
(c) Polyhedral aggregates and spherical voids
(d) Polyhedral aggregates and ellipsoidal voids
3. Combined NumericalStatistical Method
3.1. Description of the Method
A combined numericalstatistical method is proposed in this paper to study the material behaviour of concrete in a statistical sense. The detailed procedure is as follows:(1)Generate a model with prescribed variables, for example, sample size, aggregate volume fraction, void content, and aggregate shape.(2)Perform a finite element simulation of the sample for given boundary conditions.(3)Compute the mean value, standard deviation, and coefficient of variation (CoV) of effective property for the considered sample size.(4)Repeat steps (1) to (3) for sufficient number of random samples to meet the required precision, and conduct statistical analysis.
This procedure is automated by running a batch file in this study.
Results from all realisations are evaluated statistically. The standard deviation within a series of samples is given bywhere is the series average result and is the result from sample .
To compare results obtained with different sample sizes quantitatively, we use the coefficient of variation given byThis expresses the fluctuation of measured property relative to its average value.
3.2. Cohesive Zone Model in ABAQUS
The cohesive zone model developed by Barenblatt [34, 35] and Dugdale [36] enables the simulation of the energy dissipation process in the fracture process zone (FPZ) during fracture, as illustrated in Figure 5. A bilinear cohesive crack model is used here to predict discrete crack propagation in concrete specimens under tension loading. A linear ascending branch is added in the softening curve to model the initially uncracked material.
The separation displacement is difficult to derive from experiments, so cohesive fracture energy and cohesive strength are usually used. Among them, only two parameters are independent, and the fracture energy can be obtained aswhere is the cohesive fracture energy, is the cohesive strength, and is the separation displacement.
In the 3D cohesive zone model, it is assumed that there exist a normal traction and two tangential tractions (shear cohesion) and across the crack surfaces, through mechanisms such as material bonding, aggregate interlocking, and surface friction in the fracture process zone. Figures 6(a) and 6(b) illustrate typical linear softening curves for and , where and are the critical relative displacements when the tractions diminish for normal traction and tangential traction components, respectively. The unloading paths are also indicated.
(a) curve in normal direction
(b) curve in tangential direction
It is worth noting that the initial tensile stiffness should be high enough to ensure displacement continuity across the cohesive interfaces before the tensile strength is reached, but not too high to cause numerical instability due to illconditioned global stiffness matrix. Reasonable initial shear stiffness, and , is also needed before the shear strength, or , is reached. The effects of initial stiffness on computational results have been investigated previously [37, 38]. The following relation is suggested in [38] as a guideline for initial stiffness selection:where and are Young’s modulus and Poisson’s ratio, is the characteristic size of elements, and is taken as 10~100 from the experience in [38].
The areas under the curves in Figures 6(a) and 6(b) (calculated by (3)) represent the normal facture energy and twice the tangential fracture energy, and , respectively:
Cohesive elements in ABAQUS based on the cohesive zone model are used here. The damage is characterised by a scalar index representing the overall damage of the crack caused by all physical mechanisms. The effective relative displacements combining the effects of , , and can be obtained aswhere is the Macaulay bracket and
The damage evolution law is given bywhere is the maximum effective relative displacement attained during the loading history. and are effective relative displacements at damage initiation and complete failure, respectively. It is obvious that the damage variable evolves monotonically from 0 to 1 upon further loading after the initiation of damage.
The damage initiation and evolution will degrade the unloading and reloading stiffness coefficients , , and , which can be calculated asThe tractions are also affected by the damage according towhere and are the traction components predicted by the elastic tractiondisplacement behaviour for the current separation without damage.
In this study, damage in the cohesive zone model is assumed to initiate when a quadratic interaction function involving the nominal stress ratios reaches a value of one
The material properties, such as density, Young’s modulus, Poisson’s ratio, tensile strength, and fracture energy, are set for continuum elements of aggregates and mortar, three different interface cohesive elements. The material heterogeneity is investigated by considering different phases in the concrete specimen with corresponding material properties.
In the fracture process zone for a 2D case, tractions exist in the normal direction and shear direction across the crack interface, and the corresponding relative displacements are and . So the 2D constitutive law could be easily deduced from 3D characterization described above by taking one of the shear tractions and displacements out.
3.3. Cohesive Interface Elements Insertion
In this study, all finite element meshing is performed with the preprocessing functionality in commercial FE package ANSYS [29]. To mesh the mesoscopic structure of concrete, different material parts (mortar and aggregates) should maintain continuity on their surfaces. Hence, the finite element boundaries are coincident with material surfaces and there are no material discontinuities within the elements. The mortar and aggregates are meshed with triangular elements (plain stress) in 2D and hexahedral elements in 3D so that more realistic crack paths can be obtained. The specially developed ANSYS parametric design language (APDL) programs in combination with the ANSYS batch processing provide a powerful tool of automatic mesh generation for a large number of samples required for statistical analysis. Following the meshing procedure, the generated mesh will automatically have shared nodes at the interfaces between two elements. If the interfaces surrounding the elements are to be modelled explicitly as potential microcrack sources, a duplicate set of nodes will be required at the interface locations. Here 4node or 6node zero inplane thickness CIEs are preinserted into the existing element interfaces in 2D or 3D. This is conducted by a purpose written and flexible inhouse computer program cohesive interface elements insertion (CIEIN). Three sets of CIEs with different tractionseparation softening laws are inserted (see Figure 7), namely, CIEAGG inside the aggregates (grey in Figures 7(b) and 7(d)), CIEMOR inside the mortar (green in Figures 7(b) and 7(d)), and CIEINT on the aggregatemortar interfaces (yellow in Figures 7(b) and 7(d)). The element and node numbers are denoted by and , respectively. The detailed numbering of elements and nodes in the initial mesh and the mesh with inserted cohesive elements shows the insertion procedure with the new nodes generated at the same positions and CIEs between the continuum elements.
(a) 2D initial mesh
(b) 2D mesh with zerothickness CIEs
(c) 3D initial mesh
(d) 3D mesh with zerothickness CIEs
4. Numerical Simulation and Results
4.1. Description of the Numerical Model
2D concrete specimens with elliptical aggregates and circular voids and 3D specimens with ellipsoidal aggregates and spherical voids (, , , 0.2 mm) were modelled, loaded under uniaxial tension. 25 mm × 25 mm concrete square in 2D and 25 mm × 25 mm × 25 mm concrete cubic in 3D are used here. All the other mesoscale parameters were fixed to be the same for 2D and 3D models. Figure 8 shows the geometry, boundary, and loading conditions of the numerical model. Horizontal displacements were prescribed to all nodes on the left and right surfaces of the specimen, with a value equal to zero for the left surface, and a uniformly distributed displacement on the right surface. Vertical displacements for the same nodes are left free, except for the nodes at the left lower corner of the specimen, which are fixed to prevent rigid body translation. A displacement controlled loading scheme was used and all the analyses were ended at a displacement mm and loading time s [3]. Considering the balance between accuracy and efficiency, the element length 1 mm was used for all the meshes in this study [3].
(a) 2D
(b) 3D
Generally, aggregates have much higher strength than mortar and interfaces in normal concrete. In this section, no cracks are allowed to initiate inside the aggregates by assuming elastic behaviour without damage in CIEAGG. The linear tension/shear softening laws described above were used to model CIEs with quadratic nominal stress initiation criterion and linear damage evolution criterion. Similar material properties extracted from [39] were used in this study, which are listed in Table 2. The fracture properties related to the shear behaviour were assumed to be the same as the normal ones due to the lack of experimental data. Both initial tensile and shear stiffness for all the cohesive elements were set to 10^{6} MPa/mm after trial and error.

ABAQUS/explicit with small time increments (typically s) was utilised to solve the highly nonlinear equation systems by the high performance computing (HPC) cluster at the University of Manchester computational shared facility (CSF). A typical run of simulations with 50 samples in this study takes about 6 hours by parallel computation with 32 CPUs.
4.2. Tensile Behaviour
Due to the statistical nature of mesostructure models, an extensive series of numerical simulations would be necessary to capture the range of behaviours. Fifty random samples were modelled in this study to ensure that the results are statistically converged.
In Figure 9, resulting mean stressdisplacement and toughnessdisplacement curves, together with the variation range for both 2D and 3D, are depicted. It can be seen that the predicted stressdisplacement curves are qualitatively similar with a clear peak and sharp initial postpeak drop followed by a long shallow tail. Numerical mean stress for every curve is obtained by summing the nodal reaction forces in the constraints and dividing by the specimen cross section area. The mean peak stress predicted by 2D and 3D modelling is 2.65 MPa and 3.49 MPa, respectively, representing an increase of 24.1%. The mean toughness (energy absorbed per unit volume up to the break point, i.e., the area under the stressstrain curve up to the break point, measured in J/m^{3}) predicted in 2D and 3D modelling is 2.47 KJ/m^{3} and 3.15 KJ/m^{3}, respectively, representing an increase of 21.6%; Figure 10 shows the standard deviationdisplacement curves for both 2D and 3D modelling. The standard deviation of peak stress is 0.11 MPa and 0.05 MPa, respectively, representing a decrease of 54.5%. These increases in peak stress or decreases in standard deviation are due to the constraint effects in the thickness direction in 3D and the less smooth cracking surfaces from 3D modelling. However, the standard deviation of 3D modelling at softening range is larger than that in 2D modelling which indicates that 3D modelling is more complicated and dependent on different samples. This is due to the 3D heterogonous distribution of aggregates and voids. It is clear that the third dimension cannot be neglected due to the nature of fracture process. The importance of the thickness effect has also been pointed out experimentally by Van Mier and Schlangen [40].
(a) Stressdisplacement curves
(b) Toughnessdisplacement curves
The mean curve, mean value, and standard deviation of the stress and toughness shown in Figures 9 and 10 are extracted from the stressdisplacement curves for 50 samples with different aggregate and void distribution. In order to evaluate the effect of sample number on simulation results as required for the proposed combined numericalstatistical method, two special loading points, namely, point A at peak stress and point B at displacement 0.03mm where largest standard deviation is observed (denoted in Figures 9(a) and 10), are investigated. The influence of sample number on the coefficient of variation (CoV) of stress at points A and B is illustrated in Figure 11. It can be seen that 50 samples are enough to get convergent values of CoV of stress with a stable fluctuation. This is an indication that the combined numericalstatistical method with 50 random realisations offers a good balance between general applicability and testing of different samples. It can be observed that COV of peak stress for different void and aggregate distribution is low (4.2% for 2D modelling and 1.5% for 3D modelling), which demonstrates that the peak stress is relatively insensitive to the void and aggregate distribution. The high CoV of stress at point B further demonstrates the necessity to conduct a statistical analysis as the results vary greatly for different aggregate and void distribution in softening range.
Figure 12 shows the statistical distribution of peak stress, together with the probability density and the best fit Gaussian Probability Density Function (PDF). The probability density curve can be used to calculate structural reliability or failure reliability against given external loadings and material properties for structural design.
(a) 2D
(b) 3D
4.3. Crack Patterns
The complex mesoscale crack propagation is realistically simulated using the proposed method, and typical crack patterns for both 2D and 3D at failure are shown in Figure 13. To clearly visualise the fracture surfaces, models with damaged cohesive elements, models without cohesive elements, and failed interfaces only are used for 2D models (see Figures 13(a)–13(c)), and models with damaged cohesive elements, models without cohesive elements, and morphology of failed surface are used for 3D models (see Figures 13(d)–13(f)). The cracks shown in Figures 13(a), 13(c), and 13(d) are represented by red CIEs with damage index ( means complete failure). A magnification factor of 20 and 200 was used for the deformed models with damaged CIEs and those without damaged CIEs, respectively. The noticeable difference of 2D and 3D modelling discussed above could also be summarized in that the 3D modeling predicts more realistic fracture surfaces in the thickness direction whereas 2D modeling only assumes plane fracture problems. This result confirms the importance of including the third dimension into the analysis.
(a) 2D model with damaged CIEs
(b) 2D models without CIEs
(c) 2D models with failed interfaces
(d) 3D model with damaged CIEs
(e) 3D model without CIEs
(f) 3D morphology of failed surface
5. Conclusions
Models of numerical concrete with random mesostructures comprising elliptical/polygonal aggregates and circular/elliptical voids in 2D and ellipsoidal/polyhedral aggregates and spherical/ellipsoidal voids in 3D have been generated in this study. Numericalstatistical analysis was carried out using a cohesive zone model to simulate damage and failure of concrete. Crack patterns are realistically simulated using the technique of preembedding cohesive interface elements. The main conclusions, based on the results obtained from the statistical analysis under a uniaxial tension loading, are as follows: (1) statistical analysis is necessary in both 2D and 3D due to the high dependence of material behaviour on different aggregate and void spatial distribution; (2) the third dimension is demonstrated to have a pronounced influence on both macroscopic mechanical properties and crack patterns in tension; (3) compared to 2D modelling, 3D modelling demonstrates a larger mean peak stress and a smaller standard deviation in the prepeak response, attributed to more uniformly distributed microcracks within the specimen; a larger standard deviation/CoV of stress in the postpeak response is attributed to the larger number of possibilities for microcrack coalescence under the constraint of randomly distributed features. It has to be pointed out that the conclusions obtained are based on mesoscale modelling with specific specimen and aggregate size. The variation on resulting mechanical behaviour is also associated with the size of specimen and aggregates, and this phenomenon may not exist when specimen is large enough, for example, at the length scale of engineering structures.
Conflict of Interests
The authors declare that there is no conflict of interests regarding the publication of this paper.
Acknowledgments
This work is supported by the EPS Faculty Ph.D. Studentship from the University of Manchester for Xiaofeng Wang. The authors would like to acknowledge the assistance given by IT Services and the use of the Computational Shared Facility (CSF) at The University of Manchester.
References
 F. Wittmann, “Structure of concrete with respect to crack formation,” in Fracture Mechanics of Concrete, pp. 43–74, Elsevier, 1983. View at: Google Scholar
 A. Caballero, 3D mesomechanical numerical analysis of concrete fracture using interface elements [Ph.D. thesis], Polytechnic University of Catalonia, 2005.
 X. F. Wang, Z. J. Yang, J. R. Yates, A. P. Jivkov, and C. Zhang, “Monte Carlo simulations of mesoscale fracture modelling of concrete with random aggregates and pores,” Construction and Building Materials, vol. 75, pp. 35–45, 2015. View at: Publisher Site  Google Scholar
 H. K. Man, Analysis of 3D scale and size effects in numerical concrete [Ph.D. thesis], ETH Zurich, Zürich, Switzerland, 2010.
 R. Sharma, P. Mahajan, and R. K. Mittal, “Elastic modulus of 3D carbon/carbon composite using imagebased finite element simulations and experiments,” Composite Structures, vol. 98, pp. 69–78, 2013. View at: Publisher Site  Google Scholar
 E. Roubin, A. Vallade, N. Benkemoun, and J.B. Colliat, “Multiscale failure of heterogeneous materials: a double kinematics enhancement for Embedded Finite Element Method,” International Journal of Solids and Structures, vol. 52, pp. 180–196, 2015. View at: Publisher Site  Google Scholar
 C. M. López, I. Carol, and A. Aguado, “Mesostructural study of concrete fracture using interface elements. I: numerical model and tensile behavior,” Materials and Structures, vol. 41, pp. 583–599, 2007. View at: Google Scholar
 C. M. López, I. Carol, and A. Aguado, “Mesostructural study of concrete fracture using interface elements. II: compression, biaxial and Brazilian test,” Materials and Structures, vol. 41, no. 3, pp. 601–620, 2008. View at: Publisher Site  Google Scholar
 X. Wang, Z. Yang, and A. P. Jivkov, “Monte Carlo simulations of mesoscale fracture of concrete with random aggregates and pores: a size effect study,” Construction and Building Materials, vol. 80, pp. 262–272, 2015. View at: Publisher Site  Google Scholar
 X. F. Xu and L. GrahamBrady, “A stochastic computational method for evaluation of global and local behavior of random elastic media,” Computer Methods in Applied Mechanics and Engineering, vol. 194, no. 42–44, pp. 4362–4385, 2005. View at: Publisher Site  Google Scholar
 T. Reichert, Development of 3D lattice models for predicting nonlinear timber joint behaviour [Ph.D. thesis], Edinburgh Napier University, 2009.
 J. P. B. Leite, V. Slowik, and H. Mihashi, “Computer simulation of fracture processes of concrete using mesolevel models of lattice structures,” Cement and Concrete Research, vol. 34, no. 6, pp. 1025–1033, 2004. View at: Publisher Site  Google Scholar
 A. P. Jivkov, M. Gunther, and K. P. Travis, “Sitebond modelling of porous quasibrittle media,” Mineralogical Magazine, vol. 76, no. 8, pp. 2969–2974, 2012. View at: Publisher Site  Google Scholar
 M. L. Daudeville, Role of coarse aggregates in the triaxial behavior of concrete: experimental and numerical analysis [Ph.D. thesis], University of Grenoble, Grenoble, France, 2013.
 M. Xie, Finite element modelling of discrete crack propagation [Ph.D. thesis], University of New Mexico, 1995.
 P. Wriggers and S. O. Moftah, “Mesoscale models for concrete: homogenisation and damage behaviour,” Finite Elements in Analysis and Design, vol. 42, no. 7, pp. 623–636, 2006. View at: Publisher Site  Google Scholar
 H.K. Man and J. G. M. van Mier, “Influence of particle density on 3D size effects in the fracture of (numerical) concrete,” Mechanics of Materials, vol. 40, no. 6, pp. 470–486, 2008. View at: Publisher Site  Google Scholar
 I. M. Gitman, H. Askes, and L. J. Sluys, “Representative volume: existence and size determination,” Engineering Fracture Mechanics, vol. 74, no. 16, pp. 2518–2534, 2007. View at: Publisher Site  Google Scholar
 M. Bailakanavar, Y. Liu, J. Fish, and Y. Zheng, “Automated modeling of random inclusion composites,” Engineering with Computers, vol. 30, pp. 609–625, 2012. View at: Publisher Site  Google Scholar
 E. N. Landis and D. T. Keane, “Xray microtomography,” Materials Characterization, vol. 61, no. 12, pp. 1305–1316, 2010. View at: Publisher Site  Google Scholar
 J. F. Unger, S. Eckardt, and C. Könke, “Modelling of cohesive crack growth in concrete structures with the extended finite element method,” Computer Methods in Applied Mechanics and Engineering, vol. 196, no. 4144, pp. 4087–4100, 2007. View at: Publisher Site  Google Scholar  Zentralblatt MATH
 M. Vořechovský and V. Sadílek, “Computational modeling of size effects in concrete specimens under uniaxial tension,” International Journal of Fracture, vol. 154, no. 12, pp. 27–49, 2008. View at: Publisher Site  Google Scholar
 A. K. H. Kwan, Z. M. Wang, and H. C. Chan, “Mesoscopic study of concrete II: nonlinear finite element analysis,” Computers & Structures, vol. 70, no. 5, pp. 545–556, 1999. View at: Publisher Site  Google Scholar
 D. A. Hordijk, “Tensile and tensile fatigue behaviour of concrete; experiments, modelling and analyses,” Heron, vol. 37, no. 1, pp. 1–79, 1992. View at: Google Scholar
 A. Hillerborg, M. Modéer, and P.E. Petersson, “Analysis of crack formation and crack growth in concrete by means of fracture mechanics and finite elements,” Cement and Concrete Research, vol. 6, no. 6, pp. 773–781, 1976. View at: Publisher Site  Google Scholar
 W. Shiu, F.V. Donzé, and L. Daudeville, “Compaction process in concrete during missile impact: a DEM analysis,” Computers and Concrete, vol. 5, no. 4, pp. 329–342, 2008. View at: Publisher Site  Google Scholar
 Z. Tu and Y. Lu, “Mesoscale modelling of concrete for static and dynamic response analysis part 1: model development and implementation,” Structural Engineering and Mechanics, vol. 37, no. 2, pp. 197–213, 2011. View at: Publisher Site  Google Scholar
 T. J. Hirsch, “Modulus of elasticity of concrete affected by elastic moduli of cement paste matrix and aggregate,” ACI Journal Proceedings, vol. 59, no. 3, pp. 427–451, 1962. View at: Google Scholar
 ANSYS, ANSYS Academic Research, Release 14.0, ANSYS, 2012.
 ABAQUS 6.7, User documentation. Dessault systems, 2007.
 D. A. Hordijk, Local approach to fatigue of concrete [Ph.D. thesis], Delft University of Technology, Delft, Netherlands, 1991.
 MATLAB. R2011a, The MathWorks Incoporation, Natick, Mass, USA, 2011.
 X. Wang, M. Zhang, and A. P. Jivkov, “Computational technology for analysis of 3D mesostructure effects on damage and failure of concrete,” International Journal of Solids and Structures. In press. View at: Google Scholar
 G. I. Barenblatt, “The mathematical theory of equilibrium cracks in brittle fracture,” Advances in Applied Mechanics, vol. 7, pp. 55–129, 1962. View at: Google Scholar
 G. I. Barenblatt, “The formation of equilibrium cracks during brittle fracture. General ideas and hypotheses. Axiallysymmetric cracks,” Journal of Applied Mathematics and Mechanics, vol. 23, pp. 622–636, 1959. View at: Publisher Site  Google Scholar  Zentralblatt MATH  MathSciNet
 D. S. Dugdale, “Yielding of steel sheets containing slits,” Journal of the Mechanics and Physics of Solids, vol. 8, no. 2, pp. 100–104, 1960. View at: Publisher Site  Google Scholar
 Z. Yang, X. Wang, D. Yin, and C. Zhang, “A nonmatching finite elementscaled boundary finite element coupled method for linear elastic crack propagation modelling,” Computers & Structures, vol. 153, pp. 126–136, 2015. View at: Publisher Site  Google Scholar
 T. Qiang, X. Kou, and W. Zhou, “Threedimensional FEM interface coupled method and its application to arch dam fracture analysis,” Chinese Journal of Rock Mechanics and Engineering, vol. 19, no. 5, pp. 562–566, 2000. View at: Google Scholar
 A. Caballero, C. M. López, and I. Carol, “3D mesostructural analysis of concrete specimens under uniaxial tension,” Computer Methods in Applied Mechanics and Engineering, vol. 195, no. 52, pp. 7182–7195, 2006. View at: Publisher Site  Google Scholar
 J. Van Mier and E. Schlangen, “On the stability of softening systems,” in Fracture of Concrete and RockRecent Developments, S. P. Shah, S. E. Swartz, and B. Barr, Eds., pp. 387–396, Elsevier Applied Science, London, UK, 1989. View at: Google Scholar
Copyright
Copyright © 2015 Xiaofeng Wang and Andrey P. Jivkov. This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.