Fracture initiation in multi-phase materials: a systematic three-dimensional approach using a FFT-based solver

Fracture initiation in multi-phase materials: a systematic three-dimensional approach using a FFT-based solver

T.W.J. de Geus M. Cottura B. Appolaire Laboratoire d’Etude des Microstructures, CNRS/Onera, BP72, 92322 Châtillon Cedex, France R.H.J. Peerlings Department of Mechanical Engineering, Eindhoven University of Technology,
P.O. Box 513, 5600 MB Eindhoven, The Netherlands
M.G.D. Geers Department of Mechanical Engineering, Eindhoven University of Technology,
P.O. Box 513, 5600 MB Eindhoven, The Netherlands
Corresponding author:  –  –

This paper studies a two-phase material with a microstructure composed of a hard brittle reinforcement phase embedded in a soft ductile matrix. It addresses the full three-dimensional nature of the microstructure and macroscopic deformation. A large ensemble of periodic microstructures is used, whereby the individual grains of the two phases are modeled using equi-sized cubes. A particular solution strategy relying on the Fast Fourier Transform is adopted, which has a high computational efficiency both in terms of speed and memory footprint, thus enabling a statistically meaningful analysis. This solution method naturally accompanies the regular microstructural model, as the Fast Fourier Transform relies on a regular grid.

Using the many considered microstructures as an ensemble, the average arrangement of phases around fracture initiation sites is objectively identified by the correlation between microstructure and fracture initiation – in three dimensions. The results show that fracture initiates where regions of the hard phase are interrupted by bands of the soft phase that are aligned with the direction of maximum shear. In such regions, the hard phase is arranged such that the area of the phase boundary perpendicular to the principal strain direction is maximum, leading to high hydrostatic tensile stresses, while not interrupting the shear bands that form in the soft phase. The local incompatibility that is present around the shear bands is responsible for a high plastic strain. By comparing the response to a two-dimensional microstructure it is observed that the response is qualitatively similar (both macroscopically and microscopically). One important difference is that the local strain partitioning between the two phases is over-predicted by the two-dimensional microstructure, leading to an overestimation of damage.

Fracture initiation in multi-phase materials: a systematic three-dimensional approach using a FFT-based solver

Corresponding author:  –  –

Keywords: micromechanics; ductile failure; damage; multi-phase materials; FFT solver

1 Introduction

Multi-phase materials are used in a wide variety of applications because they often encompass good engineering properties: electrical, thermal, mechanical, etc. By combining two or more distinct phases at the level of the microstructure, a combination of properties is obtained that cannot always be achieved otherwise. This paper focuses on the class of two-phase materials that display both strength and ductility using a microstructure in which a soft and ductile matrix is reinforced using a hard and brittle phase.

The elasto-plastic response of such materials is reasonably well understood. Simple rules of mixture yield satisfactory predictions for many purposes [1, 2], but also mean-field [3, 4] and unit cell models have been applied [5, 6, 7, 8]. Many of them have been extended to incorporate microstructural damage to predict the loss of material integrity observed at the macro-scale. However, in most cases, the role of the microstructure, and in particular of the phase distribution, is only partially considered and not fully understood.

Different researchers have performed studies that aimed at characterizing the role of the spatial distribution of the phases on the initiation and growth of fracture at the level of individual grains [9, 10, 11, 12]. In these papers it was found that damage preferentially nucleated in or in-between reinforcement particles that are neighboring along the tensile axis – see also the experimental work by Lewandowski et al. [13]. A positive hydrostatic stress was one of the configurational characteristics driving damage nucleation and growth.

Recently, De Geus et al. [12] have demonstrated that the role of the microstructure can be rationalized by considering the average phase distribution around damage sites through a conditional spatial average. The strength of this approach is that the results need little interpretation as the correlation between morphology and fracture initiation is quantified. Besides confirming the role of neighboring reinforcement particles in the tensile direction, the critical microstructure features a pattern composed of matrix phase directly adjacent to the damage site in the direction perpendicular to the tensile direction, and in bands in the directions of maximum shear. Later, it was shown that this phase distribution is also critical for fracture initiation by cleavage in the reinforcement phase [14].

Most published models are restricted to two spatial dimensions. Likewise experimental observations generally are restricted to cross-sections. In spite of the valuable insights gained, the impact of the geometrical approximation remains unclear. Some studies were performed in three dimensions [10, 11, 15], but a systematic analysis including many statistical variations is often unfeasible. Ramazani et al. [15] observed a clear difference in the macroscopic mechanical properties predicted by two- and three-dimensional micromechanical models. For a different physical process, Xu et al. [16] compared the accuracy of two-dimensional representative volume elements (RVEs) to full three-dimensional RVEs in the context of humidity induced electrochemical degradation of fuel cells with a Nickel – Yttria-stabilized Zirconia microstructure. They concluded that the yield stresses are of the same order of magnitude in the two-dimensional model as in the three-dimensional model, but the former cannot capture the thermoelectrical properties accurately as its percolation is inherently limited. These distinct mechanisms may also have a strong influence on fracture initiation for the class of materials considered here, as the plastic strain localizes in a path of soft phase that is connected in the direction of shear, prior to fracture initiation.

A systematic analysis of the local mechanical response of a three-dimensional microstructure, with statistical relevance, is still beyond the capabilities of conventional computational solid mechanics, based on the finite elements method. Different methods have been developed that enable such an analysis. Among them, the Voronoi tesselation finite elements method of Moorthy and Ghosh [17] is worth mentioning. This method employs a microstructure that is generated using a Voronoi tessellation. Each Voronoi cell forms an element for which the response is calculated using an assumed stress distribution. Although the resulting number of degrees of freedom is limited, the extension to three dimensions is still non-trivial [18].

Besides the different strategies for improving the efficiency of the calculations based on the finite elements method, full-field calculations based on the Fast Fourier Transform (FFT), as first proposed by Moulinec and Suquet [19], have proved to be a serious alternative for periodic inhomogeneous media when the number of degrees of freedom becomes significant, as it is the case in 3-D. This method takes advantage of the fact that equivalent linear elastic problems can be formulated and recast into integral forms, where the relevant Green functions are convolved with some appropriate source fields. The convolution product is expressed as a simple product in Fourier space and is computationally evaluated by means of an efficient FFT-library. Even though the application of the Fast Fourier Transform restricts the discretization to a regular grid, this method has the great advantage that the computational costs, both in terms of CPU time and memory footprint, scale approximately linearly with the number of grid points. The method can thus handle large systems and appears to be the best suited for 3-D microstructures.

The current paper aims to extend the analysis by De Geus et al. [12] to a three-dimensional setting. By comparing different three-dimensional deformation paths for a large number of three-dimensional microstructures, the role of the microstructural morphology on fracture initiation is systematically identified. The results are compared to those obtained from two-dimensional microstructures, complementing the analyses of fracture initiation mechanisms but also allowing to identify the limitations of such two-dimensional models. This is achieved by means of the algorithm proposed by Moulinec and Suquet [19] to solve mechanical equilibrium in heterogeneous materials. Besides its low memory footprint, the simplifying but reasonable assumption that the elastic moduli are homogeneous renders the method incomparably fast and the most suitable for our statistical analysis.

This paper is structured as follows. The microstructural model and the numerical implementation are discussed in Sections 2 and 3 respectively. In Section 4, the two- and three-dimensional microstructures are compared in terms of macroscopic and microscopic elasto-plastic response and fracture initiation. For the three-dimensional microstructures, the effects of different deformation paths are investigated in Section 5. The paper ends with several concluding remarks in Section 6.


second order tensor
fourth order tensor
dyadic tensor product ()
single tensor contraction ()
double tensor contraction ()
ensemble average
volume (unit cell) average

2 Model

2.1 Microstructure

The microstructure is represented by an ensemble of three-dimensional periodic unit cells. Each unit cell comprises cubes arranged in a regular grid, see for a typical example Figure 1(b). These cubes represent the individual grains or particles of the two phases considered – a soft and a hard phase with respectively low and high yield stresses. They are referred to as grains below, to avoid confusion. It must be emphasized that the model does not account for any crystalline orientation such that the boundaries between grains of the same phase are mechanically meaningless. The introduced idealization enables a clearly defined statistical analysis in which the average spatial distribution of phases around a fracture initiation site can be transparently quantified. It must be stressed that to provide a good estimate of the strain and stress fields, the numerical discretization is finer than the grid of grains, as discussed in Section 3. However, consistent with the adopted idealization, only grain averaged quantities are considered.

The phases are randomly distributed in each unit cell according to a predetermined volume fraction of the hard phase. To achieve this, each individual grain in each of the unit cells is assigned the properties of the hard phase if a random number in the interval is smaller than the chosen . All other grains are assigned the properties of the soft phase. Throughout this paper, is used. The obtained ensemble average volume fraction of the hard phase, , is within of the target value . The hard phase volume fraction in the individual unit cells varies between around this value.

Figure 1: Typical unit cell from the ensemble of unit cells representing the two-phase material (b). The numerical discretization of a single grain; each grain is discretized in the same way (a).

2.2 Constitutive model

Different constitutive models provide a satisfactory description of the micromechanical response at the considered scale and level of approximation, that of an aggregate of grains. Here, both phases are modeled using small strain elasto-visco-plasticity, whereby the parameters for the viscous part are taken to effectively approximate a rate-independent model. It has been verified that this model is sufficiently accurate for the purpose of this study by comparing the response of the selected model to rate-independent models in small and finite strain; the results of this comparison are not included here for conciseness.

In the constitutive model, the linear strain tensor is additively split in an elastic and a plastic part:


with the elastic strain tensor and the plastic strain tensor.

Assuming linear elasticity, the stress tensor is determined from the elastic strain tensor by Hooke’s law:


with the Cauchy stress tensor. is the elastic stiffness, which is assumed to be homogeneous (independent of space) and isotropic:


where the bulk modulus and the shear modulus depend on the Young’s modulus and Poisson’s ratio in the usual way: and . The tensors and respectively are the second-order identity tensor and the fourth-order symmetric identity tensor.

Plasticity is modeled using standard plasticity, by adopting the following yield function:


with Von Mises’s equivalent stress – defined as , with the stress deviator – and the yield stress. The latter depends on the accumulated equivalent plastic strain


using a linear hardening relation:


with the initial yield stress and the hardening modulus.

Associative flow is assumed to model the evolution of the plastic strain rate:


where the plastic multiplier is obtained through Norton’s law, i.e.


with the rate dependence exponent and a material constant.

The two phases have the same elastic moduli and differ only through their plastic behavior: the hard phase has a yield stress and hardening modulus which are twice those of the soft phase. The material parameters are summarized in Table 1. Note that the rate-dependence exponent of adequately approximates a rate-independent elasto-plastic model.

2.3 Applied deformation

As a result of the assumed periodicity of the microstructure in all spatial directions, only the volume averaged deformation needs to be prescribed. In this way, local (microscopic) fluctuations in stress and strain are permitted throughout the cell and along its boundary. In the following, two different average strain paths are considered: planar and axisymmetric shear. Both strain paths are isochoric and result in a vanishing macroscopic stress triaxiality. Note that the triaxiality in the individual grains is determined by equilibrium and compatibility, and generally does not vanish.

The applied macroscopic strain tensor, , for planar shear reads


in which is the macroscopic equivalent strain. Note that all directions are treated in the same way, fluctuations are thus allowed in the -direction equivalently to the - and -direction. Locally the deformation is therefore three-dimensional.

For axisymmetric shear, the strain tensor reads


where retains its meaning of the macroscopic equivalent strain.

For both cases, the microstructure is deformed up to a value of . This value is large enough to be in a relevant regime for fracture initiation, but small enough to stay within the limits of the small strain approximation.

2.4 Fracture initiation criterion

The initiation of fracture at the local grain level is modeled using a damage descriptor. This descriptor reveals fracture initiation in the grains, but does not cause any material degradation. For the soft phase the Johnson-Cook model is used, which is commonly applied to study ductile fracture initiation [20]. In this model the rate of accumulated plastic strain is compared to a critical plastic strain and then integrated over the deformation history, i.e.


The accumulated plastic strain rate, , is the time derivative of (5) and is therefore non-negative; the critical plastic strain is a function of the (local) stress triaxiality as follows:


The the stress triaxiality is the ratio of the hydrostatic stress and the equivalent stress ; the material parameters , , and are summarized in Table 1. Fracture initiation occurs in the grain once the descriptor . It is assumed that there is no fracture initiation in the hard phase.

soft 0.3 100 1 0.2 1.7 0.05
hard 0.3 100 1 - - -
Table 1: Material parameters of the constitutive model and the fracture initiation criterion [e.g. 21, 22, 23]. Note that the phases are elastically identical, and that the hard phase is assumed not to fracture.

2.5 Fracture initiation hot-spot

The relation between the microstructure and initiation of fracture is studied by calculating the characteristic phase distribution around the fracture initiation sites. To obtain a statistically representative result, all fracture initiation sites in the entire ensemble of the random unit cells are taken into account at once. The resulting average pattern indicates what is the probability to find hard and soft phase in the vicinity of any fracture initiation site. This method, as introduced in [12], is discussed in detail based on a single unit cell. The extension to the ensemble average is straightforward and therefore omitted.

A phase indicator, , is introduced to define the phase for each grain as a three-dimensional scalar field, i.e.


where corresponds to the position in the unit cell. For the cubic microstructures considered in this paper, the index is the voxel index of each grain. Similarly, fracture initiation is indicated by which is unitary if the damage descriptor , and equals zero in all other cases.

The (spatially) average(d) phase distribution around the fracture initiation sites (where ) is given by the following weighted average phase probability function:


where is the relative position with respect to the initiation fracture sites, and where in the sums loops over all grains in the unit cell, taking the periodicity into account. Note that (14) effectively quantifies the cross-correlation between and , and that it may be evaluated as a convolution (e.g. with FFT). It may be interpreted as the probability of finding the hard phase as a function of the position relative to the fracture initiation site. As the microstructure consists of only two phases, the probability of finding the soft phase is also contained in this result: if at a certain relative position (the hard phase volume fraction), then hard phase at this relative position promotes fracture initiation, while if , soft phase at that relative position will promote fracture initiation.

3 Spectral solver for mechanical equilibrium

3.1 Introduction

The numerical scheme introduced by Moulinec and Suquet [19] is used in this paper. In this scheme, mechanical equilibrium is solved iteratively whereby a reference elastic medium is used to project an intermediate solution to a new, compatible, intermediate solution until equilibrium is satisfied. This method does not require the assembly of a tangent matrix or the solution of a linear system, and therefore its memory footprint scales approximately with the number of grid-points instead of its square. The projection itself is performed using the Discrete Fourier Transform which enhances the computational efficiency. Although it has been shown that the convergence rate of this particular algorithm depends critically on the choice of the reference medium, e.g. [24], this issue is avoided in the present study because the microstructure is considered as elastically homogeneous. For the same reason, oscillations in the solution (e.g. the well-known Gibbs phenomenon) are also avoided here.

Below, we recall briefly the algorithm of Moulinec and Suquet [19], before explaining how the specific constitutive model enters this algorithm.

3.2 Solution using auxiliary problem

For this method, the local strains are the unknowns. The goal is therefore to find a field of periodic strains, , that satisfies mechanical equilibrium, i.e.


An auxiliary problem is introduced for which the stress is the response of a homogeneous linear elastic body, with stiffness , subjected to a polarization stress , i.e.


The solution to this problem is given by the following convolution, known as the periodic Lippmann-Schwinger equation:


with the periodic Green operator associated with and the average strain tensor. The convolution can be evaluated in Fourier space using a simple tensor contraction, as follows:


where the denotes that quantities have been transformed to Fourier space, and denotes the wave-vector dependence in that space.

The polarization stress is generally composed of two contributions: (i) one ensuing from the inhomogeneity of elastic constants and (ii) one involving all eigenstrains such as the plastic strain :


Since both and are not known yet, the solution has to be found iteratively. The rate of convergence is thereby determined by the choice of the reference medium (through ).

3.3 Problem specific solution

In the present work, the microstructure is elastically homogeneous, i.e. independent of , and the polarization stress becomes:


where denotes the current time increment. Hence, iterations can be entirely avoided by using an explicit time integration scheme for assessing the plastic strain at a new time step:


substituting equations (48)


Finally, knowing the polarization stress, its Fourier transform is inserted into (18) to determine the Fourier transform of the strain field at the new time increment, . The inverse Fourier transform then provides the corresponding strain field in real space, , which were the unknowns.

3.4 Applied discretization

The Fourier transform and its inverse are evaluated using the Discrete Fourier Transform, for which efficient open source implementations are readily available. Each of the grains of the microstructure is discretized using grid points (see Figure 1(a)). It has been verified that this discretization leads to a converged solution in terms of the grain averaged tensors used in this paper. The refinement study confirmed that the approximation is smaller than 1% relative error in stress and (plastic) strain measured in the individual grains.

To minimize the approximation of the explicit time integration, the deformation is applied in small steps. It has been verified that this temporal discretization is small enough to ensure stability and convergence of the scheme. Note that even though many steps are applied, the computations at each time step are fast, even for the large three-dimensional unit cells.

4 Two-dimensional vs. three-dimensional microstructure

The effect of the three-dimensional (3-D) nature of the considered microstructures is investigated by comparing them to their two-dimensional (2-D) counterpart for which each cross-section in -direction is the same. In this comparison both the 2-D and 3-D microstructures are subjected to macroscopic planar shear (see Equation 9). Since the number of considered morphologies is smaller for the 2-D microstructures, the size of the ensemble is increased to random unit cells.

4.1 Macroscopic response

The computed macroscopic equivalent stress response, , as a function of the applied equivalent strain, , is shown in Figure 2. For the 3-D microstructures, in black, and the 2-D microstructures, in blue, the upper and lower bounds of the different unit cells are plotted. The constitutive response of the individual phases is also shown using a dark blue curve for the soft matrix and red for the hard reinforcement phase. As observed, all curves coincide in the elastic regime, as the microstructure is elastically homogeneous. In the plastic regime, the response of the microstructure is a non-linear combination of the response of the individual phases. The spread between the responses of the different unit cells is larger for the 2-D model than for the 3-D model. The fact that the scatter in volume fraction between the 2-D unit-cells is somewhat larger (it is in the range ) thereby only partly explains the difference. A clear difference is also observed in the plastic response directly after yielding. The response of the 2-D model is well below that of the 3-D model. The difference is explained by examining the local plastic response in the next section.

Figure 2: The macroscopic equivalent stress, , as a function of the applied equivalent strain, . The range of the responses of the unit cells is shown for the 2-D microstructures (blue) and the 3-D microstructures (black), both for applied planar shear (Eq. 9). The constitutive responses of the hard phase (red) and the soft phase (dark blue) are shown as reference.

4.2 Microscopic response

The local accumulated plastic strain, , is shown in Figure 3 for the 2-D microstructure and in Figure 4 for the 3-D microstructure; in both cases for an applied equivalent strain . The figures include a three-dimensional view of the unit cell and planar views of the visible (outer) planes (for the 2-D microstructure only one view is relevant). Note that because of the periodicity of the microstructure, the outer planes are in fact cross-sections of the material.

For both 2-D and 3-D microstructures, it is observed that the highest values of occur in the soft phase (with an average of and respectively). In most of the hard grains is considerably lower (with an average of and respectively). From the cross-section along the -direction (Figures 3(b) and 4(c)), it is observed that the highest values of are localized in bands of connected soft phase at degree angles with respect to the -axis, in particular if these bands are surrounded by hard phase. For the 3-D microstructure, the cross-sections along the -plane and the -plane (Figure 4(b,d)) reveal that these shear bands are extended along the -axis.

Plastic strain thus localizes in shear bands at degree angles with respect to the - and -axis. This coincides with the direction of maximum shear of the applied macroscopic strain. It can be concluded that the highest plastic strain occurs where the soft phase is connected in the direction of shear, and surrounded by the hard phase.

Figure 3: Local accumulated plastic strain, , in the 2-D microstructure, resulting from an applied macroscopic planar shear strain with equivalent value . The hard grains are indicated using a black outline. The different views: (a) a three-dimensional view, and (b) the cross-section along the -plane.
Figure 4: Local accumulated plastic strain, , in the 3-D microstructure (of Figure 1); cf. Figure 3. The different views: (a) a three-dimensional view, (b–d) each of the visible cross-sections.

These observations, made on a single unit cell, can be extended to the entire ensemble. In Figure 5 the probability density, , of the accumulated plastic strain, , is shown for different strain increments. The figure includes both the 2-D (Figure 5(a–b)) and 3-D microstructures (Figure 5(c–d)); the hard phase (left column) and the soft phase (right column) are shown separately. The value of , along the vertical axes, is normalized such that the area underneath the curve is unity. The results indicate that, throughout the deformation history, the strain is more strongly partitioned between the soft and the hard phase for the 2-D microstructure than for the 3-D microstructure (i.e. larger plastic strains, mainly in the soft phase). This explains why the macroscopic hardening mostly originates from the soft phase for the 2-D microstructure. Furthermore, the plastic strain distribution is more heterogeneous for the 2-D microstructure, because any local phase arrangement that promotes plastic straining extends infinitely in the out-of-plane direction (i.e. bands of soft phase at ). In the 3-D case, the sub-surface microstructure restricts the deformation, resulting in a less heterogeneous distribution of plastic strain.

Figure 5: Probability density, , of the accumulated plastic strain, , at different strain increments. The rows correspond to: planar shear on a 2-D microstructure (a–b), planar shear on a 3-D microstructure (c–d), and axisymmetric shear on a 3-D microstructure (e–f) at different strain increments. The columns refer to hard phase (left) and to soft phase (right).

A similar observation holds for the damage descriptor in the soft phase, shown in Figure 6(a–b). For the 2-D microstructures, more grains have a higher value of . More importantly, also more grains display a value of , indicating fracture initiation. In other words, fracture initiation is ‘promoted’ in a two-dimensional microstructure, for the same reasons as observed above.

Figure 6: Probability density, , of the ductile damage descriptor, , in the soft phase at different strain increments; for (a–b) planar shear on a two- and three-dimensional microstructure respectively, and (c) axisymmetric shear on a three-dimensional microstructure.

4.3 Hot-spot for fracture initiation

The fracture initiation hot-spot is statistically identified by considering all unit cells. As discussed in Section 2.5, this leads to the probability of finding the hard phase at a certain position relative to fracture initiation: . The relative voxel coordinates are thereby aligned with the spatial coordinates . The result is shown in Figure 7 for the 2-D microstructure and in Figure 8 for the 3-D microstructure. Both of these figures include a cross-section along the -plane. For the 3-D microstructure a threshold is applied in order to visualize the most important morphological patterns. The probability scale-bar is taken symmetric about the hard phase volume fraction to avoid an optical bias. The results may be interpreted as follows: (red) indicates an elevated probability of finding the hard phase at a position relative to the fracture initiation sites and (blue) an elevated probability of finding the soft phase there.

For the 2-D microstructure, Figure 7 shows that, as assumed at the onset, fracture initiates in the soft phase which is indicated by the blue grain at . Directly next to fracture initiation in -direction hard phase is found in most cases, since . Perpendicular, in -direction, the soft phase is found with . Further away, in -direction regions of the hard phase are found on both sides. The regions of the hard phase are interrupted by regions of the soft phase under degree angles with respect to the -axis. This result is consistent with [12] even though a different constitutive model and damage indicator were used there.

For the 3-D microstructure, the cross-section along the -plane in Figure 8(a) closely resembles the result of the 2-D microstructure. The probability is higher of finding hard phase () left and right, and soft phase () above and below. Figure 8(c) shows that fracture initiates where a ‘plate’ of hard phase is interrupted by two ‘channels’ of soft phase. Soft phase under degree angles is only recovered along the in-plane direction.

Figure 7: Fracture initiation hot-spot resulting from the 2-D microstructures loaded in planar shear. The color indicates the probability of finding the hard phase, , at a position relative to fracture initiation.
Figure 8: Fracture initiation hot-spot resulting from the 3-D microstructures loaded in planar shear; cf. Figure 7. From left to right: (a) a cross-section parallel to the -plane, (b) the same cross-section with an applied threshold , and (c) the isolated threshold result.

5 Two-dimensional vs. three-dimensional deformation

The two different deformation paths from Section 2.3 are compared next, both applied to 3-D microstructures.

5.1 Macroscopic response

The macroscopic responses are shown in Figure 9. The black curves correspond to applied planar shear and the green to applied axisymmetric shear, in between the responses of the individual phases (blue for soft and red for hard). There is little scatter, below 1%, for the different microstructures in the ensemble, i.e. the upper and lower bounds of all unit cells almost coincide. At the onset of yielding, slightly more hardening is observed under axisymmetric shear than for planar shear, as the green curves lie just above the black curves for (the difference is so small that it is difficult to see).

Figure 9: The macroscopic equivalent stress response, , normalized by the initial yield stress of the soft phase, , as a function of the applied equivalent strain of all microstructures, for: planar shear (black) and axisymmetric shear (green). The bounds representing the soft (blue) and hard phase (red) are also depicted.

5.2 Microscopic plastic response

In addition to the response to planar shear, in Figure 4, the response to axisymmetric shear is shown in Figure 10. It is observed that the plastic strain response is more homogeneous for this case. Less extreme values are observed in the soft phase, and also the average of is smaller than for planar shear. At the same time the plastic strain is higher in the hard phase, where the average value is . In each of the three cross-sections, small regions of elevated occur in the soft phase, at degree angles with respect to the three axes. These localized regions are located where soft phase is surrounded by hard phase. In addition, isolated grains are observed in the soft phase where is high. These regions are part of a previously described shear band along the in-plane direction.

Figure 10: Local accumulated plastic strain, , resulting from an applied macroscopic axisymmetric shear strain with equivalent value . The microstructure and cross-sections are identical to Figure 4.

The probability density, , of the accumulated plastic strain, is shown in Figures 5(e–f) for axisymmetric shear. In comparing to the earlier results for planar shear (cf. to Figures 5(c–d)) it is observed that, throughout the deformation history, the plastic strain is more heterogeneous for planar shear than for axisymmetric shear. For the soft phase, Figure 5(d,f), a larger number of grains display a high plastic strain, and also the maximum is higher. For the hard phase, Figure 5(c,e), it is observed that the average plastic strain is higher for axisymmetric shear than for planar shear. Furthermore, the distribution is more skewed towards low values of for planar shear, indicating that most of the hard grains reveal a low .

5.3 Ductile fracture initiation

Figure 11 shows the damage descriptor, , in the microstructure of Figure 1 for an applied equivalent strain of for both load cases. The top row, Figure 11(a–c), shows three different cross-sections parallel to the -plane for the applied planar shear case. The bottom row, Figure 11(d–f), shows the same cross-sections for the axisymmetric shear case. The different cross-sections are taken at respectively one and two grains in negative -direction with respect to the reference view on the left. The hard grains are colored gray, as damage in the hard phase is not considered.

Figure 11: Damage descriptor for an applied equivalent strain of for planar shear (a–c) and axisymmetric shear (d–f). From left to right the depth of the cross-section increases by one grain. The hard phase grains have a black outline, the gray value in these grains denote that fracture initiation of the hard phase is not considered.

For axisymmetric shear it is observed that the damage, , is lower in many of the grains when compared with the planar shear case (cf. Figures 11(a,d) for a single unit cell and Figures 6(b,c) for the ensemble). This is a direct consequence of the less heterogeneous distribution of the plastic strain. The purpose of the damage descriptor is however to identify fracture initiation, i.e. the focus should be on the sites where . The number of these fracture initiation sites is approximately the same (cf. Figures 6(b,c)). A large part of the fracture initiation sites coincide for the two load cases (e.g. the grain highlighted with a white circle), but also clear differences exist. One example is highlighted with a white square, for which the damage is high for axisymmetric shear, but not for planar shear (cf. Figures 11(a,d)). For this case, the morphology consisting of a soft grain flanked by hard grains on both sides in -direction and soft grains on both sides in -direction is not extended into the subsurface (Figure 11(e)). This shows that spatial phase distribution, directly adjacent to but also further away from the fracture initiation site, in combination with the deformation, has to be critical for fracture to initiate.

5.4 Hot-spot for fracture initiation

The fracture initiation hot-spot for axisymmetric shear is shown in Figure 12, with the same views as used in Figure 8. In the case of axisymmetric shear, the cross-section in Figure 12(a) reveals a pattern that is qualitatively the same as above. Adjacent to the soft grain in the center, in -direction, hard phase is identified (), and in -direction soft phase is found in most of the cases (). The influence of the morphology next to the fracture initiation site is thereby more pronounced than for planar shear. Further away, regions of hard phase are found in -direction, which are interrupted by regions of the soft phase at degrees with respect to the -axis. The latter are however less prominent than for planar shear. Extending the observations to three dimensions in Figure 12(b–c), the regions of hard phase in -direction tend to be axisymmetric around the -axis. Regions of the soft phase are found in cones at around the -axis. Note that for all neighboring voxels at the resulting indicates soft phase, but for some it is below the threshold value.

In earlier work, De Geus et al. [12] used a simple mechanical analysis to confirm that the combination of a band of hard phase in the direction of applied tension, in combination with a phase boundary perpendicular to it, results in a hydrostatic tensile stress state. The regions of the soft phase in the directions of principal shear induce a large plastic strain. The combination of these morphological features promotes fracture initiation. This also explains the results for the three-dimensional microstructure and deformation. For planar shear, in Figure 8, the plate of hard phase corresponds to a phase boundary whose area perpendicular to the tensile axis is maximized in such a way that the soft regions at degrees with respect to the -axis are uninterrupted. The principal shear acts only in the -plane and therefore soft phase is only recovered in that plane. For axisymmetric shear, the principal shear acts on all planes at degree angles around fracture initiation, explaining the resulting morphology of soft phase. The observed pattern of hard phase again maximizes the area of the phase boundary perpendicular to the axis of applied tension avoiding interruption of the shear bands.

Figure 12: The fracture initiation hot-spot resulting form applied planar shear. The views are identical to those used in Figure 8.

6 Concluding remarks

A two-phase microstructure is studied for which the two phases have different mechanical properties (yield stresses), resulting in a sharp mechanical contrast. Consequently, fracture initiation in the soft phase is strongly influenced by the spatial distribution of the phases. A large ensemble of random three-dimensional microstructures is analyzed to systematically characterize the configuration that is most critical for fracture initiation. To this end, all fracture initiation sites are considered at once by calculating the average phase distribution around fracture initiation sites.

Fracture preferentially initiates when the following conditions apply:

  1. regions of soft phase are aligned with the shear directions;

  2. phase boundaries perpendicular to the direction of principal strain, i.e. hard phase on both sides, aligned with the tensile axis;

  3. regions of hard phase are located next to the shear bands.

Using a two-dimensional microstructure, the same conclusions are obtained. Indeed, the average phase distribution is a cross-section of the result obtained in 3-D. However, the strain is more strongly partitioned between the phases as the deformation is not restricted or promoted by the sub-surface microstructure. As a result, fracture initiation is more frequently observed in 2-D at the same overall strain.


This research was carried out under project number M22.2.11424 in the framework of the research program of the Materials innovation institute M2i (


plus 0.3ex


  • Koo et al. [1980] J.Y. Koo, M.J. Young, and G. Thomas. On the law of mixtures in dual-phase steels. Metall. Trans. A, 11(5):852–854, 1980. doi: 10.1007/BF02661217.
  • Davies [1978] R.G. Davies. The deformation behavior of a vanadium-strengthened dual phase steel. Metall. Trans. A, 9(1):41–52, 1978. doi: 10.1007/BF02647169.
  • Maire et al. [1997] E. Maire, D.S. Wilkinson, J.D. Embury, and R. Fougères. Role of damage on the flow and fracture of particulate reinforced alloys and metal matrix composites. Acta Mater., 45(12):5261–5274, 1997. doi: 10.1016/S1359-6454(97)00147-X.
  • González and LLorca [2000] C. González and J. LLorca. A self-consistent approach to the elasto-plastic behaviour of two-phase materials including damage. J. Mech. Phys. Solids, 48(4):675–692, 2000. doi: 10.1016/S0022-5096(99)00057-5.
  • Brockenbrough and Suresh [1990] J.R. Brockenbrough and S. Suresh. Plastic deformation of continuous fiber-reinforced metal-matrix composites: Effects of fiber shape and distribution. Scr. Metall. Mater., 24(2):325–330, 1990. doi: 10.1016/0956-716X(90)90264-H.
  • Brockenbrough and Zok [1995] J.R. Brockenbrough and F.W. Zok. On the role of particle cracking in flow and fracture of metal matrix composites. Acta Metall. Mater., 43(1):11–20, 1995. doi: 10.1016/0956-7151(95)90256-2.
  • Povirk [1995] G.L. Povirk. Incorporation of microstructural information into models of two-phase materials. Acta Metall. Mater., 43(8):3199–3206, 1995. doi: 10.1016/0956-7151(94)00487-3.
  • Schröder et al. [2010] J. Schröder, D. Balzani, and D. Brands. Approximation of random microstructures by periodic statistically similar representative volume elements based on lineal-path functions. Arch. Appl. Mech., 81(7):975–997, 2010. doi: 10.1007/s00419-010-0462-3.
  • Kumar et al. [2006] H. Kumar, C.L. Briant, and W.A. Curtin. Using microstructure reconstruction to model mechanical behavior in complex microstructures. Mech. Mater., 38(8-10):818–832, 2006. doi: 10.1016/j.mechmat.2005.06.030.
  • LLorca and Segurado [2004] J. LLorca and J. Segurado. Three-dimensional multiparticle cell simulations of deformation and damage in sphere-reinforced composites. Mater. Sci. Eng. A, 365(1-2):267–274, 2004. doi: 10.1016/j.msea.2003.09.035.
  • Segurado and LLorca [2006] J. Segurado and J. LLorca. Computational micromechanics of composites: The effect of particle spatial distribution. Mech. Mater., 38(8-10):873–883, 2006. doi: 10.1016/j.mechmat.2005.06.026.
  • de Geus et al. [2015] T.W.J. de Geus, R.H.J. Peerlings, and M.G.D. Geers. Microstructural topology effects on the onset of ductile failure in multi-phase materials - A systematic computational approach. Int. J. Solids Struct., 67-68:326–339, 2015. doi: 10.1016/j.ijsolstr.2015.04.035. arXiv: 1604.02858.
  • Lewandowski et al. [1989] J.J. Lewandowski, C. Liu, and W.H.J. Hunt. Effects of matrix microstructure and particle distribution on fracture of an aluminum metal matrix composite. Mater. Sci. Eng. A, 107:241–255, 1989. doi: 10.1016/0921-5093(89)90392-4.
  • de Geus et al. [2016] T.W.J. de Geus, R.H.J. Peerlings, and M.G.D. Geers. Competing damage mechanisms in a two-phase microstructure: how microstructure and loading conditions determine the onset of fracture. Int. J. Solids Struct., Accepted, 2016. doi: 10.1016/j.ijsolstr.2016.03.029. arXiv: 1603.05841.
  • Ramazani et al. [2013] A. Ramazani, P.T. Pinard, S. Richter, A. Schwedt, and U. Prahl. Characterisation of microstructure and modelling of flow behaviour of bainite-aided dual-phase steel. Comput. Mater. Sci., 2013. doi: 10.1016/j.commatsci.2013.05.017.
  • Xu et al. [2013] W. Xu, X. Sun, D. Li, S. Ryu, and M.A. Khaleel. Mechanism-based representative volume elements (RVEs) for predicting property degradations in multiphase materials. Comput. Mater. Sci., 68:152–159, 2013. doi: 10.1016/j.commatsci.2012.10.026.
  • Moorthy and Ghosh [1998] S. Moorthy and S. Ghosh. A Voronoi Cell finite element model for particle cracking in elastic-plastic composite materials. Comput. Methods Appl. Mech. Eng., 151(3-4):377–400, 1998. doi: 10.1016/S0045-7825(97)00160-6.
  • Ghosh and Moorthy [2004] S. Ghosh and S. Moorthy. Three dimensional Voronoi cell finite element model for microstructures with ellipsoidal heterogeneties. Comput. Mech., 34(6):510–531, 2004. doi: 10.1007/s00466-004-0598-5.
  • Moulinec and Suquet [1998] H. Moulinec and P. Suquet. A numerical method for computing the overall response of nonlinear composites with complex microstructure. Comput. Methods Appl. Mech. Eng., 157(1-2):69–94, 1998. doi: 10.1016/S0045-7825(97)00218-1.
  • Johnson and Cook [1985] G.R. Johnson and W.H. Cook. Fracture characteristics of three metals subjected to various strains, strain rates, temperatures and pressures. Eng. Fract. Mech., 21(1):31–48, 1985. doi: 10.1016/0013-7944(85)90052-9.
  • Sun et al. [2009] X. Sun, K.S. Choi, W.N. Liu, and M.A. Khaleel. Predicting failure modes and ductility of dual phase steels using plastic strain localization. Int. J. Plast., 25(10):1888–1909, 2009. doi: 10.1016/j.ijplas.2008.12.012.
  • Al-Abbasi and Nemes [2003] F.M. Al-Abbasi and J.A. Nemes. Micromechanical modeling of dual phase steels. Int. J. Mech. Sci., 45(9):1449–1465, 2003. doi: 10.1016/j.ijmecsci.2003.10.007.
  • Vajragupta et al. [2012] N. Vajragupta, V. Uthaisangsuk, B. Schmaling, S. Münstermann, A. Hartmaier, and W. Bleck. A micromechanical damage simulation of dual phase steels using XFEM. Comput. Mater. Sci., 54:271–279, 2012. doi: 10.1016/j.commatsci.2011.10.035.
  • Brisard and Dormieux [2010] S. Brisard and L. Dormieux. FFT-based methods for the mechanics of composites: A general variational framework. Comput. Mater. Sci., 49(3):663–671, 2010. doi: 10.1016/j.commatsci.2010.06.009.
Comments 0
Request Comment
You are adding the first comment!
How to quickly get a good reply:
  • Give credit where it’s due by listing out the positive aspects of a paper before getting into which changes should be made.
  • Be specific in your critique, and provide supporting evidence with appropriate references to substantiate general statements.
  • Your comment should inspire ideas to flow and help the author improves the paper.

The better we are at sharing our knowledge with each other, the faster we move forward.
The feedback must be of minimum 40 characters and the title a minimum of 5 characters
Add comment
Loading ...
This is a comment super asjknd jkasnjk adsnkj
The feedback must be of minumum 40 characters
The feedback must be of minumum 40 characters

You are asking your first question!
How to quickly get a good answer:
  • Keep your question short and to the point
  • Check for grammar or spelling errors.
  • Phrase it like a question
Test description