<figure><img src="image_1.png" /></figure>

HAL
open science

# Multiscale analysis of brittle failure in heterogeneous materials

Elie Eid, Rian Seghir, Julien Réthoré

<figure><img src="image_2.png" /></figure>

## To cite this version:

Elie Eid, Rian Seghir, Julien Réthoré. Multiscale analysis of brittle failure in heterogeneous materials. Journal of the Mechanics and Physics of Solids, 2021, 10.1016/j.jmps.2020.104204. hal-02978612

HAL Id: hal-02978612

https://hal.science/hal-02978612

Submitted on 26 Oct 2020

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_3.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
# Multiscale analysis of brittle failure in heterogeneous materials

EID Elie¹*, SEGHIR Rian¹, RÉTHORÉ Julien¹

GeM, Institut de Recherche en Génie Civil et Mécanique UMR 6183 CNRS

## Abstract

This paper presents a versatile model-free approach for linking the damage in highly heterogeneous materials at multiple scales. The proposed scheme evolves from phase-field modelling at the microscopic scale to simulate brittle failure, towards the estimation of the effective elastic, toughness and strength properties of the material at the mesoscopic scale, via a model-free coarse-graining technique. On one side, it's shown that, in comparison with the classical homogenisation approaches, the considered upscaling method: (i) requires no RVE (representative volume element), (ii) can be applied when the statistical homogeneity of the material ceases to exist and (iii) when sharp localisations are present. On the other, (iv) the quasi-brittle behaviour of the material is justified without any assumption on the model at the mesoscopic scale. Most prominently this paper shows that (v) the consideration of an effective homogeneous continuum to substitute a microscopically heterogeneous one is dictated by the use of a much larger regularisation parameter than what has been classically established.

<footnote>*corresponding author</footnote>

<footnote>Email addresses: elie.eid@ec-nantes.fr (EID Elie), rian.seghir@ec-nantes.fr (SEGHIR Rian), julien.rethore@ec-nantes.fr (RÉTHORÉ Julien)</footnote>
*Keywords:* heterogeneous materials, quasi-brittle failure, coarse-graining, damage modelling, multiscale

# 1. Introduction

The desire for better-performing materials has long been established, and expansion of the boundaries of the material property space has already been achieved in multiple ways, i.e., whether by manipulating the chemistry, through developing new alloys and polymers, or by manipulating the microstructure through thermomechanical processing [1]. Innovative materials for aerospace and automotive applications are recently requiring the improvement of the mechanical properties while reducing the structure's weight. Accordingly, engineers have shown interest in controlling the architecture of the materials through thoughtfully designing them in a certain fashion, to acquire improved mechanical properties over their constituents. However, the use of those highly heterogeneous materials is bridged by some limitations, de facto, treating the heterogeneities for accurately simulating the complex behaviour of such structures requires huge computational resources. Thus, of course, it's appealing to describe a simpler nature of those materials.

A variety of upscaling methods were proposed to reveal the relations between the microstructural heterogeneities from one side and the behaviour at higher scales from the other [2]. The underlying principle of existing classical homogenisation techniques lies on the description of a structure with the help of a much smaller specimen, known as the representative volume element (RVE). This implicitly assumes the presence of two separated scales: (i) the microscopic scale that is small enough to capture the heterogeneities in the
material, and (ii) the overall scale of the structure where the effects of the heterogeneities are expected to be smeared out, and on which effective material properties are considered [3]. The classical (first-order) homogenisation method is based on the construction of a boundary value problem on the RVE that allows the determination of the effective material properties at the higher scales [2]. In this case, the RVE should be big enough to statistically capture the heterogeneities and be constitutively valid, yet small enough to be considered as a volume element of continuum mechanics. More details on the classical homogenisation method can be found in the literature, and the respective work outlines the following limitations: upscaled deformation modes of an RVE for a first-order homogenisation are linear; the first-order methods cannot take into account the size effects, nor large gradients of deformation, nor localisation, i.e., also, in case of large gradients, even materials with small microstructure cannot be accurately modelled [4]. Moreover, the first-order schemes do not work for softening materials. To surpass these problems, an extension to higher-order approaches has also been addressed [4, 5]. The solution of the microscopic boundary-value problem in the case of higher-order computational homogenisation is effortless, yet allows for an enriched upscaled continuum with higher-order strain and stress fields [4]. Although the higher-order techniques are able to treat softening materials, they present their limitations: localisation bands beyond a quadratic nature for the displacements cannot be resolved, i.e., softening materials in the presence of sharp localisation regions from the presence of a crack and/or high heterogeneities [6].

Aside from determining the effective material properties, and with the
increasing demand on architecture (highly heterogeneous) materials, there's a relevant need for incorporating small-scale mechanisms of deformation and damage to essentially assess reliability and lifetime of structures within reasonable computations. As damage localises in narrow regions in a considered continuum, the length scale that determines the variation of the defect falls below the considered scale of the mechanical fields (RVE) leading thus to what is known as gradient effects [7]. Gradient theories emerging from the multiscale nature of the mechanical framework are based on the enrichment of the classical continuum description with additional terms; those allow taking those gradient effects into account [8]. When the constitutive equations at the higher scales are difficult to write, general methods based on concurrent finite element simulations (FE2) can be applied [9]. FE2 methods do not require any constitutive equations because all non-linearities come directly from the homogenisation of microscopic quantities after applying localisation rules to determine local solutions. Interests are presently concentrating on the development of a continuous-discontinuous homogenisation scheme, to allow the assessment of the presence of both micro and macro cracks, and where localisation bands are incorporated at the macroscale [10].

Recently, work has been done on deriving a homogenised cohesive law at the macroscale from computations of crack propagation in a microscopic sample [11]. In [12], X-FEM approaches are used to incorporate the discontinuity at the macroscale, but as previous methods, this technique relies heavily on the principle of separation of scales as well as the presence of an RVE and it remains a homogenisation technique where a small part of the domain is considered to extract the full response of the structure; plus,
both those methods rely on concurrent multilevel finite element (FE2) which
is computationally expensive, and requires the difficult task of writing a
consistent homogenisation scheme to link the scales [9]. Moreover, the need
for enrichment of the description of the multiscale problem is directly per-
ceived [13, 14, 15]. Although those methods are theoretically prominent, their
experimental applicability remains questionable. In [16], a different approach
has been proposed, where the effective toughness of the heterogeneous media
was directly evaluated a priori (without concurrent computations). Recently,
[17] followed the work of [16] to identify the different parameters of a damage
model at the mesoscale by fitting a typical force-displacement response on a
heterogeneous structure. Yet, as the effective material properties are deter-
mined macroscopically from force-displacement responses, it's believed that
micro-cracks and their influence on the structural responses fail to be taken
into consideration.

The above mentioned methods stand as long as the separation of scales is
prominent, or as long as an RVE can be well-defined, which naturally leads
to a homogeneous description of the microstructures at the macroscopic scale.
Nonetheless, when the microstructure's heterogeneities and/or the damage
distribution are not statistically homogeneous, the effective description is
expected to be dependent on the position in space and becomes influenced
by simultaneous interactions between the damage and the microstructure.
This restrains the above mentioned methods from accurately transferring
the information from the microscopic scale to the macroscopic scale. Plus,
to the author's best knowledge, there's no available formulation allowing a
consistent transfer of such information between the scales.
Accordingly, we believe that a consistent micro-meso analysis on the dam-
age of highly heterogeneous materials' (architected materials between them)
requires proceeding through multiple intermediate scales to ensure the proper
modelling between the scales. By setting the intermediate mesoscales, different
continuum descriptions are sought going from the microscale consistently all
the way to the larger macroscale. For this purpose, we follow the bottom-up
approach of analysis in which information at the microscale is considered
to inform the larger scales. We use a model-free coarse-graining technique
[18] - a widely used technique in molecular dynamics studies - where pseudo-
molecular systems are built to reproduce physically consistent behaviour of
all-atom with easier and faster computations. The coarse-graining technique
[18] is adapted to evaluate continuum mechanics at different intermediate
mesoscales solely from the gathered data at the microscopic scale and by
a manipulation of the inviolable conservation laws. The effective fields of
mechanical properties are thus established. Ultimately, we seek to find the
conditions under which one may describe the crack propagation process in het-
erogeneous materials by substituting it by an effective homogeneous medium
with effective elastic moduli, fracture toughness and strength. The proposed
scheme relies on simulating (i) the failure process of architected materials
(modelled as periodic and quasi-periodic microstructures) at the microscopic
scale by explicitly taking the heterogeneities into account. The acquired infor-
mation is then (ii) upscaled to mesoscopic scales by the means of the proposed
model-free coarse-graining technique; before finally (iii) analysing the obtained
information at the considered scale. This method allows the construction
of consistent density, displacement, strain and stress fields at larger scales
<figure><img src="image_4.png" /><figcaption>Figure 1: The scales of interest for the failure problem of microscopically heterogeneous materials: the microscopic scale (a) the intermediate mesoscopic scales with properties obtained upon the length scales considered (b) and the structural macroscopic scale (c).</figcaption></figure>

based on the actual physics in question at the scale of the heterogeneities. Without any a priori on the material’s behaviour, the herein proposed scheme provides a genuine evaluation of the effective material and failure properties at the considered scales. The different scales of interest are presented in Section 2.1. The paper is organised as follows: the computational method considered for the simulations at the microscopic scales is briefly introduced in Section 2.2 and the proposed coarse-graining method is provided in Section 2.3. Application of the scheme for the analysis of mesoscopic damage on typical microstructures is then constructed and discussed in Section 3; the effective material and failure properties are assessed in a lower-cost and more straightforward manner.

## 2. Computational Approach

### 2.1. *General Statements*

Let us briefly recall the multiple scales of interest at which damage problems can be tackled: the microscopic level is the level at which the
material's architecture is prominent. In this study, reproducible lattices of
pores are considered, e.g., Figure 1(a). Continuum mechanics apply; brittle
failure occurs and can be simulated by the linear elastic fracture mechanics
or its approximations (Phase-field Modelling [19], Eigen erosion [20]...). At
the macroscale, the material is seen as a homogeneous bulk (Figure 1(c)),
where linear elastic fracture mechanics elements can be implemented, and
microstructural features give rise to resistance-curve behaviour. Effective
macroscopic physical parameters shall be determined. In between those scales
are the mesoscales. At the mesoscales, the literature suggests considering a
homogeneous material (Figure 1(b)), and the effects of the microstructural
heterogeneities are implemented in the modelling with the introduction of a
process zone of size related to a certain length parameter to which damage
spreads. Phase-Field Modelling, Eigen Erosion and Thick Level Set applied
to quasi-brittle failure are examples of models considered at this scale. The
homogeneity of the fields at the intermediate mesoscales is put into question
in this study and the bottom-up approach is considered via a model-free
upscaling technique for the analysis to answer the arisen questions.

## 2.2. Phase-field modelling

In this section, the method considered for the micromechanical simulation
of crack propagation is presented. With the expected complex crack networks,
techniques like X-FEM [21] that require the predefinition of a crack and/or
cohesive element methods [22] that only allow separations on mesh's boundary
are discarded, and the variational approach [19] known as phase-field method
is considered for building the micro-mechanical numerical experiments. The
robustness and versatility of the approach regarding its independence of the
finite element mesh were decisive in the choice of the phase-field model for the micro-mechanical simulations. Of course, more recently, a Thick Level-Set (TLS) [23] approach was developed. In the TLS approach, the damage zone is separated from the undamaged zone and the damage variable is linked to the level set function which itself is a parameter of the model. But with fewer assumptions on the model, the phase-field framework is believed to be more convenient for our study.

In this section, the phase-field approach as implemented in [24] is briefly presented; the choice regarding the method's parameters is based on the study led in [25].

Assuming small strains, the phase-field introduces the following energy functional for a cracked body in a regularized framework:

$$E(u, \Gamma(\alpha)) = E_u(u, \Gamma(\alpha)) + E_s(\Gamma(\alpha)) \quad (1)$$

Where $u$ and $\Gamma$ are the variables representing the displacement and the crack surface respectively. The crack surface being a function of a continuous damage variable $\alpha$. $\alpha$ describes the material damage state: it takes the value 0 in the intact region of the material and $0 < \alpha \le 1$ to represent the crack. $E(u, \Gamma(\alpha))$ is the strain energy stored in the cracked body, $E_s(\Gamma(\alpha))$ is the energy required to create the crack according to Griffith Criterion - known as the fracture energy.

The sharp crack is smeared-out by a regularisation parameter $l_c$ [26], and the fracture energy is thus written as a function of the regularised crack density
function $\gamma(\alpha, \nabla\alpha)$.

$$E(u, \alpha) = E_u(u, \Gamma(\alpha)) + E_s(\Gamma(\alpha)) = \int_{\Omega} W_u(\varepsilon(u), \alpha) d\Omega + g_c \int_{\Omega} \gamma(\alpha, \nabla\alpha) d\Omega \quad (2)$$

It's recalled that this approach uses the continuous field $\alpha$ to describe the discontinuities coming from the presence of a crack by the means of a crack density function $\gamma$. The term $W_u(\varepsilon(u), \alpha)$ represents the strain energy density in the cracked body, $\varepsilon$ is the displacement symmetric gradient and $g_c$ is the fracture toughness. The regularisation parameter $l_c$ choice has also been previously assessed ([27, 28]). And it was shown that this phase-field modelling approach converges to the classical brittle failure when the regularisation parameter $l_c$ approaches 0. More details on the choice of the parameters considered for the finite element implementation are given in Section 3.1

## 2.3. Coarse-graining

When seeking to find continuum mechanics for the fracture process, out of information gathered at the microscopic scale, the challenge is to find an appropriate technique that takes into consideration the real solicitation of the material as well as the singularities coming from crack propagation and the presence of the heterogeneities. The micromechanical fields coming from phase-field simulations are upscaled by adapting the method from [29] and [18]: a physically consistent upscaling coarse-graining method that allows going from discrete probability density into an upscaled continuum (Figure 2).

The parameters of this method are:

* the convolution or the coarse-graining function $\phi$ that can be any
sufficiently regular function with a local support: a variety of forms were studied and similar results were obtained; in this paper, the normalised Gaussian distribution (Figure 2(b)) of zero mean $\mu$ and a standard deviation $\sigma$ is considered. It takes the following form:

$$ \phi(x, \sigma, \mu) = \frac{1}{\sigma\sqrt{2\pi}}e^{\frac{-(x-\mu)^2}{2(\sigma)^2}} \quad (3) $$

* the width of the convolution function $l_{CG} = w/2$. In 3, $w = 2 \times 3\sigma$: it's the most important parameter that defines the different length-scales at which the problem is inspected. The normalised Gaussian function can be rewritten as follows:

$$ \phi_{l_{CG}}(x, l_{CG}) = \frac{1}{\frac{l_{CG}}{3}\sqrt{2\pi}}e^{\frac{-x^2}{2(\frac{l_{CG}}{3})^2}} \quad (4) $$

* the discretization $H$ considered for the coarser mesh (the support for the coarse-graining): this parameter defines the resolution of the coarse-grained fields (Figure 2(c)) and does not affect their distribution. Identical results are obtained from investigations on multiple discretization sizes validating thus the mesh objectivity of the proposed upscaling technique.

In [18], a system of particles indexed by $e$ is considered (Figure 2(a)), with known masses $m_e(t)$ and centres of masses $\underline{r}_e(t)$ at time $t$. The coarse-grained mass density at position $r$ and time $t$ is given by:

$$ \rho(\underline{r}, t) = \sum_{i} m_e \phi[\underline{r} - \underline{r}_e(t)] \quad (5) $$
<figure><img src="image_5.png" /><figcaption>Figure 2: The coarse-graining function (b) sweeps over the different points in the domain (c). Information from the particles system (a) is smoothed and continuous fields are computed.</figcaption></figure>

Unlike in [18], continuum data is considered at the fine-scale from the micromechanical simulations. Let $\Omega^0$ be the domain of interest in the microstructure, a discretization of $\Omega^0$ into finite elements serves as a support for the coarse-graining computations. The coarse-grained mass density $R(x, t)$ at position $x$ in $\Omega^0$, at time $t$, is defined as the convolution between the microscopic density function $\rho$ and the predefined coarse-graining function $\phi$:

$$R(\underline{x}, t) = \int_{\Omega^0} \rho(\underline{x} - \underline{x}', t)\phi(\underline{x}', t)d\underline{x}' \quad (6)$$

For the sake of simplicity, the following notation is considered to replace the convolution:

$$\langle \rho(\underline{x}, t) \rangle_{\phi} = \int_{\Omega^0} \rho(\underline{x} - \underline{x}', t)\phi(\underline{x}', t)d\underline{x}' \quad (7)$$

and the coarse-grained mass density $R(x, t)$, at position $x$ and time $t$, would be as follows:

$$R(x, t) = \langle \rho(x, t) \rangle_{\phi} \quad (8)$$

From this spatial/temporal definition of the coarse-grained mass density, and by imposing the mechanical conservation laws at both the microscopic and
coarse-grained scale, expressions for the impulsions, velocities, displacements and stresses are obtained at different positions $\underline{x}$ and times $t$, at the coarse-grained scale. We start by recalling the conservation laws written at the microscopic scale; $i$ and $j$ denote the different directions in the considered space:

• Balance of Mass

$$ \frac{\partial \rho}{\partial t} + \frac{\partial (\rho v_i)}{\partial x_i} = 0 \qquad (9) $$

• Balance of Momentum

$$ \frac{\partial}{\partial t} \rho v_i + \frac{\partial}{\partial x_i} \rho v_i v_j = \frac{\partial}{\partial x_j} \sigma_{ij} \qquad (10) $$

### 2.3.1. Balance of mass

A simple manipulation of (9) allows the computation of the expression for the velocity at the coarse-grained scale. Computing the convolution of both sides of the equation, one can obtain:

$$ \langle \frac{\partial \rho}{\partial t} \rangle_{\phi} = - \langle \frac{\partial (\rho v_i)}{\partial x_i} \rangle_{\phi} \qquad (11) $$

The left side of the equation denotes the time derivative of the coarse-grained density $R(x, t)$:

$$ \frac{\partial R}{\partial t} = - \langle \frac{\partial}{\partial x_i} \rho v_i \rangle_{\phi} \qquad (12) $$

Using the basic rule of the derivation of convolution, one can write:

$$ \frac{\partial R}{\partial t} + \frac{\partial}{\partial x_i} \langle \rho v_i \rangle_{\phi} = 0 \qquad (13) $$
Writing the balance of mass at the coarse-grained scale, with $RV_i$ denoting the impulsion $P_i$ at the coarse-grained scale:

$$ \frac{\partial R}{\partial t} + \frac{\partial}{\partial x_i} RV_i = 0 \qquad (14) $$

and by identification between (2.3.1) and (2.3.1), we can conclude that

$$ RV_i = \langle \rho v_i \rangle_{\phi} \qquad (15) $$

Identifying the coarse-grained impulsion, $P_i = RV_i$, and the microscopic impulsion, $p_i = \rho v_i$, one can see that the coarse-grained impulsion is equal to the coarse-graining of the microscopic impulsion, which is not the case for the velocity field. The velocity at the coarse-grained scale is the ratio between the upscaled impulsion and the coarse-graining mass density:

$$ V_i = \frac{\langle \rho v_i \rangle_{\phi}}{R} = \frac{\langle p_i \rangle_{\phi}}{R} = \frac{P_i}{R} \qquad (16) $$

In this study, continuum mechanics is assumed to hold in all length scales involved; the derivation is restricted to small displacement gradients and the discussion is confined to two- dimensional quasi-static problems on perfectly solid materials. Therefore, the coarse-grained displacement $U_i$ and velocity $V_i$ fields have similar expressions, from (2.3.1):

$$ U_i = \frac{\langle \rho u_i \rangle_{\phi}}{R} \qquad (17) $$

Next, it is natural to proceed with a strain calculation based on the coarse-grained displacements:

$$ \underline{E} = \frac{1}{2}(\underline{\nabla}U + \underline{\nabla}^T U) \qquad (18) $$
### 2.3.2. Balance of linear momentum

At the microscopic scale, the balance of linear momentum states:

$$ \frac{\partial}{\partial t} \rho v_i + \frac{\partial}{\partial x_i} \rho v_i v_j = \frac{\partial}{\partial x_i} \sigma_{ij} \quad (19) $$

at the mesoscopic scale, a similar expression is expected with coarse-grained mechanical fields, to be written as:

$$ \frac{\partial}{\partial t} R V_i + \frac{\partial}{\partial x_i} R V_i V_j = \frac{\partial}{\partial x_i} S_{ij} \quad (20) $$

From the time derivative of the coarse-grained impulsion $P_i = R V_i = \langle \rho v_i \rangle_\phi$ (2.3.1), and using the basic rule of derivation, the expression of the stresses at the coarse-grained scale is determined:

$$ \frac{\partial P_i}{\partial t} = \langle \frac{\partial}{\partial t} \rho v_i \rangle_\phi \quad (21) $$

from the balance of momentum at the microscopic scale (10), we can write (2.3.2) as:

$$ \frac{\partial P_i}{\partial t} = \langle \frac{\partial}{\partial x_i} (\sigma_{ij} - \rho v_i v_j) \rangle_\phi \quad (22) $$

It's here interesting to introduce what is called 'fluctuating velocity' $v'_i = v_i - V_i$. This velocity does not add any impulsion to the system, and the coarse-grained fluctuation impulsion vanishes as:

$$ \langle \rho v_i' \rangle_\phi = \int \rho (v_i(\underline{x} - \underline{x}', t) - V_i(\underline{x}, t)) \phi(\underline{x}) d\underline{x}' = \langle \rho v_i \rangle_\phi - V_i \langle \rho \rangle_\phi = P_i - R V_i = 0 \quad (23) $$

Once, $v_i$ is replaced in (2.3.2) by $v'_i + V_i$, the following equation can be written:

$$ \frac{\partial P_i}{\partial t} + \frac{\partial}{\partial x_i} R V_i V_j = \frac{\partial}{\partial x_i} \langle \sigma_{ij} - \rho v_i' v_j' \rangle_\phi \quad (24) $$
Now writing the coarse-grained balance of linear momentum as a function of the coarse-grained variables, and by identification between (2.3.2) and (2.3.2) the expression of the stress at the coarse-grained scale is obtained:

$$S_{ij} = \langle \sigma_{ij} - \rho v'_{i} v'_{j} \rangle_{\phi} \qquad (25)$$

In quasi-statics, the dynamic terms $v'_{i}v'_{j}$ will be neglected and the stress at the coarse-grained stress field $\underline{S}$ scale is found to be equivalent to the convolution of the microscopic stress field with the coarse-graining function: $\underline{S} = \langle \underline{\sigma} \rangle_{\phi}$.

Finally, from equations (2.3.1), (2.3.1) and (2.3.2), displacement, strain and stress fields are constructed out of micromechanical simulations. As seen, this upscaling technique requires no condition on the geometrical aspect of the microstructure, nor on the micromechanical fields, nor puts any a priori on the behaviour at the mesoscopic scales. It's indeed applicable on arbitrary heterogeneous (e.g., non-periodic) materials even when sharp localisation is present.

A general scheme of the analysis based on this upscaling technique is illustrated by the following steps:
1. Generating the geometries with suggested microstructures.

2. Acquiring the microscopic displacement, stress and strain fields before failure.

3. Acquiring the microscopic displacement, stress and strain fields from the micro-mechanical simulations of failure.

4. Coarse-graining the mechanical responses obtained from step 2 and 3 with different $l_{CG}$.

5. Analysing the constructed database on heterogeneous microstructures and determining the effective elastic properties ($\kappa, \mu, E, \nu$) from steps 2 and 4, the failure toughness $G_d$ and strength $\sigma_f$ properties from steps 3 and 4.

### 3. Application on Mesoscale Analysis of damage

In this section, we present detailed information about an application of the proposed scheme where we present the typical microstructures of interest and their corresponding micromechanical simulation followed by an in-depth investigation of the results.

#### 3.1. Typical Microstructures

Due to their wide presence and uses, periodic and quasi-periodic microstructures are considered in the study, with the idea that any material would behave between a perfectly periodic one and a quasi-periodic material presenting long-range heterogeneities. For this purpose, hexagonal and kite& dart Penrose paving are studied. The materials' symmetry order (6 and 5-fold symmetry respectively) should lead to elastic isotropic equivalent
<table><caption>Table 1: Typical microstructures and their properties studied: The periodic microstructure presents a unique characteristic length (the distance between the holes); while the quasi-periodic shows two [30]</caption><thead><tr><th colspan="3">Typical microstructures properties</th></tr><tr><th>Geometry</th><th>Characteristic length(s) (μm)</th><th>Hole radius (μm)</th></tr></thead><tbody><tr><td>Periodic</td><td>3000</td><td>750</td></tr><tr><td>Quasi-Periodic - Type 1</td><td>2270 and 3670</td><td>750</td></tr><tr><td>Quasi-Periodic - Type 2</td><td>2270 and 3670</td><td>750</td></tr></tbody></table>

<figure><img src="image_6.png" /><figcaption>Figure 3: Typical microstructures considered for the simulations of periodic material (a) and two types of quasi-periodic microstructures: holes distribution at the nodal positions of a kite& dart paving (Type 1) (b) and on the centroids of the kite& dart paving (Type 2) (c).</figcaption></figure>
<table><caption>Table 2: Material properties considered in the phase-field simulation</caption><tr><td>ρ</td><td>1200kg/m³</td></tr><tr><td>E</td><td>3GPa</td></tr><tr><td>ν</td><td>0.35</td></tr><tr><td>g<sub>c</sub></td><td>250J/m²</td></tr><tr><td>l<sub>c</sub></td><td>400μm</td></tr></table>

media - a property to be discussed in the upcoming sections. We recall that for the hexagonal distribution, there's only one characteristic length of the microstructure that is the distance between holes, denoted by $d$; while for the Penrose kite& dart paving, there are two characteristic lengths. As mentioned previously, the microstructure is modelled by the drilled holes inside a bulk material which mechanical properties are presented in Table 2. $ρ$ is the material density, $E$ and $ν$ are respectively the Young Modulus and Poisson ratio, $l_c$ is the internal length of the phase-field model and $g_c$ the fracture toughness. In order to meticulously compare the microstructures, the same hole radii $r_h = 750μm$ are taken, and the mean distances between the holes is fixed to $d = 3000μm$ corresponding to a volume fraction of 78% for the periodic microstructure and 75% for the quasi-periodic ones. The generated microstructures are shown in Figure 3, and their respective geometrical aspects are displayed in Table 1. Two types of quasi-periodic microstructures are considered, both based on the kite& dart Penrose paving. Type 1 corresponds to the holes drilled at the nodes of the paving, while type 2 suggests drilling holes at the centroids of the kites and darts in the paving, leading thus to same-yet-different perfectly controlled microstructures each
<figure><img src="image_7.png" /><figcaption>Figure 4: Full (a) and close-up (b) diagrams of a typical finite element mesh of a TDCB specimen with the considered coarse-graining support mesh region (in light blue) allowing the study of large scales $l_{CG}$ up to 10 times the mean distance between holes $d$.</figcaption></figure>

one presenting specific quasi-periodic patterns. For the micro-mechanical simulations, the microstructure is put at the core of a Tapered Double Cantilevered Beam (TDCB) fracture geometry to provide crack growth stability from the tapered profile of the specimen [31, 32]; the microstructure is surrounded by a homogeneous bulk material, the dimensions are put forth in Figure 4. Displacement boundary conditions are applied for the phase-field micro-mechanical simulation. In fact, the stability in crack growth provided by the TDCB specimen is believed to match the stability provided by the application of a surfing boundary condition [16]. In both cases, the crack evolves as it pleases inside the microstructure. In this study, the former - more straightforward - approach is adopted. We recall that we confine
ourselves to the quasi-static evolution of cracks neglecting thus dynamics effects. The numerical discretization $h = 200\mu m$ is thoughtfully adapted to the heterogeneities sizes and placements inside the domain and based on [25], the internal parameter of the phase-field model $l_c$ is set to $400\mu m$. Both lengths are much smaller than the structure's heterogeneities leading to mesh-independent crack initiation and propagation.

The influence of the microstructures on the crack propagation is prominent and put forth in Figure 3. For the periodic material, a simple linear crack path is obtained suggesting thus the presence of weak planes [30]. While for the quasi-periodic materials, each type proposes different-more tortuous paths suggesting more complexity of the damaging process. The influence of the microstructure is visibly reflected on both paths. A detailed crack path analysis is found in Section 3.4. After acquiring the simulation results, a domain of interest $\Omega^0$ is considered for the study and then subdivided into finite constant-stress-constant-strain elements of size $H = 1mm$ without loss of generality (Figure 4(b)). The considered coarse-graining function then sweeps over the mesh support constructing a database of mechanically and physically consistent fields at multiple larger scales to be analysed without any a priori on the model at each scale and by considering the physics happening at the scale of heterogeneities. $\Omega^0$ is at the core of the microstructure which allows the use of large $l_{CG}$ and therefore large observation scales.

Exploiting the adapted coarse-graining method, one is able to investigate the microstructure at different transitional scales by building density, displacement, strain and stress fields at each scale.
<figure><img src="image_8.png" /><figcaption>Figure 5: Effective density fields of type 1 quasi-periodic microstructure for three different scales.</figcaption></figure>

## 3.2. Density

As established in the method, the density field $R(x)$ is first computed. Effective density fields of the quasi-periodic type 1 microstructure at different coarse-graining scales $l_{CG}$ are presented in Figure 5. Results actually show the ability of the method to construct heterogeneous/homogeneous continuous density fields depending on the scale of interest without any a priori on the effective field's homogeneity. To analyse the "homogenization" of the effective density as the coarse-graining scale $l_{CG}$ is increased, we plot the evolution of the mean and the coefficient of variation of $R$ for the three considered microstructures (Figure 6 (a)). It's considered that the homogeneity of a field is attained once its corresponding coefficient of variation $COV$ drops below the $\le 1\%$ threshold. The arithmetic mean is calculated as the sum of the sampled values (whether at the nodes or the Gaussian points of the coarse elements in $\Omega^0$) divided by the total number of samples (nodes and Gaussian points respectively). The standard deviation is found by taking the square root of the average of the squared differences of the values from their mean value. The coefficient of variation is defined as the ratio of the standard
<figure><img src="image_9.png" /><figcaption>Figure 6: Density conservation through mesoscales (a) and the evolution of the corresponding coefficient of variation $COV_R$ - defining the heterogeneity of the effective density - evolution with $l_{CG}$ (b)</figcaption></figure>

deviation to the arithmetic mean.

As the implemented coarse-graining method is based on the inviolable conservation laws - mass continuity between them -, it can be seen that the mean effective density in the studied domain is conserved through the scales and only the homogeneity of the field is altered. The heterogeneities of the effective density of the material are smeared-out much faster when considering a periodic microstructure at $l_{CG}/d = 1$, while the quasi-periodic microstructures require higher coarse-graining scales for the density heterogeneities to smear-out at $l_{CG}/d = 4$ (Figure 6 (b)).

### 3.3. *Elastic Properties*

Once the density fields are computed, manipulating the balance of mass at the fine and the coarser scale leads to the computation of the effective displacement fields that can be differentiated to determine the strain fields. The balance of linear momentum allows the evaluation of the effective stress
<figure><img src="image_10.png" /><figcaption>Figure 7: Fields of the $C_{11}$ component of the effective elasticity tensor of type 1 quasi-periodic microstructure at three different scales.</figcaption></figure>

First, we aim to determine the behaviour of the material prior to damping. For this purpose, we put the microstructure at the core of a rectangular specimen on which two tensile and one shear test simulations are conducted. From the microscopic mechanical fields, coarse-grained displacements, strains and stresses can be evaluated for different $l_{CG}$. Coarse-grained strain and stress couples in $\Omega_0$ (obtained from the three tests) provide an evaluation of the nine components $C_{ij}$ of the effective stiffness tensor $C$, representative of the elastic behaviour, at each material point for the considered scales. It was observed that in the specimen coordinate system, the shear-extension coupling terms vanish; the reduced expression for the effective elasticity tensor written in (26) is thus adopted:

$$ \begin{bmatrix} S_{11} \\ S_{22} \\ S_{12} \end{bmatrix} = \begin{bmatrix} C_{11} & C_{12} & 0 \\ C_{21} & C_{22} & 0 \\ 0 & 0 & C_{66} \end{bmatrix} \begin{bmatrix} E_{11} \\ E_{22} \\ 2E_{12} \end{bmatrix} \quad (26) $$

Fields of the $C_{11}$ component of the effective elasticity tensor of the quasi-periodic type 1 microstructure at three different coarse-graining scales $l_{CG}$ are
presented in Figure 7; the heterogeneity of $C_{11}$ is shown to persist for large regularisation scales and thus the influence of the distribution of holes on the effective stiffness fields. From the computed effective elasticity tensors, both the anisotropy and homogeneity can be evaluated. As mentioned previously, the symmetry order of the studied periodic and quasi-periodic microstructures (6 and 5-fold symmetry respectively) are expected to lead to an equivalent isotropic response. From here, we aim to determine the scale from which the symmetry orders actually governs the elastic isotropy. To do so, the two-dimensional elastic anisotropy index $a_r$ - defined in [33] - is analysed at each material point for each observation scale. An explicit expression of $a_r$ as a function of $C_{ij}$ and $S_{Cij}$ can be written as follows:

$$ a_r = \left( \left[ \frac{1}{4}(C_{11} + C_{22} + 2C_{12})(S_{C11} + S_{C22} + 2S_{C12}) - 1 \right]^2 + 2\left[ \frac{1}{16}(C_{11} + C_{22} - 2C_{12} + 4C_{66})(S_{C11} + S_{C22} - 2S_{C12} + S_{C66}) - 1 \right]^2 \right)^{\frac{1}{2}} \quad (27) $$

Where $S_{Cij}$ represent the components of the compliance tensor $S_C$ defined as $S_C = C^{-1}$. $a_r$ takes the value of 0 in the case of perfect isotropy. Otherwise, $a_r$ increases as the anisotropy strengthens. The main advantage of $a_r$ over other anisotropy indices (Kube [34], Zener [35],...) is its direct applicability in 2D, for any symmetry type, and its direct evaluation from the elasticity tensor. In this study, $a_r$ is computed from the evaluated elasticity tensors $C$ and $S_C$ at each material point of the domain for each coarse-graining scale $l_{CG}$. The distribution of these indices is shown in Figure 8. It's observed that for the periodic material, all the points present an isotropic response
<figure><img src="image_11.png" /><figcaption>Figure 8: Material points distribution based on the elastic anisotropy index $a_r$ computed at different coarse-graining scales and for the considered microstructures.</figcaption></figure>

(with an error lower than 1%) once $l_{CG}$ is larger than the characteristic length $d$ of the microstructure. Although the material symmetry order of the quasi-periodic materials considered suggests elastic isotropic behaviour, it's quite clear that this isotropy of the elasticity tensor cannot be reached at small scales but in fact, suggests a required scale of observation of at least 3 times the characteristic length $d$ of the microstructure where the anisotropy indices of all the material points fall within 1% to 0 corresponding to the isotropy of the elasticity tensor in the whole domain.
Next, the homogeneity of the elastic moduli is put forth at different $l_{CG}$. Only the observation scales at which the effective elasticity in the whole domain is isotropic (within 1% error) are considered in the study (Figure 8), i.e., $l_{CG}/d \ge 1$ for the periodic microstructure and $l_{CG}/d \ge 3$ for the two quasi-periodic microstructures. From the isotropic elasticity tensors fields, effective Young modulus (E) and Poisson ratio ($\nu$) are computed for different $l_{CG}$ on each element in the domain $\Omega_0$. Statistical analysis is conducted on the evaluated effective elastic moduli to identify their possible "homogenisation" for each microstructure. The mean and standard deviation of the two material properties are plotted in Figures 9 and 10 as well as the evolution of the coefficient of variation of each property with the coarse-graining scale. The Poisson ratio converges and homogenises at considerably low length scales of observation. For the periodic microstructure, the effective Young modulus is homogenised at considerably small observation scales ($l_{CG}/d = 1.5$). The average effective Young modulus for the different microstructures rapidly converges, contrary to its variation, especially for the quasi-periodic microstructures where the effective Young modulus field requires larger length scales, at least 7 times the characteristic length $d$, to homogenise. Next, we proceed with the multi-scale damage analysis, starting with the influence of the observation scale on the crack path towards the effect on the strength and toughness fields as the crack advances inside the microstructures.

### 3.4. Crack Path

Displacement and stress fields from the phase-field simulations of fracture on the TDCB are upscaled, and coarse-grained displacements, strains and stresses with damage consideration are obtained. As the crack advances in
<figure><img src="image_12.png" /><figcaption>Figure 9: The effective Young modulus computed at different $l_{CG}$ for the three microstructures (a) and the evolution of the corresponding coefficient of variation $COV_E$ - defining the heterogeneity of the effective Young modulus - with $l_{CG}$ (b)</figcaption></figure>

<figure><img src="image_13.png" /><figcaption>Figure 10: The effective Poisson ratio computed at different $l_{CG}$ for the three microstructures (a) and the evolution of the corresponding coefficient of variation $COV_ν$ - defining the heterogeneity of the Poisson ratio - with $l_{CG}$ (b)</figcaption></figure>
the microstructure, interactions between the crack and the structure of the
material are observed, especially in the quasi-periodic microstructures. The
crack can be expected to follow the path that would allow the maximum dis-
sipation of energy. For type 1 quasi-periodic distribution, "resilient patterns"
(red circles in Figure 11(b)) impose the deviation of the crack. For the type-2
quasi-periodic distribution, kinking of the crack is present, and due to the high
amount of elastic energy stored in the specimen before kinking (red circles in
Figure 11(c)), failure becomes unstable. The question arising here regards
the ability to replace the complex crack path (at the microscopic scale) by an
equivalent failure band (at larger scales of observations) represented by an
effective straight crack path. To answer that question, the effective crack path
at the coarse-grained scale is determined. Without loss of generality, we define
the effective crack tip at time *t* as the local zone where the maximum stress
occurs at this time *t*. By considering different criteria (Rankine, Maximum
Volumetric Stress, ...) to evaluate the position of maximum stress, identical
results were obtained; that is of course due to the stress singularity at the
crack tip. The choice of the Rankine criterion naturally emerges for the
description of brittle and quasi-brittle failure. From here, the effective "meso-
scopic" crack tip position is defined as the local zone where the maximum
Rankine stress is reached for each microscopic crack tip position. The length
scale introduced via the coarse-graining *l*<sub>CG</sub> method suggests softening of the
material (Figure 14 - more details about the stress-strain response is found in
Section 3.5)- and an equivalency to a process zone is present suggesting that
the critical stress before softening happens at a distance to the true crack
tip (the crack tip determined from the micromechanical simulations). The
<figure><img src="image_14.png" /><figcaption>Figure 11: Mesoscopic crack paths for the type 1 (b) and type 2 (c) microstructures at different $l_{CG}$ and the corresponding crack tortuosity evolution with $l_{CG}$ (a).</figcaption></figure>

evolution of the crack paths and their corresponding tortuosity at different scales in the quasi-periodic microstructures are shown in Figure 11. Counter intuitively, one can observe prominent tortuosity of the crack path even at large scales. In fact, for the type 1 quasi-periodic microstructure, the crack tortuosity only drops of less than 6% at $l_{CG}/d = 10$ while the type 2 crack path tortuosity drops around 9% from its microscopic value (Figure 11(a)). The conservation of the tortuosity across the considered scales drops the idea of the consideration of an effective straight failure band to replace complex crack paths. FFT analysis on the crack paths at different observation scales offers an insight on the amplitudes and wavelengths that contribute to the effective crack deflection. Later, the interaction between the wavelengths driving the crack deflection and the wavelengths driving the variations of the critical stresses and energy dissipation along the crack path are confronted. Figure 12 shows the FFT analysis results on the crack path inside the two quasi-periodic microstructures. As expected, small wavelengths $\lambda/d \le 3$ are smeared-out as the coarse-graining scale increases and more weight on the larger wavelengths is observed $\lambda/d \ge 10$. The type 1 response shows decaying
<figure><img src="image_15.png" /><figcaption>Figure 12: FFT analysis on the mesoscopic crack path inside the quasi-periodic type 1 (a), and quasi-periodic type 2 (b) microstructures</figcaption></figure>

amplitudes at wavelengths $\lambda/d = 2-3$ with a drop of about 62% (Figure 12(a)). The contribution of the wavelengths around $\lambda/d = 5$ corresponding to the distance between zones of "resilient patterns" (Figure 11(b)) is conserved through the scales and is responsible for 80% of the crack deflection. For the type 2 quasi-periodic microstructure, one clearly observe the absence of uniquely conserved high amplitude wavelengths across the scales (Figure 12(b)), except for $\lambda/d = 2-3$ that actually corresponds to the kinking spots. The crack path, in comparison with the type 1 -except for the kinking spots- does not conserve the same wavelengths suggesting a more easily smoothed crack path. Next, we focus on the influence of the microscopic heterogeneity on the resistance and toughness fields at different coarse-graining scales, and we confront the wavelengths controlling the heterogeneities with the wavelengths present in the crack paths.
<figure><img src="image_16.png" /><figcaption>Figure 13: Upscaled vertical displacement $U_y$ at different $l_{CG}$. The discontinuity of the displacement field is regularized by the coarse-graining function and dictated by its width</figcaption></figure>

3.5. *Fracture strength*

Displacement and stress fields from the phase-field simulations of fracture on the TDCB specimens with the microstructures at their cores are upscaled, and coarse-grained displacements, strains and stresses with damage consideration are obtained. Figure 13 shows the vertical displacement field at different coarse-graining scales. The jump on the crack faces are smoothed as the observation scale increases and the sharpness of the crack at these scales is expected to be smeared out. A plot of the stress-strain relation computed at the coarse-grained scales of a TDCB test corresponds to a stress-strain response of a quasi-brittle behaviour (Figure 14). Here, the history of the hydrostatic stresses as a function of hydrostatic strains is plotted for different points in $\Omega_0$.

As shown in Figure 14, the material undergoes a linear elastic trend followed by a non-linear region before reaching the critical stress where softening occurs. For small length scales, the increased width of the pack in the linear elastic region suggests heterogeneity of the modulus, to confirm the previous results regarding the homogenization of the elastic properties
<figure><img src="image_17.png" /><figcaption>Figure 14: Stress-strain behaviour of elements along the crack path and in a neighbourhood of the crack at two length scales $l_{CG}$</figcaption></figure>

of the materials. Non-uniqueness of the stress-strain response is observed when studying the response in a neighbourhood around the crack. The curve reaching the maximum stress corresponds to the response of the elements along the crack path. As the distance to the crack path $d_c$ increases, the stress-strain response of the corresponding elements reach lower maximum stress states followed by some softening. When the elements are much further to the crack path, they remain undamaged and show typical linear elastic stress-strain response. For instance Figure 14 shows the stress-strain response of a material in the neighbourhood of a crack at a specific abscissa in the domain, at $l_{CG}/d = 1$ which corresponds to $l_{CG} = d = 3mm$. As the coarse-graining mesh size $H$ is equal to $1mm$, we find 6 different stress-strain responses of the elements corresponding to the discretization of the damageable zone. The same trend is found at different coarse-graining scales $l_{CG}$ (Figure 14(a)) where a unique relation between the local stress and local strain states does not exist. Instead, as we move further from the crack path, and due to the regularizing nature of the coarse-graining, softening persists yet starts at lower critical stress levels leading thus to non-unique relations between the
stresses and strains. Non-locality of the fields is thus probable.

Moreover, when following the stress-strain response of the elements along the crack path, especially at low scales, and as the crack propagates through different patterns in the microstructure, the critical stress state reached before softening is found to differ from one position to another along the path. We define the maximum Rankine stress reached at each position of the crack tip as the *fracture strength* or simply *strength* denoted $\sigma_f$, without abuse of language. The evolution of the critical Rankine stress at different scales of observation is plotted in Figure 15 for the periodic (a), quasi-periodic type 1 (b) and quasi-periodic type 2 (c). The below figures thus show the significant influence of the local differences inside the microstructure - distribution of holes - on the fracture strength of the effective continuum. As the crack gets trapped inside the holes, much higher loading is required for the crack re-initiation, this phenomenon is mainly observed on the periodic geometry at the smallest scales, while for the quasi-periodic geometries, not only crack trapping influence the critical stress state reached, but also the crack deflection and deviation around the special "resilient patterns".

The strength fields at different coarse-graining scales $l_{CG}$ are studied and both the average critical stress and its coefficient of variation over the crack path for each length scale $l_{CG}$ are plotted in Figure 16. It's clear that the critical stress $\sigma_f$ decreases when $l_{CG}$ increases. Moreover, a tendency stating $\sigma_f \sim \frac{1}{\sqrt{l_{CG}}}$ is found. Similar expressions relating the tensile strength $\sigma_f$ to the characteristic length $l_c$, E and $g_c$ are found in [25, 36, 37, 38] for gradient and non-local damage models. Kindly refer to Section 3.6 for more details about this finding.
<figure><img src="image_18.png" /><figcaption>Figure 15: Fracture strength $\sigma_f$ evolution along the crack path for the periodic (a), quasi-periodic type 1 (b) and quasi-periodic type 2 (c) microstructures.</figcaption></figure>
In Figure 16, one observes higher strength for the periodic microstructure in comparison with the two quasi-periodic microstructures. A study on the evolution of the "homogeneity" of the strength is also conducted. Here, the coefficient of variation is evaluated by considering the squared differences of the values from the normalised values of a homogeneous response - allowing the study of the deviation of the strength surpassing the geometrical influences of the TDCB geometry and the loading conditions. Unlike the elasticity, the effective strength of the microstructures remains highly heterogeneous even at large $l_{CG}$. From Figure 16(b), it's clear that the coefficient of variation stabilizes near slightly increases at large scales of observations independently of the microstructure. The coefficient of variation of the strength of type 1 and the type 2 microstructures stagnates at about 2% and 3% respectively. The periodic's $COV_{\sigma_f}$ converges to 0.6%. Counter-intuitively, the coefficient of variation of the effective strength of the type 1 microstructure is lower than that of type 2 even though the path is more complex, but this may be caused by the presence of kinking in specific places leading to a huge increase of the loadings before the crack saps in comparison to the rest of domain where the crack path is smooth and straight.

Looking at Figure 15, it's hard to quantify both the microstructural effects and the influence of $l_{CG}$ on the fracture strength evolution in the material. For this purpose, the fracture evolution is studied in the frequency domain, and FFT analysis allows to display the wavelengths and amplitudes to better depict the interactions of the microstructures and the coarse-graining impacting the strength. Figure 17 compares the wavelength spectrum of the FFT analysis transformed from the fracture strength's deviation $D_{\sigma_f}$. The
<figure><img src="image_19.png" /><figcaption>Figure 16: Evolution of the mean fracture strength $\sigma_f$ as a function of $l_{CG}$ for the three microstructures (a) and the evolution of the corresponding coefficient of variation $COV_{\sigma_f}$ defining the heterogeneity of the effective strength field of the continuum (b)</figcaption></figure>

amplitude is computed as the ratio of the amplitude of each wavelength to the RMS (root mean square) amplitude of the input signal. The periodic microstructure analysis shows a dominant peak occurring at a wavelength $\lambda = d$ corresponding to the length scale of the microstructure. The relatively low (30%) drop of the contribution from the small scale to the largest considered scale of $\lambda = d$ suggests that the variation of the strength $\sigma_f$ of the material is controlled mainly by the distance $d$ between the holes even on large observation scales for a periodic microstructure. This observation could explain the stagnation of $COV_{\sigma_f}$ across the scales. As seen previously on the crack path analysis, more weight is put on the larger wavelengths as the coarse-graining scales increase. The size of the holes does not present any influence on the strength variations. Regarding the type 1 quasi-periodic microstructure in Figure 17(b), one clearly observes the drop of the contributions of the small length scales and the rise of the contributions of the larger wavelengths as
the coarse-graining scales increase. The peaks on the wavelengths around $\lambda/d = 2 - 3$ persist with increased influence as the scale of observation enlarges. Comparing the FFT analysis of the crack path and the strength offers insights on the link between the wavelengths controlling the crack inside the microstructure and the effective strength. As the crack deflects mainly every $\lambda/d = 5$, it would have travelled two "resilient patterns" zones. This reflects the periodicity of the crack path that is twice the periodicity of the mechanical response. Once again, we can see that the main variation of the strength of the material is directly controlled by the distribution of holes. For the type 2 quasi-periodic microstructure (Figure 17(c)), even at the smallest scales, the heterogeneity of the strength comes from the long-range variations. At the smallest scale, 50% of the deviation comes from the larger wavelengths, which is due to the kinking that happen for the crack inside this microstructure at wide distances in the domain. Moreover, the small peak present at the small scale at $\lambda/d = 1$ corresponds to the crack travelling from one hole to another in a straight path between the kinking zones. The other small peak at $\lambda/d = 2 - 3$ (also met in type 1 microstructure) decays as the scale is increased. As the coarse-graining scale enlarges, the small wavelengths contributions are smeared out leading to the extreme rise of the influence of the large scales onto the strength variations.

### 3.6. *Fracture toughness*

The difference in the loading history of the material points inside the microstructure - whether distributed along their crack path or in the neighbourhood of a crack (Figure 14) - is directly related to the total amount of energy absorbed by the material points until fracture. Here, the focus is
<figure><img src="image_20.png" /><figcaption>Figure 17: FFT analysis on the fracture strength $\sigma_f$ for the periodic microstructure (a), quasi-periodic type 1 (b) and quasi-periodic type 2 (c) microstructures.</figcaption></figure>
on the energy absorbed by the material points along the crack path. The effective toughness $G_d$ is defined as the energy to total failure evaluated as:

$$G_d = \int_{t_0}^{t_f} \underline{S} : \dot{\underline{E}} dt \qquad (28)$$

allowing thus the measure of the energy dissipated per unit volume from the start of the loading ($t = t_0$) until the fracture of the specimen ($t = t_f$). $\underline{S}$ is the mesoscopic stress tensor and $\dot{\underline{E}}$ is the mesoscopic strain rate tensor. The evolution of $G_d$ at different scales of observation is plotted in figure 18 for the periodic (a), quasi-periodic type 1 (b) and quasi-periodic type 2 (c). The below figures thus show the significant influence of the local differences inside the microstructure -holes distribution- on the overall dissipated energy along the crack path of the obtained continuum. As the crack gets trapped inside the holes, much higher energy is required for the crack re-initiation, this phenomenon is mainly observed on the periodic geometry at the smallest scales, while for the quasi-periodic geometries, not only crack trapping influence the energy dissipation, but also the crack deflection and deviation around the special "resilient patterns".

A similar analysis to the one in Section 3.5 is conducted: a study of the evolution of this toughness parameter $G_d$ at different coarse-graining scales followed by an FFT analysis to better understand the relationship between the microstructure and the effective toughness.

The average toughness $G_d$ is inversely proportional to the coarse-graining scale $l_{CG}$ (Figure 19(a)), and we find $G_d \sim \frac{1}{l_{CG}}$. Indeed, coarse-graining admits that the displacement, stress and strain fields on a material point actually depend on the state variables distribution in a neighbourhood of the point under consideration. The size of the neighbourhood is depicted by
<figure><img src="image_21.png" /><figcaption>Figure 18: Fracture toughness $G_d$ evolution along the crack path for the periodic (a), quasi-periodic type 1 (b) and quasi-periodic type 2 (c) microstructures.</figcaption></figure>
$l_{CG}$. Here, we find $G_d \sim \frac{1}{l_{CG}}$ and $\sigma_f \sim \frac{1}{\sqrt{l_{CG}}}$. As previously stated, similar expressions relating the tensile strength $\sigma_f$ to the characteristic length $l_c$, $g_c$ and $E$ can be found in [25, 36, 37, 38]. From [25], the critical value of the tensile strength in uniaxial traction is given by:

$$ \sigma_f = \frac{9}{16} \sqrt{\frac{E g_c}{3 l_c}} \qquad (29) $$

This suggests the following: for a fixed toughness $g_c$, the same relation between $\sigma_f$ and $l_{CG}$ holds between $\sigma_f$ and $l_c$: $\sigma_f \sim \frac{1}{\sqrt{l_c}}$. Ultimately, for an experimentally determined $\sigma_f$, the relationship $G_d \sim \frac{g_c}{l_c}$ stands. We recall that the relations found in the literature hold only for uniaxial traction. We mention that no consensus on the relation linking $E$, $\sigma_f$, $l_c$ and $g_c$ taking into account the different loading conditions, and/or specimen geometries can be found in the literature. To emphasize, we recall that all the results presented in this paper are found without any a priori on the material behaviour at the coarse-grained scales.

A study on the evolution of the "homogeneity" of the effective toughness is conducted. In Figure 19, we plot the evolution of the effective fracture toughness $G_d$ as a function of $l_{CG}$ for the three microstructures (a) and the evolution of the corresponding coefficient of variation $COV_{G_d}$ defining the heterogeneity of the effective toughness field (b).

As compared to the strength $\sigma_f$, one notices that the heterogeneity of the toughness field $G_d$ is higher than the one of the strength field at the same scales. Yet no stability of the coefficient of variation of this quantity is observed at the large scales. As long as the crack path is straight, both the strength and the toughness evolve in the same manner at different coarse-graining scale. Regarding the quasi-periodic microstructures, one observes
inversion of "homogeneity" of the toughness in comparison with the strength,
i.e., when comparing the plots in both Figure (a) and Figure (b), the blue
and green curves corresponding to the coefficients of variation of σ<sub>f</sub> and G<sub>d</sub>
are inverted; COV<sub>σ<sub>f</sub></sub> is greater for the quasi-periodic type 2 microstructure
as compared to quasi-periodic type 1, while COV<sub>G<sub>d</sub></sub> is smaller. In fact, the
crack might deflect to maximize its energy dissipation. This raises questions
on the drift from LEFM where the notions of critical stress intensity factors
and the critical energy release rate are somehow two faces of the same coin.
Again, FFT analysis transformed from the effective toughness deviation D<sub>G<sub>d</sub></sub>
is presented in Figure 20.

For all three microstructures, small wavelengths λ/d ≤ 3 are smeared-
out as the coarse-graining scale increases and more weight on the larger
wavelengths is observed λ/d ≥ 10. In comparison with the strength responses,
we can see that the larger wavelengths contributions to the heterogeneity
of the toughness fields significantly increase as l<sub>CG</sub> increases (a rise of more
than 100%) for all the microstructures). The peak on λ/d = 2 − 3 observed
at the smallest scale for the type 1 microstructure (in both the strength
and the toughness responses) is decreasing as the coarse-graining increases
to give way to the larger wavelengths (see Figure 20). As compared to
the FFT analysis of the fracture strength, the dominant peaks observed at
smaller scales are no longer present through the observation scales and this
suggests the following: as the regularization via coarse-graining increases, the
influence of the microstructure on G<sub>d</sub> is smeared-out and the local phenomena
intervening in the energy dissipation process are smoothed, which can lead to
an actual homogenisation of this parameter in comparison with the strength
<figure><img src="image_22.png" /><figcaption>Figure 19: Evolution of the effective fracture toughness $G_d$ as a function of $l_{CG}$ for the three microstructures (a) and the evolution of the corresponding coefficient of variation $COV_{G_d}$ defining the heterogeneity of the effective fracture toughness field of the continuum (b)</figcaption></figure>

$\sigma_f$ where the microstructural influence persists even for larger coarse-graining scales.

4. Discussion and Concluding Remarks

This paper proposes a model-free coarse-graining method that does not require specific boundary conditions (as opposed to classical homogenisation schemes) and is indeed applicable to non-periodic structures presenting high strain localisation. Additionally, a model-free analysis of failure at the mesoscopic scale is performed. Phase-field simulation of failure on different microstructures was considered for the micromechanical simulations. The obtained data are upscaled at different mesoscopic scales $l_{CG}$ via the proposed coarse-graining technique that is solely based on the conservation laws. Density, displacement, stresses and strain fields at the mesoscopic scales are constructed and analysis on the established database is led.
<figure><img src="image_23.png" /><figcaption>Figure 20: FFT analysis on the fracture toughness $G_d$ for the periodic microstructure (a), quasi-periodic type 1 (b) and quasi-periodic type 2 (c) microstructures.</figcaption></figure>
**Quasi-brittle behaviour** The coarse-graining method introduces a length scale which implies softening of the material - and thus an equivalency to a softening behaviour with a process zone, without any a priori on the behaviour at the larger scales. Plotting the stresses against the strains shows a typical response of quasi-brittle materials where a linear elastic region is followed by a non-linear region before softening. The notion of strength is thus notable.

**Non-local effects** Without any assumption on the material behaviour, the absence of a unique behaviour law that links the local variables, i.e., local strains and local stresses, is illustrated. The stress-strain history of the elements at different distances to the crack path is different suggesting thus non-locality of the behaviour.

Once the behaviour at larger scales was established from the microscopic data, failure analysis was led on different heterogeneous materials, and for this purpose, three microstructures were considered: the periodic hexagonal distribution of holes and two types of kite& dart Penrose distribution of holes: at the nodal positions (type 1) and on the centroids of the kites and darts (type 2). A study on the effective material properties was led followed by a failure analysis to further understand interactions between the cracks, the scales and the microstructures.

**Density** The equivalent densities at different coarse-graining scales are conserved - an expected results since the essence of the coarse-graining method is the conservation principles. A homogeneous density field ($COV < 1\%$) is obtained at $l_{CG}/d = 1$ for the periodic and $l_{CG}/d = 4$ for the quasi-periodic microstructures.
**Elastic properties** From the micromechanical simulations before failure, coarse-grained elastic stresses and strains allowed the computation of elasticity fields. Although the anisotropy of a microstructure is depicted by its symmetry order, the proposed scheme shows that a relatively large observation scale is needed for the symmetry order to reveal its influence. From the scale $l_{CG}/d = 1$, isotropy converges for the periodic microstructure while the isotropy of the quasi-periodic structure requires $l_{CG}/d = 3$. The heterogeneities found at the microscopic scale influence the equivalent elastic properties at the mesoscopic scale. We show that in order to consider a homogeneous isotropic elastic equivalent medium, the required $l_{CG}$ exceeds the values considered in the literature and in fact is much larger when considering non-periodic microstructures with long-range heterogeneities. For the periodic distribution, an equivalent isotropic homogeneous elastic medium is found starting $l_{CG}/d = 1.5$, while for both quasi-periodic distributions, an $l_{CG}/d = 7$ is required to obtain an equivalent elastic homogeneous medium.

**Crack path** Phase-field simulations on the considered microstructures showed highly complex crack networks especially for the quasi-periodic microstructures. The effective crack path was analysed at different observation scales. We define the crack tip at each time step as the position at which the equivalent Rankine stress is maximum. The evolution of the crack paths shows prominent tortuosity across the observation scales. FFT analysis on the crack paths offered insight on the amplitudes and wavelengths that contribute to the effective crack paths. For the type 1 quasi-periodic microstructure, peaks at similar wavelengths for different coarse-graining scales persist suggesting thus that the microstructural distribution of holes underlines the conserva-
tion of the tortuosity and prohibits the consideration of a smooth effective
failure band. For the type 2 microstructure, one clearly observes the absence
of uniquely conserved high amplitude wavelengths across the scales which
suggests a more readily smoothed crack path between the kinking spots.

**Fracture strength** The notion of strength emerging from the obtained coarse-grained stress-strain response was analysed. We define the fracture strength as the critical equivalent Rankine stress state reached at each effective crack tip position. The effective fracture strength is found to vary from one position of the crack tip to another along the crack path. This is explained by different phenomena including the trapping, re-nucleation, and deflection of the cracks after advancing inside the microstructure. This property is found to be the hardest to smear-out; in fact, the influence of the microstructure persists in all three microstructures even for relatively large coarse-graining scales. The stagnation of the coefficient of variation of the strength for the considered microstructures suggests the inability of the consideration of a homogeneous strength field for an effective medium that takes into account the heterogeneities at the smaller scales and that has a dominant role in stress concentration, crack initiation and general propagation. To better understand the microstructure's influence coupled with the scales, *FFT* of the crack paths and the fracture strengths are confronted. Even at large coarse-graining scales, i.e., where smoothing of local phenomena is expected, the strength remains highly influenced by the relatively small scales of the microstructures and large wavelengths do not have a significant influence on this quantity. This, with the evolution of the coefficient of variation of *σ*<sub>f</sub> leads to believe in the impossibility to completely smear-out the heterogeneities involved in
the fracture strength field of locally heterogeneous material.

**Fracture Toughness** We also proposed a definition of the fracture toughness. Following the evolution of the coefficient of variation for the different microstructures, coarse-graining shows good ability to smear-out the microstructural effects on the equivalent toughness with no stagnation of the coefficient of variation at the large scales. FFT analysis of the fracture toughness evolution in the domain shows the large wavelengths contributions suggesting more easily smoothed parameter in comparison with the fracture strength, whereas the small wavelength amplitudes related to the microstructures are smeared-out.

This leads us to believe that considering a homogeneous model for simulating quasi-brittle failure in highly heterogeneous materials requires the consideration of extremely large scales to always be on the safe side when considering equivalent isotropic homogeneous properties. We also show the inevitability of the consideration of a non-homogeneous material in which the influence of substructures is preserved at the mesoscopic scales.

**Open issues and extensions** Within the proposed framework, we are able to perform a micro-meso analysis of quasi-brittle failure on periodic and quasi-periodic materials. The examples presented herein are limited to quasi-static crack propagation with imposed displacement. The absence of branching and multi cracks limited the analysis to a single crack propagation. The consideration of the phase-field model at the microscopic scale to build the ground for the multi-scale analysis presented some limitations, especially regarding the model parameter and the calculation times. Yet, once the micromechanical data are obtained, the proposed coarse-graining method can
be easily and efficiently implemented and is readily paralleled. The arbitrary
mesoscopic discretization choice *H* did not influence on the overall results
and only alters the resolution of the coarse-grained fields. The maximum
observation scale *l*<sub>CG</sub> is bounded by the considered TDCB geometry. Con-
fronting macroscopic load-displacement, with the obtained mesoscopic results
at the larger scales and projecting the obtained displacement fields on the
Williams' series [39] should give more insight on the multiscale behaviour. We
are currently exploiting the framework for a meso-macro analysis. Moreover,
we are working on a multi-scale dynamic analysis of crack propagation inside
the considered microstructures.

Acknowledgements

We acknowledge the financial support of Pays de La Loire, the city of
Nantes and Ecole Centrale de Nantes through grant ConnectTalent IDS. The
support of the French Research National Agency programme through grant
ANR-16-CE30-0007-01 is also acknowledged.

References

[1] N. A. Fleck, V. S. Deshpande, M. F. Ashby, Micro-architectured materials: past, present and future, Proceedings of the Royal Society A: Mathematical, Physical and Engineering Sciences 466 (2121) (2010) 2495-2516. doi:10.1098/rspa.2010.0215.

URL https://royalsocietypublishing.org/doi/10.1098/rspa.2010.0215
[2] F. Ongaro, Estimation of the effective properties of two-dimensional cellular materials: a review, Theoretical and Applied Mechanics Letters 8 (4) (2018) 209-230. doi:10.1016/j.taml.2018.04.010.
URL https://linkinghub.elsevier.com/retrieve/pii/S2095034918301843

[3] J. Michel, H. Moulinec, P. Suquet, Effective properties of composite materials with periodic microstructure: a computational approach, Computer Methods in Applied Mechanics and Engineering 172 (1-4) (1999) 109-143. doi:10.1016/S0045-7825(98)00227-8.
URL https://linkinghub.elsevier.com/retrieve/pii/S0045782598002278

[4] M. G. Geers, V. Kouznetsova, W. A. Brekelmans, Gradient-enhanced computational homogenization for the micro-macro scale transition, Le Journal de Physique IV 11 (PR5) (2001) Pr5-145-Pr5-152. doi:10.1051/jp4:2001518.
URL http://www.edpsciences.org/10.1051/jp4:2001518

[5] V. Kouznetsova, M. G. D. Geers, W. A. M. Brekelmans, Multi-scale constitutive modelling of heterogeneous materials with a gradient-enhanced computational homogenization scheme, International Journal for Numerical Methods in Engineering 54 (8) (2002) 1235-1260. doi:10.1002/nme.541.
URL http://doi.wiley.com/10.1002/nme.541

[6] M. Geers, V. Kouznetsova, W. Brekelmans, Multi-scale computational homogenization: Trends and challenges, Journal of
Computational and Applied Mathematics 234 (7) (2010) 2175–2182.
doi:10.1016/j.cam.2009.08.077.

URL https://linkinghub.elsevier.com/retrieve/pii/
S0377042709005536

[7] G. Z. Voyiadjis, R. K. Abu Al-Rub, A. N. Palazotto, Thermodynamic framework for coupling of non-local viscoplasticity and non-local anisotropic viscodamage for dynamic localization problems using gradient theory, International Journal of Plasticity 20 (6) (2004) 981–1038. doi:10.1016/j.ijplas.2003.10.002.

URL https://linkinghub.elsevier.com/retrieve/pii/
S0749641903001414

[8] Y. Yang, A. Misra, Higher-Order Stress-Strain Theory for Damage Modeling Implemented in an Element-free Galerkin Formulation (2010) 37.

[9] F. Feyel, A multilevel finite element method (FE2) to describe the response of highly non-linear structures using generalized continua, Computer Methods in Applied Mechanics and Engineering 192 (28-30) (2003) 3233–3244. doi:10.1016/S0045-7825(03)00348-7.

URL https://linkinghub.elsevier.com/retrieve/pii/
S0045782503003487

[10] M. Geers, V. Kouznetsova, W. Brekelmans, Multi-scale computational homogenization: Trends and challenges, Journal of Computational and Applied Mathematics 234 (7) (2010) 2175–2182.
doi:10.1016/j.cam.2009.08.077.
URL https://linkinghub.elsevier.com/retrieve/pii/
S0377042709005536

[11] V. P. Nguyen, O. Lloberas-Valls, M. Stroeven, L. J. Sluys, Homogenization-based multiscale crack modelling: From micro-diffusive damage to macro-cracks, Computer Methods in Applied Mechanics and Engineering 200 (9-12) (2011) 1220-1236.
doi:10.1016/j.cma.2010.10.013.
URL https://linkinghub.elsevier.com/retrieve/pii/
S0045782510002914

[12] S. Loehnert, T. Belytschko, A multiscale projection method for macro/microcrack simulations, International Journal for Numerical Methods in Engineering 71 (12) (2007) 1466-1482. doi:10.1002/nme.2001.
URL http://doi.wiley.com/10.1002/nme.2001

[13] K. Matous, M. Kulkarni, P. Geubelle, Multiscale cohesive failure modeling of heterogeneous adhesives, Journal of the Mechanics and Physics of Solids 56 (4) (2008) 1511-1533. doi:10.1016/j.jmps.2007.08.005.
URL https://linkinghub.elsevier.com/retrieve/pii/
S0022509607001718

[14] M. G. Kulkarni, K. Matouš, P. H. Geubelle, Coupled multi-scale cohesive modeling of failure in heterogeneous adhesives, International Journal for Numerical Methods in Engineering 84 (8) (2010) 916-946. doi: 10.1002/nme.2923.
URL http://doi.wiley.com/10.1002/nme.2923
[15] C. Hirschberger, S. Ricker, P. Steinmann, N. Sukumar, Computational multiscale modelling of heterogeneous material layers, Engineering Fracture Mechanics 76 (6) (2009) 793-812.
doi:10.1016/j.engfracmech.2008.10.018.
URL https://linkinghub.elsevier.com/retrieve/pii/S0013794408002993

[16] M. Hossain, C.-J. Hsueh, B. Bourdin, K. Bhattacharya, Effective toughness of heterogeneous media, Journal of the Mechanics and Physics of Solids 71 (2014) 15-32. doi:10.1016/j.jmps.2014.06.002.
URL https://linkinghub.elsevier.com/retrieve/pii/S0022509614001215

[17] N. Nguyen, J. Yvonnet, J. Réthoré, A. B. Tran, Identification of fracture models based on phase field for crack propagation in heterogeneous lattices in a context of non-separated scales, Computational Mechanics 63 (5) (2019) 1047-1068. doi:10.1007/s00466-018-1636-z.
URL http://link.springer.com/10.1007/s00466-018-1636-z

[18] I. Goldhirsch, C. Goldenberg, On the microscopic foundations of elasticity, The European Physical Journal E 9 (3) (2002) 245-251. doi:10.1140/epje/i2002-10073-5.
URL http://link.springer.com/10.1140/epje/i2002-10073-5

[19] G. Francfort, J.-J. Marigo, Revisiting brittle fracture as an energy minimization problem, Journal of the Mechanics and Physics of Solids 46 (8) (1998) 1319-1342. doi:10.1016/S0022-5096(98)00034-9.
URL https://linkinghub.elsevier.com/retrieve/pii/
S0022509698000349

[20] A. Pandolfi, M. Ortiz, An eigenerosion approach to brittle fracture: AN EIGENEROSION APPROACH TO BRITTLE FRACTURE, International Journal for Numerical Methods in Engineering 92 (8) (2012) 694-714. doi:10.1002/nme.4352.
URL http://doi.wiley.com/10.1002/nme.4352

[21] N. Moes, J. Dolbow, T. Belytschko, A finite element method for crack growth without remeshing (1999) 20.

[22] X.-P. Xu, A. Needleman, Numerical simulations of fast crack growth in brittle solids, Journal of the Mechanics and Physics of Solids 42 (9) (1994) 1397-1434. doi:10.1016/0022-5096(94)90003-5.
URL https://linkinghub.elsevier.com/retrieve/pii/
0022509694900035

[23] N. Moës, C. Stolz, P.-E. Bernard, N. Chevaugeon, A level set based model for damage growth: The thick level set approach, International Journal for Numerical Methods in Engineering 86 (3) (2011) 358-380.
doi:10.1002/nme.3069.
URL http://doi.wiley.com/10.1002/nme.3069

[24] T. T. Nguyen, J. Réthoré, M.-C. Baietto, Phase field modelling of anisotropic crack propagation, European Journal of Mechanics - A/Solids 65 (2017) 279-288. doi:10.1016/j.euromechsol.2017.05.002.
URL https://linkinghub.elsevier.com/retrieve/pii/
S0997753816303321

[25] T. T. Nguyen, J. Yvonnet, M. Bornert, C. Chateau, K. Sab, R. Romani, R. Le Roy, On the choice of parameters in the phase field method for simulating crack initiation with experimental validation, International Journal of Fracture 197 (2) (2016) 213–226. doi:10.1007/s10704-016-0082-1.
URL http://link.springer.com/10.1007/s10704-016-0082-1

[26] C. Miehe, M. Hofacker, F. Welschinger, A phase field model for rate-independent crack propagation: Robust algorithmic implementation based on operator splits, Computer Methods in Applied Mechanics and Engineering 199 (45-48) (2010) 2765–2778.
doi:10.1016/j.cma.2010.04.011.
URL https://linkinghub.elsevier.com/retrieve/pii/
S0045782510001283

[27] M. J. Borden, C. V. Verhoosel, M. A. Scott, T. J. Hughes, C. M. Landis, A phase-field description of dynamic brittle fracture, Computer Methods in Applied Mechanics and Engineering 217-220 (2012) 77–95.
doi:10.1016/j.cma.2012.01.008.
URL https://linkinghub.elsevier.com/retrieve/pii/
S0045782512000199

[28] L. Ambrosio, V. M. Tortorelli, Approximation of functional depending on jumps by elliptic functional via t-convergence, Communications on Pure and Applied Mathematics 43 (8) (1990) 999–1036. doi:10.1002/
cpa.3160430805.

URL http://doi.wiley.com/10.1002/cpa.3160430805

[29] B. J. Glasser, I. Goldhirsch, Scale dependence, correlations, and fluctuations of stresses in rapid granular flows, Physics of Fluids 13 (2001) 407-420. doi:10.1063/1.1338543.

[30] A. Glacet, Study of quasi-periodic architectured materials: Vibrations, dynamic fracture and homogenization, Ph.D. thesis, L'Institut National des Sciences Appliquées de Lyon, Université de Lyon.

[31] E. N. Brown, Use of the tapered double-cantilever beam geometry for fracture toughness measurements and its application to the quantification of self-healing, The Journal of Strain Analysis for Engineering Design 46 (3) (2011) 167-186. doi:10.1177/0309324710396018.

URL http://journals.sagepub.com/doi/10.1177/0309324710396018

[32] T. M. Grabois, J. Neggers, L. Ponson, F. Hild, R. D. Toledo Filho, On the validation of integrated DIC with tapered double cantilever beam tests, Engineering Fracture Mechanics 191 (2018) 311-323. doi:10.1016/j.engfracmech.2017.12.015.

URL https://linkinghub.elsevier.com/retrieve/pii/S0013794417309232

[33] R. Li, Q. Shao, E. Gao, Z. Liu, Elastic anisotropy measure for two-dimensional crystals, Extreme Mechanics Letters 34 (2020) 100615. doi:10.1016/j.eml.2019.100615.
URL https://linkinghub.elsevier.com/retrieve/pii/
S2352431619302834

[34] C. M. Kube, Elastic anisotropy of crystals, AIP Advances 6 (9) (2016) 095209. doi:10.1063/1.4962996.
URL http://aip.scitation.org/doi/10.1063/1.4962996

[35] C. Zener, Contributions to the Theory of Beta-Phase Alloys, Physical Review 71 (12) (1947) 846-851. doi:10.1103/PhysRev.71.846.
URL https://link.aps.org/doi/10.1103/PhysRev.71.846

[36] A. Benallal, J.-J. Marigo, Bifurcation and stability issues in gradient theories with softening, Modelling and Simulation in Materials Science and Engineering 15 (1) (2007) S283-S295. doi:10.1088/0965-0393/15/1/S22.
URL https://iopscience.iop.org/article/10.1088/0965-0393/15/1/S22

[37] K. Pham, J.-J. Marigo, C. Maurini, The issues of the uniqueness and the stability of the homogeneous response in uniaxial tests with gradient damage models, Journal of the Mechanics and Physics of Solids 59 (6) (2011) 1163-1190. doi:10.1016/j.jmps.2011.03.010.
URL https://linkinghub.elsevier.com/retrieve/pii/S002250961100055X

[38] H. Amor, J.-J. Marigo, C. Maurini, Regularized formulation of the variational brittle fracture with unilateral contact: Numerical experiments, Journal of the Mechanics and Physics of Solids 57 (8) (2009) 1209-1229. doi:10.1016/j.jmps.2009.04.011.
[39] J. Réthoré, Automatic crack tip detection and stress intensity factors estimation of curved cracks from digital images: AUTOMATIC CRACK TIP DETECTION AND SIF ESTIMATION OF CURVED CRACKS, International Journal for Numerical Methods in Engineering 103 (7) (2015) 516-534. doi:10.1002/nme.4905.
URL http://doi.wiley.com/10.1002/nme.4905