Crack Extension and Possibility of Debonding in Encapsulation-Based Self-Healing Materials

The breakage of capsules upon crack propagation is crucial for achieving crack healing in encapsulation-based self-healing materials. A mesomechanical model was developed in this study to simulate the process of crack propagation in a matrix and the potential of debonding. The model used the extended finite element method (XFEM) combined with a cohesive zone model (CZM) in a two-dimensional (2D) configuration. The configuration consisted of an infinite matrix with an embedded crack and a capsule nearby, all subjected to a uniaxial remote tensile load. A parametric study was performed to investigate the effect of geometry, elastic parameters and fracture properties on the fracture response of the system. The results indicated that the effect of the capsule wall on the fracture behavior of the matrix is insignificant for tc/Rc ≤ 0.05. The matrix strength influenced the ultimate crack length, while the Young’s modulus ratio Ec/Em only affected the rate of crack propagation. The potential for capsule breakage or debonding was dependent on the comparative strength between capsule and interface (Sc/Sint), provided the crack could reach the capsule. The critical value of Sc,cr/Sint,cr was obtained using this model for materials design.


Introduction
As an alternative to conventional passive repair, bio-inspired self-healing materials are of great interest to material scientists because they have the built-in capability to repair structural damage either autogenously or with minimal help from an external stimulus.
A variety of innovative strategies have been explored and conducted to attempt to endow materials with self-healing abilities [1][2][3][4][5][6][7][8][9][10][11][12][13][14][15][16][17][18][19][20]. Encapsulation-based self-healing materials are one of the most promising techniques in this area, and they are receiving considerable attention since being initially demonstrated in the pioneering work of White et al. [6]. This process essentially involves crack propagation rupturing the embedded microcapsules, and the incorporated healing agent released into the crack face through certain physical mechanisms such as capillary action [6]. The polymerization of the healing agent, which can be triggered by the embedded catalyst, can narrow and even close the detached crack surfaces.
Therefore, the breakage of capsules takes a primary role in the realization of crack healing. In other words, the capsules need to be prudently designed with the appropriate geometry and mechanical properties for the capsules to break when a crack reaches them. To achieve this, the stress transferred from a crack to the capsules should be high enough to rupture the capsules. Unfortunately, a crack may propagate along the capsule-matrix interface and result in premature debonding at a weak interface. This is particularly problematic for inorganic hosts such as cement-based materials considering encapsulation techniques mostly apply to polymers. In contrast, a polymer matrix is normally chemically compatible with the capsule materials.
The full process of capsule breakage consists of two steps: (i) crack propagation in a matrix prior to reaching a capsule; and (ii) capsule rupture or debonding, depending on the interface properties. The former problem (i) of crack propagation in either homogeneous or heterogeneous solids has been widely investigated and described thoroughly in the literature. Likewise, the effects of an arbitrarily shaped hole near a crack tip on the stress field and stress intensity factor have been analyzed numerically, although this configuration is not a precise representation of the presence of a capsule [21][22][23][24]. A more accurate configuration for this case was proposed in which the stress concentration in a solid with an embedded capsule was evaluated for various interfacial conditions [25]. The onset of an interfacial crack due to a lack of debonding was clearly understood, although no initiation of the crack in the host was included in the study. In a similar manner, the stress concentration has been numerically analyzed but with an additional crack placed in the matrix with its tip right below the capsule [26]. Although this problem is partially understood for some specific cases, these findings provide good insight into further and comprehensive understanding about the complete event of crack propagation in a matrix and the potential of debonding as a result.
In this study, both the crack propagation in a matrix and the capsule-matrix interaction were investigated numerically using the extended finite element method (XFEM) combined with a cohesive zone model (CZM). The combined methods of XFEM and CZM were applied to fracture simulations of asphalt mixtures with heterogeneous microstructures in which a viscoelastic constitutive model was adopted for cracking in a matrix [27]. By taking advantage of XFEM, which does not require mesh to conform to geometric discontinuities, mesh refinement was not required in the region around the crack tip to capture the singular asymptotic fields adequately. The crack propagation in the matrix was exhaustively modeled by changing the mechanics of the matrix. However, it is difficult to simulate the interaction between the crack in a host and the capsule-matrix interface using only XFEM. Therefore, to determine the effect of an interfacial material on potential debonding, the cohesive zone model, achieved by pre-inserted zero-thickness cohesive elements (CIEs), was used to understand the stress transfer between matrix and capsule. It is noteworthy that the configuration proposed in this study also resembles the debonding of cylindrical capsules in a 2D domain. In that sense, this analysis also applies to material design for that case.

Numerical Setup
Consider a two-dimensional, infinite, isotropic, and homogeneous geometry of length 2L containing two circular holes of outer radius R c and a crack of 2l 0 in length with its center at the origin in between the two holes, as depicted in Figure 1a. A capsule with wall thickness t c should be inserted and bonded to each hole. This layout should represent the potential for the crack to freely propagate through the capsule on either side. Furthermore, it is assumed that the capsule size is small enough, when compared to the entire domain, to resemble the hypothetical condition of a holed plate with an embedded crack under a remote tensile stress ∑; thus, the edges of the linear crack 2l 0 are free from traction. Taking advantage of the symmetrical configuration, this analysis has been performed on half of the geometry, as shown in Figure 1b, in which the selected geometry and boundary conditions are schematically illustrated. Specifically, the full length of the domain was set to 2L = 10 mm and the outer radius to R c = 150 mm. Three various wall thicknesses ranging in different scales of magnitude were evaluated for t c /R c : 0.006, 0.05 and 0.12, which correspond with the true dimensions of microencapsulation-based self-healing specimens in lab tests [28][29][30]. A magnified view of the vicinity around the capsule, as highlighted in the red dashed box in Figure 1b, is schematically shown in Figure 2. The detailed composition of the matrix region is in grey and the capsule is in blue-green. The attached boundary surfaces of the capsule and matrix are highlighted in green and red, respectively. The solid region of the matrix was assigned a mesh with four-node isoperimetric plane strain elements. The free meshing method with mesh control was used such that the mesh was gradually refined approaching the crack and capsule, taking into account the considerably larger region of the matrix in comparison with the capsule. A mesh grid similar to that used in [25] has been introduced considering the present configuration is scaled down with respect to that in [25]. Mesh sensitivity analysis was conducted to ensure the model's accuracy, balancing computational accuracy and time. As a result, the mesh dimension progressively decreased from 50 µm to 10 µm at the left symmetrical edge, and the meshes on the other three edges (the top, bottom and right) were all 100 µm in size. There were ultimately 51,032 meshes in the matrix overall. The uniform four-node isoperimetric plane strain elements were used with 220 circumferential elements and four elements along the reference thickness (i.e., / c c t R  0.05) due to the radial symmetry. The debonding of the interface between the capsule and matrix was modeled using fournode zero-thickness cohesive elements (CIEs) with the given elastic and fracture properties, which are available in the finite element package ABAQUS Version 6.5 or higher. Figure 3 illustrates the final mesh grid used in this work. It should be noted that only a quarter of the plate can be modeled instead with nodes relaxing method for the crack for simplicity considering the axial symmetry of the configuration. In this way, the crack path has to be prior known and to align with element boundaries, which are not applicable for the case of arbitrary cracking due to mostly, heterogeneity of matrices, randomly distributed capsules, etc. The XFEM method does not require the mesh to match the geometry of the discontinuities and therefore no considerable mesh refinement is needed in the neighborhood of the crack tip to capture the singular asymptotic fields adequately. A magnified view of the vicinity around the capsule, as highlighted in the red dashed box in Figure 1b, is schematically shown in Figure 2. The detailed composition of the matrix region is in grey and the capsule is in blue-green. The attached boundary surfaces of the capsule and matrix are highlighted in green and red, respectively. The solid region of the matrix was assigned a mesh with four-node isoperimetric plane strain elements. The free meshing method with mesh control was used such that the mesh was gradually refined approaching the crack and capsule, taking into account the considerably larger region of the matrix in comparison with the capsule. A mesh grid similar to that used in [25] has been introduced considering the present configuration is scaled down with respect to that in [25]. Mesh sensitivity analysis was conducted to ensure the model's accuracy, balancing computational accuracy and time. As a result, the mesh dimension progressively decreased from 50 µm to 10 µm at the left symmetrical edge, and the meshes on the other three edges (the top, bottom and right) were all 100 µm in size. There were ultimately 51,032 meshes in the matrix overall. The uniform four-node isoperimetric plane strain elements were used with 220 circumferential elements and four elements along the reference thickness (i.e., t c /R c = 0.05) due to the radial symmetry. The debonding of the interface between the capsule and matrix was modeled using four-node zero-thickness cohesive elements (CIEs) with the given elastic and fracture properties, which are available in the finite element package ABAQUS Version 6.5 or higher. Figure 3 illustrates the final mesh grid used in this work. It should be noted that only a quarter of the plate can be modeled instead with nodes relaxing method for the crack for simplicity considering the axial symmetry of the configuration. In this way, the crack path has to be prior known and to align with element boundaries, which are not applicable for the case of arbitrary cracking due to mostly, heterogeneity of matrices, randomly distributed capsules, etc. The XFEM method does not require the mesh to match the geometry of the discontinuities and therefore no considerable mesh refinement is needed in the neighborhood of the crack tip to capture the singular asymptotic fields adequately.

Cohesive Crack Model for Capsule-Matrix Interaction
The surface interaction between the capsule and the matrix was modeled using the cohesive crack model [31][32][33], in which energy dissipation occurs in the fracture process zone (FPZ) during fracture. This nonlinear behavior is characterized by a traction-separation law, which is a typical bilinear response, as schematically shown in Figure 4. Here, the subscripts c, m, and int are short for capsule, matrix and interface, respectively. It is assumed that the tractions exist in the normal

Cohesive Crack Model for Capsule-Matrix Interaction
The surface interaction between the capsule and the matrix was modeled using the cohesive crack model [31][32][33], in which energy dissipation occurs in the fracture process zone (FPZ) during fracture. This nonlinear behavior is characterized by a traction-separation law, which is a typical bilinear response, as schematically shown in Figure 4. Here, the subscripts c, m, and int are short for capsule, matrix and interface, respectively. It is assumed that the tractions exist in the normal

Cohesive Crack Model for Capsule-Matrix Interaction
The surface interaction between the capsule and the matrix was modeled using the cohesive crack model [31][32][33], in which energy dissipation occurs in the fracture process zone (FPZ) during fracture. This nonlinear behavior is characterized by a traction-separation law, which is a typical bilinear response, as schematically shown in Figure 4. Here, the subscripts c, m, and int are short for capsule, matrix and interface, respectively. It is assumed that the tractions exist in the normal direction t n and shear direction t s across the crack surface for a 2D case, and the corresponding relative displacements are the crack opening displacement δ n and the crack sliding displacement δ s . The traction and deformation remain linearly related prior to the onset of damage, at which point the traction reaches a strength of t 0 n or t 0 s . Then, the traction decreases linearly and monotonically as a function of the corresponding deformation δ n or δ s . The irreversible degradation of the material is characterized by the progressive reduction in stiffness, defined by the slope k n or k s of the traction-displacement curve, as marked in Figure 4. The covered area is defined as the fracture energy G f I or G f II , which corresponds to t n = 0 or t s = 0 and δ n = δ  A scalar index D, which is a function of the effective displacement eff  (the combined effects of n  and s  ), characterizes the overall damage of the crack as given by the following: where max int,eff  is the maximum effective displacement attained during the loading history. The effective displacement is: and int, int, int, int, , 0 where  is the Macaulay bracket. This definition is based upon the fact that the materials do not undergo damage under pure compression.
Based on the initial stiffness  A scalar index D, which is a function of the effective displacement δ e f f (the combined effects of δ n and δ s ), characterizes the overall damage of the crack as given by the following: where δ max int,e f f is the maximum effective displacement attained during the loading history. δ 0 int,e f f and δ f int,e f f are the effective displacements at damage initiation and failure corresponding to δ 0 int,n and δ 0 int,s , and δ f int,n and δ f int,s . The value of D varies from 0 to 1, representing an intact and a fully cracked material, respectively.
The effective displacement is: and where · is the Macaulay bracket. This definition is based upon the fact that the materials do not undergo damage under pure compression. Based on the initial stiffness k 0 int,n and k 0 int,s of the intact material, the unloading and reloading stiffness k int,n and k int,s can then be determined according to: The tractions t int,n and t int,s are written as functions of the corresponding displacements δ int,n and δ int,s respectively, as Concerning the damage initiation, a quadratic nominal stress law is assumed to combine both the effects of normal traction and tangential traction as given by: where t 0 int,n and t 0 int,s refer to the maximum value of normal traction and tangential traction of the interface, respectively.

XFEM-Based Cohesive Behavior for Crack Propagation in Matrix
An extended finite element method was used to model the propagation of a pre-defined crack in a matrix prior to reaching a capsule. This method is based on the concept of partition of unity, which allows local special enriched functions to be incorporated into a finite element approximation to represent discontinuities in a crack [34][35][36]. The enrichment functions typically consist of the near-tip asymptotic functions that representthe singularity around the crack tip and a discontinuous function that describes the abrupt change in displacement across the crack faces [34][35][36]. The approximation for a displacement vector function with the partition of unity enrichment is given by [34][35][36]: where N I (x) is the usual nodal shape function, u I is the usual nodal displacement vector associated with the continuous part of the finite element solution, a I is the nodal enriched degree of freedom vector, H(x) is the associated discontinuous jump function across the crack surfaces as written as Equation (10) below, b α I is the product of the nodal enriched degree of freedom vector and F α (x) is the associated elastic asymptotic crack-tip function that can be determined according to Equation (11).
where x is a sample (Gauss) point, x * is the point on the crack closest to x, and n is the unit outward normal to the crack at x * . In a polar coordinate system (r, θ) with its origin at the crack tip, F α (x) is: where θ = 0 is tangent to the crack at the tip. The damage and failure of an enriched element is based on the cohesive response that consists of a damage initiation criterion and a damage evolution law. A linear elastic response is assumed for the initially undamaged material. Once damage is initiated, the elements will degrade progressively according to a given damage evolution law. Figure 5 schematically depicts a typical and commonly used linear traction-separation response, where the traction t m = t 0 m starts to decrease to zero corresponding to the maximum displacement δ f m at which the damage initiation criterion is met. For the sake of simplicity, it is assumed that the enriched elements do not undergo damage under pure compression and that the response is mode-independent.
where 0   is tangent to the crack at the tip.
The damage and failure of an enriched element is based on the cohesive response that consists of a damage initiation criterion and a damage evolution law. A linear elastic response is assumed for the initially undamaged material. Once damage is initiated, the elements will degrade progressively according to a given damage evolution law. For the sake of simplicity, it is assumed that the enriched elements do not undergo damage under pure compression and that the response is mode-independent. The maximum principal stress criterion is used for damage initiation in the model, as given by: Here, m D can be determined in the same way as int D using Equations (1)-(3).

Input Parameters for Material Properties
The material parameters used in this study are listed in Table 1. Both the effects of elasticity and fracture were examined. The quantity given by the ratio c E / m E is hereafter referred to as the elastic The maximum principal stress criterion is used for damage initiation in the model, as given by: Here, σ 0 m and σ max m represent the maximum allowable principal stress and the maximum principal stress at present, respectively.
Similarly, the tractions can be determined by the damage variable D m , and the strengths, t 0 m,n for normal opening and t 0 m,s for shear sliding of the matrix can be determined: D m can be determined in the same way as D int using Equations (1)-(3).

Input Parameters for Material Properties
The material parameters used in this study are listed in Table 1. Both the effects of elasticity and fracture were examined. The quantity given by the ratio E c /E m is hereafter referred to as the elastic ratio to illustrate the influence of the elastic mismatch on cracking and on debonding potential. Various matrix strengths were selected to elucidate crack propagation through the matrix prior to the onset of debonding at an interface. The host strengths correspond to a wide range of materials, including cement-based materials (mean value of 3.5 MPa [37]) and polymer matrices such as epoxy resins (39 ± 4 MPa in [38]). Various interface strengths were modeled, ranging from weak bonding (0.1 MPa) to perfect bonding (10 times of matrix strength, i.e., 35 MPa) between the capsule and matrix. Regarding the fracture properties of matrix, the literatures provide tentative values of G f m for normal concrete (0.15 N/mm in [37]) and for cement-based materials with a certain additives, such as the well-known fiber-reinforced cement (2.1 N/mm in [39]). Taking into account , as presented in Figure 6. For both the matrix and capsule, the shear fracture properties were simply assumed to be the same as the normal fracture properties: 15 N/mm, and t 0 int /t 0 m = 10 for perfect bonding or t 0 int /t 0 m = 1 for imperfect bonding, respectively. The Poisson's ratio ν was kept at 0.18 throughout.
for normal concrete (0.15 N/mm in [37]) and for cement-based materials with a certain additives, such as the well-known fiber-reinforced cement (2.1 N/mm in [39]).
Taking into account performed by changing either: (a) the strength , as presented in Figure 6. For both the matrix and capsule, the shear fracture properties were simply assumed to be the same as the normal fracture properties:  Figure 7 shows the typical maximum in-plane stress principle in the vicinity of crack tip at interval of progressive crack extension until it reaches to the capsule wall for the reference with perfect bonding (i.e., t 0 int /t 0 m =10). To track every step of crack extension, the output was recorded at time interval of 0.002 of σ/ ∑. The stress concentration is observed at the crack tip (see Figure 7a) prior to its propagation when the maximum stress increases to 3.5 MPa (see Figure 7b). The crack has to propagate across an entire element at a time to avoid the need to model the stress singularity. Unloading occurs to the stress concentrated elements when the crack propagates and the stress concentration shifts to the new crack tip, e.g., crack extension as shown in Figure 7c,d. Depending on the respective mechanical properties, the crack will either propagate along the interface or rupture the wall when the crack reaches to the capsule wall (see Figure 7e). Here, the crack cannot propagate in either way due to a quite high interface strength (i.e., t 0 int /t 0 m = 10) and the elastic response of capsule wall as assumed. Therefore, the stress accumulates between the interface and the capsule wall. In other words, the crack can propagate through the wall if its strength is less than 6.3 MPa as shown in Figure 7f once the interface is sufficiently strong to avoid debonding to occur prior to the wall breakage.   Figure 8 shows the normalized ultimate crack extension with respect to its initial length, denoted by ∆l/l 0 = (l − l 0 )/l 0 and labeled on the left y-axis, for various initial crack sizes, from 50 µm up to 325 µm, as a function of the matrix strength. The maximum crack growth up to the capsule wall (i.e., l = 350 µm) is presented on the right y-axis. The value of ∆l/l 0 initially exhibits a sharp decrease with increasing matrix strength for t 0 m ≤ 10 MPa. The rate of decrease eventually slows down, with the exception of the crack 325 µm in length, which remains at its maximum value (i.e., 0.077) until the matrix strength increases to 20 MPa. For the remaining crack sizes, ∆l/l 0 is less than its maximum value regardless of the matrix strength. These results imply that the initial crack length should be no less than 325 µm and the matrix strength should be no higher than 20 MPa for the crack to propagate to the capsule at the present load and boundary conditions. Thus, the initial crack length was fixed at 325 µm in the subsequent study to keep the focus on crack propagation through the matrix, resulting in either capsule rupture or capsule bypass.   Figure 9b. The fracture energy ranged from values representing normal, plain cement [37] to values for fiber-reinforced cement to ensure universal coverage of materials [39]. However, the fracture energies selected did not show an obvious influence on the crack propagation. This can be explained by the choice of the criteria for crack initiation, which is stress controlled in this model, as given by Equation (12), while the fracture energy strongly affects the damage evolution in course of displacement, as characterized by the damage index defined in Equation (1). Besides, the crack extension is much smaller than its full length, i.e., 25/(325 + 25) = 0.071. Likewise, the effect of fracture energy on a capsule with constant strength on the likelihood of debonding is insignificant until 0.2-0.4 of the relative crack extension is reached with respect to the circumference of the full interface [25].  Figure 9a shows the detailed crack propagation process, as normalized by its initial length, l 0 = 325 i.e., µm, for a load applied to matrices with various strengths. The crack extended to the capsule wall for cases in which ∆l/l 0 = 0.077 (0.055 or 0.022) provided the matrix strength was no more than 20 MPa, 30 MPa, and 40 MPa, respectively. A higher load increase was necessary to initiate crack propagation for an increase of matrix strength.

Crack Extension in Matrix for Perfect Bonding
The influence of the fracture energy of the matrix with constant strength t 0 m = 3.5 MPa on crack extension is presented in Figure 9b. The fracture energy ranged from values representing normal, plain cement [37] to values for fiber-reinforced cement to ensure universal coverage of materials [39]. However, the fracture energies selected did not show an obvious influence on the crack propagation. This can be explained by the choice of the criteria for crack initiation, which is stress controlled in this model, as given by Equation (12), while the fracture energy strongly affects the damage evolution in course of displacement, as characterized by the damage index defined in Equation (1). Besides, the crack extension is much smaller than its full length, i.e., 25/(325 + 25) = 0.071. Likewise, the effect of fracture energy on a capsule with constant strength on the likelihood of debonding is insignificant until 0.2-0.4 of the relative crack extension is reached with respect to the circumference of the full interface [25]. model, as given by Equation (12), while the fracture energy strongly affects the damage evolution in course of displacement, as characterized by the damage index defined in Equation (1). Besides, the crack extension is much smaller than its full length, i.e., 25/(325 + 25) = 0.071. Likewise, the effect of fracture energy on a capsule with constant strength on the likelihood of debonding is insignificant until 0.2-0.4 of the relative crack extension is reached with respect to the circumference of the full interface [25].  Figure 10. As discussed in [25], an appropriate range of values for the capsule elasticity were used, The effect of the Young's modulus ratio E c /E m on the crack extension in the matrix is shown in Figure 10. As discussed in [25], an appropriate range of values for the capsule elasticity were used, ranging from polymer elasticity (≥1 GPa) [40] to the elasticity of stiffer materials, such as ceramic or conventional glass (≤70 GPa) [41]. The Young's modulus was kept at 25 GPa for the matrix material as a reference. The increase of E c /E m delays the crack initiation and growth but does not change the ultimate crack length. In other words, the crack was always able to propagate to the capsule with the Young's modulus selected in the present study. Compared with the aforementioned results, it is proposed that crack extension in a matrix is influenced significantly by the strength of the matrix, which governs the ultimate crack length, whereas the Young's modulus ratio E c /E m only affects the rate of crack propagation. It should be noted that the presence of defects (inclusions or voids) strongly affects the fracture mechanism because of the interaction between the defects and the cracks. A crack is either attracted or repelled by inclusions or voids depending on the comparative rigidity of the defects to the matrix. The final crack trajectory is a result of the combined effect of the inclusions or voids. In the present model, the void left by the capsule attracts the crack but the crack can be either attracted or deflected depending on the stiffness of the capsule wall.   The nodal stress values of the interface, matrix and capsule, labeled "_A" and "_B" for the upper and lower elements, respectively, as shown in Figure 11a, are collected and presented in Figure 11b. The other elastic and fracture indexes were all kept the same as reference values. The stress progressively increases with applied load, and the interface exhibits a quite low stress compared to that of the matrix and capsule. One striking feature is that, for the same phase, the nodal stress curves _A and _B are coincident for the elements until a certain load, i.e., /   = 0.132, as marked by a The nodal stress values of the interface, matrix and capsule, labeled "_A" and "_B" for the upper and lower elements, respectively, as shown in Figure 11a, are collected and presented in Figure 11b. The other elastic and fracture indexes were all kept the same as reference values. The stress progressively increases with applied load, and the interface exhibits a quite low stress compared to that of the matrix and capsule. One striking feature is that, for the same phase, the nodal stress curves _A and _B are coincident for the elements until a certain load, i.e., σ/ ∑ = 0.132, as marked by a dashed line in Figure 11b, after which the stress curves diverge. The stress of the matrix and capsule are very close until 0.132 of σ/ ∑, while the interface stays at a comparatively low stress. This demonstrates that the capsule and matrix are attached quite well, as there is a complete stress transition between the two. The reason for the abrupt change of stress at σ/ ∑ = 0.132 is that the crack tip reaches the boundary of the matrix but cannot propagate through the CIEs of interface using XFEM and hereafter the crack is characterized by the damage of CIEs. The stress at which the curves split are taken as the critical stress S cr , corresponding to the crack reaching the interface.  Figure 12 shows the effect of the Young's modulus ratio on the critical stress of the interface, which is normalized by the remote load applied to a matrix with strength of either 3.5 MPa (Figure 12a Figure 12 shows the effect of the Young's modulus ratio on the critical stress of the interface, which is normalized by the remote load applied to a matrix with strength of either 3.5 MPa (Figure 12a) or 10 MPa (Figure 12b). The higher thickness of the capsule wall results in a higher critical stress, which increases with E c /E m because of concentrated stress. The effect of E c /E m on the critical stress becomes more significant with an increase of t c /R c . Similarly, the curves for the matrix with t 0 m = 10 MPa exhibit the same trend, but the critical stress increases by approximately three times when t 0 m = 3.5 MPa. Thus, the critical stress of the interface S int,cr , as one of the crucial requirements for preventing detachment of the matrix from the capsule, can be determined for materials design. The other requirement is the critical stress of the capsule S c,cr because the onset of debonding depends on the comparative strength between the capsule and interface (S c /S int ), which will be discussed further. For instance, given t 0 m = 3.5 MPa, E c /E m = 1 and t c /R c = 0.05, the critical strength S int,cr is 0.33 × ∑.
for preventing detachment of the matrix from the capsule, can be determined for materials design.
The other requirement is the critical stress of the capsule     (Figure 13b). The value of S cr / ∑ increases almost linearly with the elastic ratio for the three levels of thickness: a value of 0.006 for t c /R c exhibits quite a low S cr / ∑ ratio for the same E c /E m value when compared with other results. However, the stress curves are nearly coincident with a further increase of t c /R c from 0.05 to 0.12. For t 0 m = 10 MPa, the curve exhibits a same trend, but the critical stress is approximately three times higher than with t 0 m = 3.5 MPa. It can be concluded that there is an optimal selection for the capsule thickness, i.e., 0.05 of t c /R c , for optimal cost saving. Other than this, the critical capsule stress S c,cr is known to ensure that the capsule can be ruptured by crack propagation, provided the interface is strong enough. For instance, the critical stress is less than one times the applied load (S c,cr / ∑ < 1) provided t 0 m = 3.5 MPa, E c /E m = 1 and t c /R c = 0.05.

S
are the critical stresses obtained in the model for the capsule and interface, respectively. Taking the reference as an example, int / c S S should be no more than 1/0.33 ≈ 3 inaccordance with Figures 12a and 13a. In general, the potential for capsule breakage is dependent on the comparative strength between the capsule and interface, provided a crack can reach the capsule with a given load, as discussed above. This can be written as [42][43][44][45]:

Debonding Due to Capsule-Matrix Interaction
S c /S int ≤ S cr,c /S cr,int where S c and S int are the design strengths of the capsule and interface, respectively. Likewise, S cr,c and S cr,int are the critical stresses obtained in the model for the capsule and interface, respectively.
Taking the reference as an example, S c /S int should be no more than 1/0.33 ≈ 3 inaccordance with Figures 12a and 13a. Figure 14 shows the progressive crack extension for the reference with imperfect bonding (i.e., t 0 int /t 0 m = 1). Stress concentration is observed in the capsule wall as the crack extends along the interface and the debonding initiates at 0.25 of σ/ ∑, as characterized by damage of CIEs, i.e., the scalar damage variable (SDEG) higher than zero, as shown in Figure 15, comparing to the condition of perfect bonding (see Figure 7) where no damage of CIEs occurs.      Figure 16 shows the debonding evolution as a function of the applied load for the reference with a variety of fracture properties of the interface. The effect of the strength t 0 int with δ 0 int kept constant on the crack growth has been presented in Figure 16a. Conversely, the results for a variety of fracture energies G 0 int with t 0 int kept constant are reported in Figure 16b. The onset of the capsule debonding shifts to a higher load for a higher interfacial strength. The curves can be divided into two groups by 0.29 of t 0 int /t 0 m . This agrees well with the results for a stress of S cr = 0.33· ∑ = 0.33 × 3 ≈ 1.0 (t 0 int = 0.29·t 0 m = 0.29 × 3.5 ≈ 1.0) as shown in Figure 12a. For the condition t 0 int /t 0 m < 0.29, the crack extension breaks into two stages, with an abrupt growth at 0.45 and 0.7 of σ/ ∑ for t 0 int /t 0 m = 0.11 and 0.17, respectively. This can be explained by debonding also initiating from other locations due to stress concentration, as shown in Figure 17. The cracks ultimately merge together, leading to complete detachment of the capsule-matrix, i.e., ∆l/(2πR c ) = 1. Regarding the t 0 int /t 0 m ≥ 0.29 condition, debonding still occurs because capsule response is assumed to be elastic in the present model. When t 0 int increases to 7.14 times that of the matrix strength, the capsule-matrix debonding does not initiate regardless of the capsule mechanical properties under given load conditions at present.  [25]. In other words, there seems to exist a limit value of the toughness of bonding after which the crack growth along the interface is hardly affected [25]. A smaller value of   same curves before that 0.2-0.4 of the circumference of the interface is detached [25]. In other words, there seems to exist a limit value of the toughness of bonding after which the crack growth along the interface is hardly affected [25]. A smaller value of t 0 int /t 0 m = 0.11 was selected for further analysis, and debonding initiated at a lower load than expected. It exhibits separate trend of crack extension after ∆l/(2πR c ) = 0.2. The critical stage shifts further to a load as high as 0.5 followed by 0.65 of σ/ ∑, as highlighted in dashed lines in Figure 16b. However, the onset of debonding does not change with fracture energy for all cases.    Figure 19 shows that the SDEG spatially distributes along the capsule-matrix interface at load intervals corresponding to the first complete debonding, as shown in Figure 16a (i.e.,   Figure 18 shows the debonding evolution for the reference but with t 0 m = 10 MPa. Again, the higher strength results in the initiation of debonding at the higher applied load. The critical condition for debonding is t 0 int /t 0 m = 0.27 according to S cr / ∑ = 0.89, as shown in Figure 12b. Likewise, debonding occurs at multiple locations for t 0 int /t 0 m > 0.27 and otherwise. The detachment between the capsule and the matrix can be perfectly avoided when the value of t 0 int /t 0 m increases to 1.5. A quite low strength of the interface, i.e., t 0 int /t 0 m = 0.04, still leads to complete debonding.    Figure 19 shows that the SDEG spatially distributes along the capsule-matrix interface at load intervals corresponding to the first complete debonding, as shown in Figure 16a (i.e.,  Figure 19 shows that the SDEG spatially distributes along the capsule-matrix interface at load intervals corresponding to the first complete debonding, as shown in Figure 16a (i.e., t 0 int /t 0 m = 0.11 and t 0 m = 3.5) and Figure 18 (i.e., t 0 int /t 0 m = 0.04 and t 0 m = 10). The damage of elements at the interface initiates when the crack tip approaches, i.e., θ = 0 • . The detachment propagates progressively along the interface with the applied load until it is close to the compression stressed zone, where the direction of the cross section is nearly parallel to the remote tensile load, i.e., θ = 90 • and 270 • . At that point, the debonding initiates from the other side by σ/ ∑ = 0.6, as highlighted in the blue dashed line in Figure 19a. The exception is at the vicinity of θ = 180 • , which remains intact due to the resultant stress state. With further load increases, the dominating stress in the vicinity of θ = 180 • transitions to tension and, thus, debonding occurs there. For the compression stressed zone, the intact section narrows with a load increase. The remaining CIEs also initiate damage with a minimum SDEG = 0.5 at σ/ ∑ = 1, and at this point, complete debonding occurs. The same trend is observed in Figure 19b, but one difference is that the SDEG of CIEs is approximately 0.8-0.9 at σ/ ∑ = 0.2 rather than 1 for t 0 int /t 0 m = 0.11 and t 0 m = 3.5, as shown in Figure 19a. The debonding propagates in the same way for various load intervals until the SDEG increases to 1. progressively along the interface with the applied load until it is close to the compression stressed zone, where the direction of the cross section is nearly parallel to the remote tensile load, i.e.,   90° and 270°. At that point, the debonding initiates from the other side by /    0.6, as highlighted in the blue dashed line in Figure 19a. The exception is at the vicinity of   180°, which remains intact due to the resultant stress state. With further load increases, the dominating stress in the vicinity of   180° transitions to tension and, thus, debonding occurs there. For the compression stressed zone, the intact section narrows with a load increase. The remaining CIEs also initiate damage with a minimum SDEG = 0.5 at /    1, and at this point, complete debonding occurs. The same trend is observed in Figure 19b, but one difference is that the SDEG of CIEs is approximately 0.8-0.9 at /    0.2 rather than 1 for 0 0 int / m t t  0.11 and 0 m t  3.5, as shown in Figure 19a. The debonding propagates in the same way for various load intervals until the SDEG increases to 1.      In Figure 21b, the debonding occurs in the vicinity of   0° for 0 0 int / m t t  1, 1.43 and 4.29, beyond which the onset of debonding is avoided.  For the capsule with a higher thickness, i.e., t c /R c = 0.12, similar results are observed, as shown in Figure 21a. No debonding is observed for an interface with strength higher than 7.14 times the strength of the matrix. Given t 0 int /t 0 m = 1, the ultimate crack growth ∆l/(2πR c ) = 0.26 is greater than the reference value of 0.19 (i.e., t 0 m = 3.5 MPa), followed by a value of 0.15 for higher matrix strength (i.e., t 0 m = 10 MPa) and of 0.03 for a decrease of Young's modulus ratio (i.e.,E c /E m = 0.4). In Figure 21b, the debonding occurs in the vicinity of θ = 0 • for t 0 int /t 0 m = 1, 1.43 and 4.29, beyond which the onset of debonding is avoided.

Conclusions
A mesomechanical model was proposed in this study to simulate the fracture behavior of encapsulation-based healing materials. The crack propagation in a matrix and the capsule-matrix interaction were investigated numerically using the extended finite element method (XFEM) combined with a cohesive zone model (CZM) in a 2D configuration consisting of an infinite matrix with an embedded crack and a capsule nearby. A remote uniaxial tensile and uniform load was applied on the system, and a typical bilinear cohesive law was used in the model. The effect of geometry, elastic parameters and fracture properties on crack extension in the matrix and the potential of debonding were investigated in the parametric study. The following conclusions can be drawn from the analysis: (i) For crack propagation in matrix prior to reaching a capsule: To ensure that the crack can propagate and reach a capsule, the initial crack length should be no less than 325 µm, and the matrix strength should be no higher than 20 MPa for the present loading conditions.
Regarding the fracture energy of matrix, the modification of its strength exhibits a more significant effect on the ultimate crack length than changing the maximum

Conclusions
A mesomechanical model was proposed in this study to simulate the fracture behavior of encapsulation-based healing materials. The crack propagation in a matrix and the capsule-matrix interaction were investigated numerically using the extended finite element method (XFEM) combined with a cohesive zone model (CZM) in a 2D configuration consisting of an infinite matrix with an embedded crack and a capsule nearby. A remote uniaxial tensile and uniform load was applied on the system, and a typical bilinear cohesive law was used in the model. The effect of geometry, elastic parameters and fracture properties on crack extension in the matrix and the potential of debonding were investigated in the parametric study. The following conclusions can be drawn from the analysis: (i) For crack propagation in matrix prior to reaching a capsule: To ensure that the crack can propagate and reach a capsule, the initial crack length should be no less than 325 µm, and the matrix strength should be no higher than 20 MPa for the present loading conditions.
Regarding the fracture energy of matrix, the modification of its strength t 0 m (δ f m = constant) exhibits a more significant effect on the ultimate crack length than changing the maximum displacement δ f m (t 0 m = constant), which is because a stress-based criterion for crack initiation was used in this model. The Young's modulus ratio E c /E m only affects the rate of crack propagation and does not change the ultimate crack length.
(ii) For the potential of debonding due to capsule-matrix interaction: The stress of the interface and capsule as the crack reaches the capsule can be known and used as the critical stress for material design.
t c /R c = 0.05 is an optimal selection for the capsule thickness value according to the critical stress field.
The potential for capsule breakage or debonding is dependent on the comparative strength of the capsule and interface (S c /S int ), provided a crack can reach the capsule with a given load. The critical value of S cr,c /S cr,int for capsule breakage can be obtained from the model for design as S c /S int ≤ S cr,c /S cr,int , whereas debonding occurs otherwise.
For the reference example, debonding can be avoided provided t 0 int /t 0 m ≥ 7.14, regardless of the capsule's mechanical properties for the present load and boundary conditions. Some limitations of the present results are the following:(i) no arbitrary cracking due to, e.g., heterogeneous and anisotropic matrix, multiple capsules with a random distribution, etc., has been included; (ii) only mode-I is assumed to be the dominant failure mode; (iii) the prescribed load is uniform and uniaxial tension; (iv) the influence of Poisson's ratio is not considered; (v) the fracture indexes are kept the same in the normal and tangential direction; and (vi) the capsule response is assumed to be linear elastic. These factors can be included in the present model framework for more general considerations and the associated study is in process.