Simulating 3D crack propagation with XFEM to investigate failure mechanism in high strength concrete

High strength concrete compared to normal strength concrete has a denser micro‐structure, i.e. tight packaging or bonding between glued aggregates and cement paste, which increases strength and durability. The disproportionate ratio in ingredients of concrete, water to cement ratio and binding agents to form cement paste are subjected to self‐desiccation during the hydration process which leads to shrinkage, that has potential to form micro‐cracks and delamination within the ITZ (Interfacial Transition Zone), which ultimately causes deterioration in health of the structure. In order to investigate this phenomenon, a computational model to simulate the formation and propagation of cracks within the micro‐structure is developed. The geometrical model for the computation and simulation of the micro‐structure is extracted by means of computer tomography (CT) scans of concrete samples. The eXtended Finite Element Method (XFEM) in combination with level set techniques is employed to simulate discrete cracks and their advancement within the micro‐structure. The criterion for the crack propagation is obtained from a gradient enhanced damage model and once it is fulfilled, the crack geometry is updated using an advance algorithm for level sets.


Introduction
Modern architecture has the demand for the alteration from conventional concrete construction with bulky cross-sections, shapes and supports to more slender, possibly free-form geometrical design without compromising the purpose of build. That urges the need for concrete which provides high strength, flexibility in forming and a long-life span. High strength concrete compared to normal strength concrete has a denser micro-structure, i.e. tight packaging or bonding between glued aggregates and cement paste, which increases strength and durability. The optimized cementitious materials, water to cement (w/c) ratio and binding agent used as ingredients of high strength concrete are subjected to self-desiccation during the hydration process which leads to shrinkage. This shrinkage process is much influenced by the w/c ratio, where a higher ratio leads to the propensity for the development of micro-cracks in between the bonding zone of aggregates and cement paste region. This happens due to surface tension between the hydrated and the un-hydrated portion of mixture inside components. This zone is called ITZ -(Interfacial Transition Zones), which can be considered as the "weak-link" part of the domain, which has significant effects on the overall toughness of concrete. In contrast to normal strength concrete high strength concrete exhibits more brittle fracture behavior under static loading once it reaches the peak load because of a more homogeneous structure. However, its mechanical properties permitting its application for more slender structures, leads to a higher susceptibility to more variable loads. Thus high strength concrete is more sensitive to oscillations which are subject of fatigue failure. Fatigue failure of high strength concrete and the influence of its composition on the failure mechanisms are subject of extensive research.
To this end for the computational modeling and simulation, the micro-structural geometry is extracted using micro computer tomography (µCT) scan data and SEM technology to capture the geometry as realistically as possible. The formation of the micro-structure is discussed in more detail in section 2. The micro-structure of concrete includes many micro-cracks within the ITZ and in the bulk material. Micro-cracks may coalesce and become a macro-crack which can be modeled as a discrete crack. The 3D eXtended Finite Element Method (XFEM) in combination with level set techniques is employed here to simulate the cracks inherent and evolved within the micro-structure of high strength concrete. The need for a crack propagation criterion is met by a gradient enhanced damage model [1] which is mesh independent and allows for tracking the material degradation behavior in the domain. Beyond a certain threshold of the gradient enhanced damage, a discrete crack propagates or nucleates. This combination and the transition from continuous damage to discrete cracks provides significant advantages. Geometrically, discrete crack propagation is handled by an implicit update of level sets [2].

Experimental part
The description of the progress of the degradation of concrete under compressive fatigue loading by using damage indicators and high-resolution imaging techniques on different scales allows establishing a connection between changes in micro-and macrostructure. The damage on the macroscale is described using damage indicators, such as strain, stiffness and acoustic emission. Computer tomography-scans, as described below, allow access to the mesoscale. SEM and light microscopy of cut and polished thin sections enable the evaluation of the microstructure. Numerical models allow an in-situ perspective of degradation during load, and different parameters not accessible or too effortful for experimental characterisation can be manipulated separately. To reconstruct a numerical model close to realistic microstructure, CT scan data of the concrete structure on the mesoscale is measured to reconstruct the geometry for a simulated framework. Later on, experimental data is used to calibrate the simulation parameters.
A high-strength concrete was investigated which consists of cement (CEM I 52,5 R HS/NA), quartz sand (0/0.5 mm), sand (0/2 mm), basalt (2 mm/8 mm), PCE plasticizer, stabilizer and water. For the fatigue tests, cylindrical specimens with d/h = 60/180 mm were used. For CT scans, smaller samples are necessary which exhibit less interference with the X-rays, allowing for higher magnifications during scanning, resulting in higher resolution. To reduce specimen size to d/h = 10/30 mm, the concrete is sieved to exclude grains above 2 mm diameter. Smaller grains increase the representativeness of the grain-matrixboundaries. The specimen had a minimum age of 56d and were storage in climate conditions to conclude crystallisation processes largely. To measure the mesoscale, an Xradia 410 Versa by Zeiss® with electric potential (ϕ) = 80 kV, electric power (P ) = 10 W, a magnification of 0.4 and exposure duration of 2 sec is used. These empirical parameters show sufficient resolution while reducing artefacts and noise. For reconstruction and post-processing, software by Zeiss® is used. Since mortar displays high contrasts, a beam-hardening algorithm is applied to reduce surface artefacts. Subsequent visualisation and processing are performed with Avizo® by ThermoFisher Scientific®. To extract the micro-structure from the reconstructed CT-Scan (shown in Fig. 1, a subvolume is extracted and filtered by a non-local-means filter to denoise scalar volume data. Afterwards, an unstructured finite-element grid composed of tetrahedra is applied to the sample.

eXtended Finite Element Method
There are various methodologies developed to handle discontinuities inside the computational domain such as standard finite element with re-meshing, the cohesive zone element, the phase field approach and the general-/extenden finite element method (GFEM/XFEM) [3]. The XFEM offers a robust platform to simulate discrete cracks in the domain by means of local enrichment functions which improve the approximation and allow for strong and weak discontinuities in the primary field quantities. The XFEM has many computational advantages. One of the most important advantages is that the mesh does not need to conform with the discontinuity and the re-meshing is not required. Instead of that sets of enrichment functions are modified within the domain. That leads to an accurate reflection of kinks and jumps in the displacement field as well as an accurate prediction of the singularity of the stress field along the crack front within elements. The enriched displacement field takes the form: Here, u I is the standard degree of freedom and a J ,b jK are the enriched degrees of freedom belonging to elements completely intersected by the crack and to crack front elements, N I (X) are the standard second order Lagrangian test functions for 10noded tetrahedral elements for the set N of standard nodes, H(X) is the modified Heaviside enrichment function for the set N * of Heaviside enriched nodes, f j (X) are the front enrichment functions that provide scope to incorporate analytical terms for a more accurate approximation of highly singular region (Crack front elements) for the set N * * of crack front enriched nodes, and q(X) is a ramp function employed in the context of the corrected XFEM [4,5].

Gradient Enhanced Damage Model
The damage model accounts for the material degradation behavior and the loss of stiffness as a consequence of micro-crack formation. Discrete crack propagation demands a crack propagation criterion. The correlation between these phenomena provides a fair reason to use damage as a fracture or crack propagation criterion [6] and build up this continuous-discontinuous approach. The classical continuum damage model experience spurious mesh dependency due to localization. Non-local damage models [1] and gradient enhanced damage models [7] are developed to circumvent that problem. Here, an implicit gradient enhanced damage model formulation is proposed where c is internal length scale parameter. This Helmholtz-type equation shows an implicit dependency of the gradientenhanced equivalent strainε on the scalar valued function of the local equivalent strainε. The right hand side of Eq. (2) is obtained by choice ofε mentioned in Eq.
(2) assuming that damage occurs always under positive strains. Here, ε i are the positive eigenvalues of the strain tensor. The history parameter κ is the maximum garient-enhanced strainε ever reached in loading history. This leads to the simple yield function f (ε, κ) and the classical Kuhn-Tucker conditions.
The state of damage κ is mapped onto a scalar valued damage variable D(κ) by means of a damage evolution law.
This damage variable is used to govern material degradation. It is also used as propagation criterion for discrete crack propagation. Similar to displacement field, the gradient-enhanced strain field can experience possible jumps across crack surface. So it is enriched similar to the displacement field for better accuracy in the vicinity of crack domain.
The complete boundary value problem is given in (6).

Crack Advancement
The level set method (LSM) [8] is used to geometrically track the advancement of cracks. It is a robust tool to define moving interfaces by solving an additional Hamilton-Jacobi equation (7) for the level set ϕ i with initial and boundary values. The XFEM and LSM are combined together firstly by [9] for 2D and extended to 3D by [10][11][12].
The crack is represented by two level sets ϕ 1 and ϕ 2 , defined as signed distance functions whose gradients are chosen to be orthogonal (∇ϕ 1 · ∇ϕ 2 = 0). The intersection of the two iso-zero level sets is considered as the location of the crack front. The level set update is governed by a crack velocity field extended from the crack front to the whole level set subdomain. The process is conducted by help of (8) for the respective extension of velocities corresponding to both level sets. Here τ is a time like parameter.
In this work, the first level set is extended implicitly by means of the global crack tracking algorithm proposed in [13] by solving (9). The second level set is updated by Eq. (7).
After some steps of propagation, the level sets are likely to loose their signed distance property which is restored by a reinitialization process. Equation (10) is solved iteratively until steady state leading to the signed distance property.
Before reinitialization of the second level set, it is reorthogonalised to the first level set by solving (11) to steady state. www.gamm-proceedings.com The whole process is tabulated below to handle the update of level sets in three dimensions.

Numerical Example
The whole procedure is integrated in an algorithm to simulate the meso-structure built from the CT-scan data. A mode 1 crack propagation is simulated under dynamic load conditions. The material data is: Young's modulus E = 34.866 GPa (Cement Mixture), 109.854 GPa (Basalt), Poisson's ration 0,2 (concrete), Basalt (0.29), density ρ = −2.27 kg/m 3 (Cement Mixture) and 2.97 kg/m 3 (Basalt). It is extracted from experimental tests on high strength concrete. The model is discretized with tetrahedral elements with quadratic shape function and the simulated outcome is shown in Fig. 2.