Characteristic Length of Meso-Structure as a Physical Fact for a Class of Brittle Materials

1. Introduction

Recently, the nonlocal macro-meso-scale consistent damage (NMMD) model (Lu and Chen 2020 [1]; Chen et al. 2021 [2]) has been applied to various static and dynamic fracture problems in the manner of articulation, accuracy and effectiveness. These applications include mixed-mode cracking problems (Ren et al. 2024 [3]; Lu et al. 2024 [4]), direct-shear fracture (Ren et al. 2023 [5]), dynamic crack branching (Lu and Chen 2021 [6]), compression fracture of concrete (Ren et al. 2024 [7]), anisotropic cracking of rocks (Xia et al. 2024 [8]) and rate-dependent failure (Zhao et al. 2024 [9]). However, the influence radius (also called the characteristic length) in NMMD model, though proved to be necessary for the mesh insensitivity of strain-softening regime, remains to be estimated indirectly with considerable arbitrariness.

The NMMD model still belongs to the family of nonlocal damage theories, which can be traced back to the remarkable works of Eringen (1981 [10]) and Bažant et al. (1984 [11]). Whatever form is used in constructing damage variable (Bažant and Jirásek 2002 [12]), the integral-type or the gradient-enhanced approach, (strictly speaking, the NMMD model should be classified as the integral-type), the common feature of these nonlocal theories is the involvement of characteristic length to avoid the mesh-size sensitivity. Nonetheless, there is no broad acceptance of consensus on the interpretative declaration of characteristic length in nonlocal models, such as nonlocal-integral damage model (Bažant and Pijaudier-Cabot 1989 [13]; Bažant and Jirásek 2002 [12]), peridynamics (Silling 2000 [14]; Silling and Lehoucq 2007 [15]), phase field model (Bourdin et al. 2000 [16]; Bourdin et al. 2008 [17]; Wu 2017 [18]; Feng et al. 2021 [19]) and gradient-enhanced model (Peerlings et al. 1996 [20]; Xue and Ren 2024 [21]). Natural questions emerge about the characteristic length employed in nonlocal models. Is it a mathematically fictitious parameter or physically-based material property? If it is a material property, then is it identifiable and measurable by the underlying physical experiment? How can it relate to ordinary material-constants?

This paper, focusing on the NMMD model and its macro-meso-scale consistent mechanism, aims to answer these above questions by experimental, numerical and theoretical investigations. It should be stressed that we limit ourselves to brittle materials like PMMA and currently do not deal with quasi-brittle and ductile fractures for the sake of preliminary exploration.

The remainder of this paper is organized as follows. Section 2 firstly states the fact that the conventional local criteria cannot adequately interpret and predict the experimentally-observed results of tests with various sizes of defect due to the lack of involving characteristic length. Sections 3 briefly introduces the basic idea of NMMD model and its application to the analysis of brittle materials like polymethyl-methacrylate (PMMA). The corresponding results and discussions are arranged in Section 4, where the theoretical calibration of the characteristic length in NMMD model is identified for the first time. Comparisons with other modified criteria involving characteristic length are also discussed in Section 5. Besides, in Section 6 theoretical calibration of characteristic length is proved to be suitable for problems of mixed-mode fracture with stress singularity. Finally, the main conclusions are drawn in Section 7, followed by an Appendix A.

2. Statement of the Problem: Crack Onset at Non-Singular Stress Concentration

Shown in Figure 1 is the dog-bone specimen made of PMMA material, in which a circular hole with the diameter $2R=\left\{0.5\mathrm{mm},1\mathrm{mm},2\mathrm{mm},4\mathrm{mm}\right\}$ is preset at the central of the specimen. These specimens, subjected to a uniaxial tension $\overline{\sigma }$ were experimentally conducted by Sapora et al. (2018 [22]) with the following dimensions: $L=100\mathrm{mm}$, $W=40\mathrm{mm}$, and thickness $B=10\mathrm{mm}$. The basic mechanical properties of PMMA material (Sapora et al. 2018 [22]) are listed in Table 1, and the experimental results of failure stress ${\sigma }_{\mathrm{f}}$, which is defined as the ultimate load divided by the minimum cross area, against the hole diameter are plotted in Figure 2. It can be seen that the failure stress of central-holed specimens strongly depends on the size of circular hole. In particular, the failure stress asymptotically and nonlinearly increases as the diameter of circular hole decreases, and it reaches the tensile strength ${\sigma }_{\mathrm{t}}=70.5\mathrm{MPa}$ as the size of circular hole tends to zero.

As is well known, referring to the polar coordinate system $\left(r,\theta \right)$, the stress distribution around the circular hole within an infinite plate under the uniaxial loading is given by (Timoshenko and Goodier 1951 [42], p.80)

$$\left\{\begin{array}{c}{\sigma }_{rr}={\displaystyle \frac{\overline{\sigma }}{2}}\left(1-R^2/r^2\right)+{\displaystyle \frac{\overline{\sigma }}{2}}\left(1-R^2/r^2\right)\left(1-3R^2/r^2\right)cos2\theta ,\hfill \\ {\sigma }_{\theta \theta }={\displaystyle \frac{\overline{\sigma }}{2}}\left(1+R^2/r^2\right)-{\displaystyle \frac{\overline{\sigma }}{2}}\left(1+3R^4/r^4\right)cos2\theta ,\hfill \\ {\sigma }_{r\theta }=-{\displaystyle \frac{\overline{\sigma }}{2}}\left(1-R^2/r^2\right)\left(1+3R^2/r^2\right)sin2\theta .\hfill \end{array}\right.$$

(1)

It is checked in Appendix that Equation (1) is still sufficiently suitable for the dog-bone specimens with the finite size.

For brittle cracking problems with non-singular stress concentration, in engineering design and analysis it was usually believed that either the maximum stress (MS) criterion

$${\sigma }_{max}={\sigma }_{\theta \theta }\left(R,\pm \pi /2\right)⩾{\sigma }_{\mathrm{t}}⟹\left\{\begin{array}{c}R=0{\sigma }_{\mathrm{f}}={\sigma }_{\mathrm{t}}\hfill \\ R>0:{\sigma }_{\mathrm{f}}={\sigma }_{\mathrm{t}}/3\hfill \end{array}\right.$$

(2)

or the net stress (NS) criterion

$${\sigma }_{\mathrm{net}}={\displaystyle \frac{\overline{\sigma }W}{W-2R}}⩾{\sigma }_{\mathrm{t}}⟹\forall R⩾0:{\sigma }_{\mathrm{f}}=\left(1-2R/W\right){\sigma }_{\mathrm{t}}$$

(3)

would provide predictions of satisfaction. The predicted results of the two commonly-used criteria are also plotted in Figure 2, however, excepting the case without hole it is clearly observed that both of them are significantly different from the mean values of experimental results.

In terms of the cases of $R>0$, the discrepancy between predictions of the two commonly-used criteria and experimental results manifests that the NS criterion predicts the upper bound of the real failure stress, while the MS criterion provides the lower bound. Similar experimental results (e.g., Carter 1992 [23], Li and Zhang 2005 [24]; Torabi et al. 2022 [25]), though not shown here, also confirm that both the MS criterion and the NS criterion, without involving some kind of characteristic length, are far from being satisfactory to handle crack onset from a hole boundary.

In the following, the NMMD model will be employed to analyze the above dog-bone specimens and provide insight into the problem of crack onset at non-singular stress concentration.

3. The NMMD Model for PMMA-like Materials

In this section, the basic idea of NMMD model and its application to the analysis of PMMA-like materials will be briefly summarized and described correspondingly. For more aspects of theoretical detail and numerical implementation of NMMD model, the interested reader can refer to Lu and Chen (2020 [1]) and Chen et al. (2021 [2]).

3.1. General Formulation of NMMD Model

Let $\mathcal{B}\subset {\mathbb{R}}^d$ ($d=1,2,3$) be the reference configuration of a solid body, and $\partial \mathcal{B}$ denote its external boundary. Given the influence radius (also called the characteristic length) $\ell \in {\mathbb{R}}^{+}$, the mesoscopic structure (see Figure 3) assigned to the solid body is defined as

$${\mathcal{M}}^{〈\ell 〉}=\left\{\left({\mathit{x}}^{\prime },\mathit{x}\right):{\mathit{x}}^{\prime }\in {\mathcal{K}}_{\mathit{x}}^{〈\ell 〉}\cap \mathcal{B},\mathit{x}\in \mathcal{B}\right\},$$

(4)

where ${\mathcal{K}}_{\mathit{x}}^{〈\ell 〉}=\left\{{\mathit{x}}^{\prime }\in {\mathbb{R}}^d:0<∥{\mathit{x}}^{\prime }-\mathit{x}∥⩽\ell \right\}$ is the neighborhood of material point $\mathit{x}\in \mathcal{B}$, and $∥·∥$ denotes the Euclidean norm. It is noted that such mesoscopic structure in essentially the simplest case of the n-point meso-structures (Chen et al. 2024 [26]).

Denote the set $\mathcal{T}=\left[0,t_{\infty }\right]$ with $t_{\infty }\in {\mathbb{R}}^{+}$ as the time of interest, and the map $\mathit{u}:\mathcal{B}\times \mathcal{T}\to {\mathbb{R}}^d$ be the displacement vector field. The deformation intensity measure of a material point-pair is expressed by the natural deformation (also called predominant elongation)

$$\lambda \left({\mathit{x}}^{\prime },\mathit{x},t\right)=\left[\mathit{u}\left({\mathit{x}}^{\prime },t\right)-\mathit{u}\left(\mathit{x},t\right)\right]·\nu ,\forall \left(\mathit{x}\prime ,\mathit{x},t\right)\in {\mathcal{M}}^{〈\ell 〉}\times \mathcal{T},$$

(5)

where $\nu =\left({\mathit{x}}^{\prime }-\mathit{x}\right)/∥{\mathit{x}}^{\prime }-\mathit{x}∥$ is the unit direction vector. The meso-scale damage can be described by the composition (Lu and Chen 2020 [1])

$$\omega \left({\mathit{x}}^{\prime },\mathit{x},t\right)=\varpi \circ \kappa \left({\mathit{x}}^{\prime },\mathit{x},t\right),\kappa =max\left[\underset{\tau \in \left[0,t\right]}{\sup }\lambda \left({\mathit{x}}^{\prime },\mathit{x},\tau \right)-{\lambda }_{\mathrm{c}},0\right],$$

(6)

where $\varpi :{\mathbb{R}}^{+}\cup \left\{0\right\}\to \left[0,1\right]$ is the meso-scale damage evolutionary rule required to be a monotone nondecreasing function, $\kappa$ is the quantity related to the deformation history of point-pair, and ${\lambda }_{\mathrm{c}}\in {\mathbb{R}}^{+}$ is the critical threshold.

Once the natural deformation of material point-pair exceeds the critical threshold, it leads to the geometric breaking of the influence domain, which can be measured by the topology damage (Lu and Chen 2020 [1])

$${\Omega }^{〈\ell 〉}\left(\mathit{x},t\right)={\int }_{{\mathcal{D}}_{\mathit{x}}^{〈\ell 〉}}\omega \left({\mathit{x}}^{\prime },\mathit{x},t\right)\mathrm{d}{\mathcal{L}}_{\phi }^d\left({\mathit{x}}^{\prime }\right),\forall \left(\mathit{x},t\right)\in \mathcal{B}\times \mathcal{T},$$

(7)

where ${\mathcal{D}}_{\mathit{x}}^{〈\ell 〉}={\mathcal{K}}_{\mathit{x}}^{〈\ell 〉}\cap \mathcal{B}$ is the influence domain, $\mathrm{d}{\mathcal{L}}_{\phi }^d=\phi \mathrm{d}{\mathcal{L}}^d$ is the weighted Lebesgue measure, and $\phi =1/{\mathcal{L}}^d\left({\mathcal{D}}_{\mathit{x}}^{〈\ell 〉}\right)$.

To establish the relationship between the geometry-based damage ${\Omega }^{〈\ell 〉}$ and the energy-based one $D^{〈\ell 〉}$, a continuous convex function $g:\left[0,1\right]\to \left[0,1\right]$ called the energetic degradation factor is employed as follows

$$D^{〈\ell 〉}\circ {\Omega }^{〈\ell 〉}+g\left({\Omega }^{〈\ell 〉}\right)=1\mathrm{s}.\mathrm{t}.D^{〈\ell 〉}⩾{\Omega }^{〈\ell 〉}.$$

(8)

Noting the definition of energy-based damage

$$D^{〈\ell 〉}=1-\psi \left(\mathit{ϵ},{\Omega }^{〈\ell 〉}\right)/{\psi }_0\left(\mathit{ϵ}\right),$$

(9)

i.e., the ratio of the current Helmholtz free energy density $\psi \left(\mathit{ϵ},{\Omega }^{〈\ell 〉}\right)$ to the elastic strain energy density ${\psi }_0\left(\mathit{ϵ}\right)$ for the appropriate strain tensor $\mathit{ϵ}$, and applying the Coleman-Noll procedure (Li et al. 2014 [43]) in continuum damage mechanics, the governing equation of NMMD model reads

$$\mathrm{div}\left[g\left({\Omega }^{〈\ell 〉}\right){\partial }_{\mathit{ϵ}}{\psi }_0\left(\mathit{ϵ}\right)\right]+\overline{\mathit{b}}=\mathit{0},\forall \left(\mathit{x},t\right)\in \mathcal{B}\times \mathcal{T}.$$

(10)

In Equation (10), the stress tensor is recognized as

$$\sigma =g\left({\Omega }^{〈\ell 〉}\right){\partial }_{\mathit{ϵ}}{\psi }_0\left(\mathit{ϵ}\right),$$

(11)

which is called the nonlocal macro-meso-scale consistent constitutive relationship, and $\overline{\mathit{b}}=\mathit{b}-\rho {\mathrm{d}}_{tt}\mathit{u}$ is the total body-force density including the inertial term $\rho {\mathrm{d}}_{tt}\mathit{u}$. The boundary conditions are subjected as

$$\sigma \mathit{n}=\overline{\mathit{t}},\forall \left(\mathit{x},t\right)\in {\partial }_{\overline{\mathit{t}}}\mathcal{B}\times \mathcal{T};\mathit{u}=\overline{\mathit{u}},\forall \left(\mathit{x},t\right)\in {\partial }_{\overline{\mathit{u}}}\mathcal{B}\times \mathcal{T},$$

(12)

where $\overline{\mathit{t}}$ and $\overline{\mathit{u}}$ are the prescribed surface force and displacement vector acting on the corresponding boundaries ${\partial }_{\overline{\mathit{t}}}\mathcal{B}$ and ${\partial }_{\overline{\mathit{u}}}\mathcal{B}$ with the outward unit normal vector $\mathit{n}$.

3.2. Application to PMMA-like Materials

The general formulation of NMMD model described in Subsection 3.1 is well suited to a wide range of applications, including brittle and ductile fracture. For PMMA-like materials investigated herein, we consider the following assumptions: (I) The solid is made up of an isotropic, homogeneous brittle material; (II) The displacement gradient $\nabla \mathit{u}$ is small, i.e., the strain tensor can be linearized as $\mathit{ϵ}=\left(\nabla \mathit{u}+{\nabla }^{\mathrm{T}}\mathit{u}\right)/2$, in which the superscript $\mathrm{T}$ denoted the transpose. In addition, the residual stress in the reference configuration vanishes; (III) The body is subjected to quasi-static loading such that the loading-rate effect and the inertial force ($\rho {\mathrm{d}}_{tt}\mathit{u}=\mathit{0}$) can be ignored.

With these assumptions, the NMMD model could be simplified. In particular, the elastic strain energy density in Equations (9) and (10) is the well-known Navier-Lamé potential function (Gurtin 1981 [44])

$${\psi }_0\left(\mathit{ϵ}\right)={\displaystyle \frac{E}{2\left(1+\mu \right)}}{∥\mathit{ϵ}∥}_{\mathrm{F}}^2+{\displaystyle \frac{E\mu }{2\left(1+\mu \right)\left(1-2\mu \right)}}{\left(\mathrm{tr}\mathit{ϵ}\right)}^2,$$

(13)

where E and $\mu$ are Young’s modulus and Poisson’s ratio respectively, ${∥·∥}_{\mathrm{F}}$ and $\mathrm{tr}\left(·\right)$ correspondingly indicate the Frobenius norm and the trace operator. The meso-scale damage evolutionary rule in Equation (6) is simultaneously reduced to the Heaviside step function, i.e.,

$$\varpi =H\left(\kappa \right)=\left\{\begin{array}{c}1,\kappa ⩾0,\hfill \\ 0,\kappa >0,\hfill \end{array}\right.$$

(14)

and for brittle material the energetic degradation factor in Equations (8) and (10) takes the following monomial

$$g\left({\Omega }^{〈\ell 〉}\right)={\left(1-{\Omega }^{〈\ell 〉}\right)}^p,p>1,$$

(15)

where the power p is the degradation parameter.

It should be noted that for the quasi-static cracking problems at each virtual time step (namely the load increment) the standard finite element procedure equipped with the local arc-length method (May and Duan 1997 [27]) is employed to obtain the complete equilibrium solutions of NMMD model, including descending branch and snap-back behaviors in the load-deflection response after attaining the peak-load.

4. Results and Discussions: Crack Onset of the Dog-Bone Specimens Without and With Hole

The dog-bone specimens described in Section 2 are revisited by the NMMD model, and the finite element mesh is shown in Figure 4. The basic mechanical properties of PMMA material listed in Table 1 are also adopted in the analysis, and all the specimens with and without hole are considered as plane stress problems.

4.1. Dog-Bone Specimen Without Hole

As can be seen in Section 3, besides the two ordinary elasticity constants (namely Young’s modulus E and Poisson’s ratio $\mu$), there are three additional parameters to be determined in the failure process analysis, i.e., the degradation parameter p, the critical threshold ${\lambda }_{\mathrm{c}}$ and the influence radius ℓ. In this subsection, we take the dog-bone specimen without hole as the first example to illustrate how these parameters influence the numerical results. Listed in Table 2 is the parameter setting for various degradation parameters (namely $p=4,5,6$), and the corresponding results of NMMD model are shown in Figure 5.

When equipped with the local arc-length method (May and Duan 1997 [27]), rather than the ordinary displacement control method, as can be seen, the NMMD model can capture the visible snap-back behavior in the post-peak regime of load-deformation curve. Meanwhile, it can be observed that the load-deformation curves obtained from NMMD model are almost the same everywhere as long as the degradation parameter is suitably large (e.g., $p⩾4$). In particular, the peak-load of all cases uniformly converges to $28.02\mathrm{kN}$, which is quite closed to the experimental result $P_{\mathrm{f}}^{exp}={\sigma }_{\mathrm{t}}\mathrm{WB}=28.2\mathrm{kN}$. These investigations highlight the ability of NMMD model to capture the whole process of brittle fracture without given pre-existing defects.

On that basis, we further investigate the calibration of critical threshold and influence radius, and the parameter setting is listed in Table 3. In these cases, the degradation parameter is consistently limited to a fixed constant (i.e., $p=5$) and the critical threshold varies in direct proportion to influence radius, namely ${\lambda }_{\mathrm{c}}\propto \ell$. The corresponding results of NMMD model are shown in Figure 6, in which the peak load of all cases almost equals to the experimental result when the proportion keeps a constant $2.38\times {10}^{-2}$. It is of interest to find that the proportion is exactly equivalent to the experimental critical strain, that is

$${\lambda }_{\mathrm{c}}/\ell ={\sigma }_{\mathrm{t}}/E={\mathit{ϵ}}_{\mathrm{t}},\forall \ell >0,$$

(16)

for cracking problems without singular stress concentration.

Besides, the load-deformation curve of all cases (though not shown here) are very similar to that as shown in Figure 5, and fairly closed to each other. Several damage profiles with regard to different values of influence radius are also shown in Figure 7. It seems that there are infinitely many pairs of critical threshold and influence radius satisfying the requirement of NMMD model in simulating the failure process of PMMA-like material.

However, in our previous analyses, e.g., Lu and Chen (2020 [1]), Ren et al. (2024 [3]) and Lu et al. (2024 [4]), the numerical results without exception point to the fact that the influence radius (i.e., characteristic length) should be deemed as a certain kind of physical or material property, rather than the mathematically fictitious parameter as argued in phase-field model (Wu 2017 [18]; Feng et al. 2021 [19]) and gradient damage model (Xue and Ren 2023 [21]). In other words, the influence radius in NMMD model should be identifiable and measurable in physical meaning. The further results and discussion are arranged in the next subsection.

4.2. Dog-Bone Specimen with Different Sizes of Circular Hole

Now let us turn to the cases with different sizes of circular hole. Shown in Figure 8 are the initial profiles before accessing to the main program of NMMD model, in which the radius of circular hole is considered to be the same size as experiments (Sapora et al. 2018 [22]), i.e., $R=\left\{0.25\mathrm{mm},0.5\mathrm{mm},1\mathrm{mm},2\mathrm{mm}\right\}$. Before analysis the topology damage at the pre-existing circular hole is set to one (i.e., the empty material in geometry). Figure 9 shows the numerical results (namely failure stress) with respect to various sizes of influence radius. As can be seen, the predictions of NMMD model would approximate that of MS criterion as the influence radius tends to zero and reach that of NS criterion when enlarging the influence radius. Among these results, the influence radius $\ell =0.38\mathrm{mm}$ is inferred from the phenomenological evidence that its result perfectly fits the experimental observations in regard to various sizes of hole.

From the above analyses, it can be seen that the characteristic length (i.e., influence radius) is identifiable and can be uniquely determined by the underlying physical experiment. With the material characteristic length, the NMMD model is capable of handling the failure stress of brittle fracture with and without pre-existing defect fairly well. Excepting the critical state of brittle failure, actually, the NMMD model could provide the simulation of the process sequences from linearly elastic deformation to complete failure, as shown in Figure 10 for the case of $R=2\mathrm{mm}$ with $\ell =0.38\mathrm{mm}$.

It is worth noting that crack has already initiated and developed to a certain extent (almost equaling to the size of influence radius; see Step I in Figure 10) when the load reaches its peak, which is consistent with many experimental observations, e.g., Romani et al. (2015 [28]) and Yin et al. (2022 [29]). Similar phenomenon also occurs in other sizes of circular hole. Figure 11 only displays their final failure morphology for the sake of non-repetition. As such, one clearly sees that the blame for incapacity of the conventional local criteria (including MS criterion and NS criterion) is attributable to the lack of material characteristic length.

4.3. Theoretical Determination of the Influence Radius

It is highlighted that the selection of influence radius in NMMD model (see Subsection 4.2) to numerically fit the experimental results is not for pleasing the eye but for aiming to assert the necessity and importance of involving characteristic length into the analysis of solid fracture. In the following, the theoretical determination of influence radius is further discussed.

Without loss of generality, the plane deformation problem is considered. Noting the fact that the crack has propagated slowly and steadily before the peak load is reached (see Figure 10), we therefore start with the cracked body. According to the recent result (Lu et al. 2024 [4]) that the crack topology represented by NMMD model always propagates along the direction of the maximum stress intensity factor of mode-I, which means the mixed mode brittle fracture can be equivalently transformed into the corresponding pure mode-I problem, hence we might as well focus on the mode-I loading, as shown in Figure 12.

Based on fracture mechanics, referring to the polar coordinate system $\left(r,\theta \right)$, the stress field of mode-I near crack tip is expressed as (Tada et al. 2000 [45], p.7)

$$\left\{\begin{array}{c}{\sigma }_{rr}={\displaystyle \frac{K_{\mathrm{I}}}{\sqrt{2\pi r}}}{\displaystyle \frac{5cos\left(\theta /2\right)-cos\left(3\theta /2\right)}{4}},\hfill \\ {\sigma }_{\theta \theta }={\displaystyle \frac{K_{\mathrm{I}}}{\sqrt{2\pi r}}}{\displaystyle \frac{3cos\left(\theta /2\right)+cos\left(3\theta /2\right)}{4}},\hfill \\ {\sigma }_{r\theta }={\displaystyle \frac{K_{\mathrm{I}}}{\sqrt{2\pi r}}}{\displaystyle \frac{sin\left(\theta /2\right)+sin\left(3\theta /2\right)}{4}},\hfill \end{array}\right.$$

(17)

where the coefficient $K_{\mathrm{I}}$ is the stress intensity factor of mode-I. Assume that the crack propagates with a small length scale $\Delta a$, as is well-known, the mode-I crack would extend along a straight line without kink-type deflection (i.e.,${\theta }_0=0^{\circ }$ ). For plane stress problems concerned in this paper, by applying the generalized Hooke’s law the corresponding strain field along the crack has

$$\left\{\begin{array}{c}{\mathit{ϵ}}_{rr}={\displaystyle \frac{K_{\mathrm{Ic}}}{E\sqrt{2\pi r}}}\left(1-\mu \right),\hfill \\ {\mathit{ϵ}}_{\theta \theta }={\displaystyle \frac{K_{\mathrm{Ic}}}{E\sqrt{2\pi r}}}\left(1-\mu \right),\hfill \\ {\mathit{ϵ}}_{r\theta }=0,\hfill \end{array}\right.$$

(18)

as $\Delta a\to 0^{+}$, where Irwin’s criterion (i.e., $K_{\mathrm{I}}=K_{\mathrm{Ic}}$) in linear elastic fracture mechanics (Anderson 2005 [46],p. 15) is employed for the extended crack.

According to Equation (5), right before crack propagation the maximum predominant elongation should occur at those material point-pairs which are perpendicular to the crack and pass through the crack tip, that is

$$\underset{\left({\mathit{x}}^{\prime },\mathit{x}\right)\in {\mathcal{M}}^{〈\ell 〉}\mathrm{s}.\mathrm{t}.\omega \left({\mathit{x}}^{\prime },\mathit{x}\right)=0}{\sup }\lambda \left({\mathit{x}}^{\prime },\mathit{x}\right)=\underset{0<s^{\prime }+s⩽\ell }{\sup }\left({\int }_0^{s^{\prime }}{\mathit{ϵ}}_{\theta \theta }\mathrm{d}r+{\int }_0^s{\mathit{ϵ}}_{\theta \theta }\mathrm{d}r\right),$$

(19)

where $s=∥\mathit{x}-\mathit{o}∥$ and $s^{\prime }=∥\mathit{x}\prime -\mathit{o}∥$ are the distances of material points $\mathit{x}$ and ${\mathit{x}}^{\prime }$ from crack tip $\mathit{o}$ (see Figure 12 for the illustration). Substituting Equation (18) into (19) and noting the evolution of meso-scale damage for brittle material, i.e., Equations (6) and (14), we have

$${\displaystyle \frac{2K_{\mathrm{Ic}}\left(1-\mu \right)\sqrt{\ell }}{E\sqrt{\pi }}}={\lambda }_{\mathrm{c}}.$$

(20)

Recalling the relationship between the critical threshold of and critical strain in Equation (16), furthermore, with combination of Equation (20) the influence radius in NMMD model can be ‌ theoretically determined by the solution, i.e.,

$$\ell ={\displaystyle \frac{4{\left(1-\mu \right)}^2}{\pi }}{\left({\displaystyle \frac{K_{\mathrm{Ic}}}{{\sigma }_{\mathrm{t}}}}\right)}^2.$$

(21)

This relationship is very similar to Irwin’s length (namely the size of crack-tip-yielding zone) in elastic-plastic fracture mechanics (Anderson 2005 [46], p. 61) with small scale yielding, that is, ${\ell }_{\mathrm{p}}={\left(K_{\mathrm{Ic}}/{\sigma }_{\mathrm{y}}\right)}^2/\pi$, where ${\sigma }_{\mathrm{y}}$ is the yield strength. It means that even in brittle material once crack growth is initiated, stresses at the crack tip would not be infinite. In this process, the irreversible breakage of meso-scale bonds, accompanied by the energetic degradation, leads to the further stress relaxation of crack-tip.

It is of interest to check that by substituting the material properties of PMMA material listed in Table 1 the theoretical influence radius ${\ell }_{\mathrm{th}}=0.378\mathrm{mm}$ in Equation (21) is very closed to the numerically identified one ${\ell }_{\mathrm{num}}=0.38\mathrm{mm}$ as discussed previously in Subsection 4.2. This result asserts again that the influence radius should not be construed to a mathematically fictitious parameter, but rather the physically-based material property.

5. Comparison with Other Modified Nonlocalized Criteria

Not coincidentally, in contrast to conventional local criteria (e.g., MS criterion and NS criterion), there are many modified nonlocalized criteria involving some kind of characteristic length, see Seweryn and Łukaszewicz (2002 [30]), Li and Zhang (2005 [24]) and Sapora et al. (2018 [22]) for the overview perspective. Among those modified criteria, the line method (LM; Novozhilov 1969 [31], Seweryn 1994 [32]), the quantized fracture mechanics (QFM; Seweryn and Łukaszewicz 2002 [30], Pugno and Ruoff 2004 [33]) and the finite fracture mechanics (FFM; Carpinteri et al. 2008 [34]) criteria are the representative works. In this section, discussions are devoted to the comparison between NMMD model and the above three modified criteria. To be consistent, examples of dog-bone specimen with a circular hole (Sapora et al. 2018 [22]) as previously discussed are still considered.

In LM criterion (Seweryn 1994 [32]), assuming that failure occurs when the averaged stress ahead of failure location over a segment of length ${\ell }_{\mathrm{LM}}$ reaches the specified strength, it states

$${\displaystyle \frac{1}{{\ell }_{\mathrm{LM}}}}{\int }_0^{{\ell }_{\mathrm{LM}}}{\sigma }_{\mathrm{LM}}\left(r\right)\mathrm{d}r={\sigma }_{\mathrm{t}},$$

(22)

and therefore,

$${\ell }_{\mathrm{LM}}={\displaystyle \frac{2}{\pi }}{\left({\displaystyle \frac{K_{\mathrm{Ic}}}{{\sigma }_{\mathrm{t}}}}\right)}^2,$$

(23)

where ${\sigma }_{\mathrm{LM}}$ is a certain kind of stress measure normal to the cracking direction.

Similarly, the QFM criterion (Pugno and Ruoff 2004 [33]) supposes that failure is achieved when the accumulated energy release rate G with regard to ${\ell }_{\mathrm{QFM}}$ (namely the crack extension) reaches the critical value $G_{\mathrm{c}}{\ell }_{\mathrm{QFM}}$, that is

$${\displaystyle \frac{1}{{\ell }_{\mathrm{QFM}}}}{\int }_a^{a+{\ell }_{\mathrm{QFM}}}G\left(s\right)\mathrm{d}s=G_{\mathrm{c}},$$

(24)

and therefore,

$${\ell }_{\mathrm{QFM}}={\displaystyle \frac{2}{c^2\pi }}{\left({\displaystyle \frac{K_{\mathrm{Ic}}}{{\sigma }_{\mathrm{t}}}}\right)}^2,c=\left\{\begin{array}{c}1\mathrm{for}\mathrm{center}\mathrm{crack},\hfill \\ 1.122\mathrm{for}\mathrm{edge}\mathrm{crack},\hfill \end{array}\right.$$

(25)

where a denotes the length of crack (already present or not).

The combination of LM criterion and QFM criterion gives birth to the FFM criterion (Carpinteri et al. 2008 [34])

$$\left\{\begin{array}{c}{\displaystyle \frac{1}{{\ell }_{\mathrm{FFM}}}}{\int }_0^{{\ell }_{\mathrm{FFM}}}{\sigma }_{\mathrm{LM}}\left(r\right)\mathrm{d}r={\sigma }_{\mathrm{t}},\hfill \\ {\displaystyle \frac{1}{{\ell }_{\mathrm{FFM}}}}{\int }_a^{a+{\ell }_{\mathrm{FFM}}}G\left(s\right)\mathrm{d}s=G_{\mathrm{c}},\hfill \end{array}\right.$$

(26)

where the characteristic length ${\ell }_{\mathrm{FFM}}$ is considered as an unknown quantity to be solved in Equation (26).

It is clear that the above three criteria are essentially nonlocalized criteria. In application to the failure analysis of dog-bone specimen with circular hole, the stress component ${\sigma }_{\mathrm{LM}}$ in LM criterion (22) and FFM criterion (26) is chosen as

$${\sigma }_{\mathrm{LM}}\left(r\right)={\sigma }_{\theta \theta }\left(r,\theta =\pm \pi \right)={\displaystyle \frac{\overline{\sigma }}{2}}\left(2+R^2/r^2+3R^4/r^4\right).$$

(27)

And with Irwin’s relationship (Anderson 2005 [46], p. 15) the energy release rate and fracture energy in QFM criterion (24) and FFM criterion (26) can be related to the stress intensity factor and fracture toughness respectively, i.e.,

$$G\left(s\right)={\displaystyle \frac{K_{\mathrm{I}}^2}{E}},G_{\mathrm{c}}={\displaystyle \frac{K_{\mathrm{Ic}}^2}{E}},$$

(28)

where the stress intensity factor is expressed as (Tada et al. 2000 [45], p.289)

$$K_{\mathrm{I}}\left(s\right)=\overline{\sigma }\sqrt{\pi s}F\left(\eta \right),\left\{\begin{array}{c}\eta =s/\left(s+R\right),\hfill \\ F\left(\eta \right)={\displaystyle \frac{1}{2}}\left(3-\eta \right)\left[1+1.243{\left(1-\eta \right)}^3\right].\hfill \end{array}\right.$$

(29)

By substituting the material properties of PMMA in Table 1, results of the above three criteria, including failure stress and the characteristic length, along with that of NMMD model and experiment (Sapora et al. 2018 [22]), are shown in Figure 13 and Figure 14.

As can be seen from Figure 13, LM criterion always provides the lower bound of the experimental failure stress for various sizes of circular hole, while the results of QFM criterion suffer the maximum percentage discrepancy (i.e., $12\%$) for the case without hole. Due to the coupling of LM and QFM criteria, the FFM criterion could effectively suppress the discrepancy at the case without hole suffered by QFM criterion, but not be of much assistance to improve the accuracy of cases with finite size of circular hole as predicted by LM criterion. Anyway, the FFM criterion pays a heavy price for overcoming the incongruence of QFM criterion, that is, the material constant stipulation of characteristic length, shown in Figure 14. Actually, as noticed that the LM, QFM and FFM criteria are all nonlocalized criteria. If it is further noticed that the energy release rate itself is a nonlocalized quantity (Lu et al. 2024 [4]), QFM and FFM are twice averaged, which is physically not convincingly reasonable. This might interpret why QFM and FFM criteria seem not as good as LM criterion. Compared with the above three modified criteria, the NMMD model with the theoretically-determined influence radius (namely the material constant $\ell \doteq 0.38\mathrm{mm}$) provides the best prediction for the mean values of experimentally-measured failure stress.

6. Revalidation: Mixed-Mode Fracture from Notch-Tip with Stress Singularity

Another experiment needed to be further investigated is the short bending beams with a pre-existing crack as reported in Mousavi et al. (2020 [35]) and Aliha et al. (2021 [36]). Figure 15 illustrates the details of geometry, loading and boundary condition. Specifically, the length, width and thickness of specimens were considered as: $L=36\mathrm{mm}$, $W=12\mathrm{mm}$ and $B=4\mathrm{mm}$. In each specimen, the pre-existing crack with length $a=6\mathrm{mm}$ makes an inclination angle $\alpha$ relative to the vertical direction, two simple-supports with the span $2S=12\mathrm{mm}$ were symmetrically placed at the bottom surface of the beam, and the vertical concentrated force P was loaded at the midpoint of the top surface. These specimens were also made of PMMA material, and mechanical properties are listed in Table 4, which are a bit different (namely the lower strength and the higher fracture toughness) from the sample in Table 1 as investigated in Sapora et al. (2018 [22]).

In the short bending beams, stress intensity factors of mode-I and mode-II are functions as follows (Mousavi et al. 2020 [35])

$$K_{\mathrm{I}}={\displaystyle \frac{P}{WB}}\sqrt{\pi a}{Y}_{\mathrm{I}}\left(a/W,S/W,\alpha \right),K_{\mathrm{II}}={\displaystyle \frac{P}{WB}}\sqrt{\pi a}{Y}_{\mathrm{II}}\left(a/W,S/W,\alpha \right),$$

(30)

where $Y_{\mathrm{I}}$ and $Y_{\mathrm{II}}$ are the geometry factors corresponding to crack length ratio $a/W$, loading span ratio $S/W$ and crack inclination angle $\alpha$. For the tested specimens, the variations of geometry factors respect to different inclination angles are shown in Figure 16. Five inclination angles (i.e., $\alpha =0^{\circ },{15}^{\circ },{26}^{\circ },{33}^{\circ },{39}^{\circ }$) were set to investigate cracking behaviors of different loading mode mixtites from pure mode-I ($\alpha =0^{\circ }$) to pure mode-II ($\alpha ={39}^{\circ }$).

By advocating Equations (16) and (21), the influence radius and the critical threshold in NMMD model are calibrated as $\ell =0.664\mathrm{mm}$ and ${\lambda }_{\mathrm{c}}=1.437\times {10}^{-2}\mathrm{mm}$ according to mechanical properties listed in Table 4. As can be seen, the influence radius for this example is slightly larger than that in the previous dog-bone specimens due to the lower strength and the higher fracture toughness. Numerical results including crack propagation path and peak load, together with experimental observations (Mousavi et al. 2020 [35]; Aliha et al. 2021 [36]) are shown in Figure 17 and Figure 18. It is found that the crack propagation paths from mode-I to mode-II loadings obtained by NMMD model are quite close to that observed in experiments. Meanwhile, the corresponding peak loads vary within the scope of experimental measurements. That is to say, the calibration procedure of NMMD model, i.e., Equations (16) and (21), though proposed for crack onset at non-singular stress concentration, is still suitable for problems of mixed mode fracture from notch-tip with stress singularity.

Such problems of short bending beam with a pre-existing crack, without a doubt, should belong to the scope of linear elastic fracture mechanics. However, as shown in Figure 19, both predictions of fracture toughness by the two celebrated local criteria, the maximum tangential stress (MTS) criterion (Erdogan and Sih 1963 [37]) and the maximum energy release rate (refers to as G) criterion (Wu 1978 [38]), cannot yet provide satisfactory consistency with the experimental results, especially for the mode-II dominated loading conditions (i.e., $K_{\mathrm{II}}/K_{\mathrm{I}}>1$ or $\alpha >{26}^{\circ }$; see Figure 16). Other commonly-used criteria, such as the principle of local symmetry criterion (Gol’dstein and Salganik 1974 [39]), the minimum strain energy density criterion (Sih 1974 [40]) and the maximum tangential strain criterion (Chang 1981 [41]), though not shown here, also encounter the similar discrepancies. As can be seen in Figure 19, the NMMD model with the theoretically-determined influence radius obviously provides better predictions of fracture toughness under mixed-mode loadings. On the other hand, comparison of kinking angle (see Figure 20) also confirms the assertion.

7. Concluding Remarks

In this work, we investigated the problems of crack onset at non-singular stress concentration and mixed-mode fracture from pre-existing notch-tip with stress singularity, in which predictions of the nonlocal macro-meso-scale consistent damage (NMMD) model, the conventional local criteria and the modified nonlocalized criteria, are discussed and compared with the corresponding experimental results. The main conclusions are summarized as follows:

(1) The conventional local criteria (e.g., maximum stress (MS) criterion and net stress (NS) criterion) cannot predict the experimentally-observed failure stress of specimens containing various sizes of circular defect with satisfaction due to the lack of involving characteristic length, whereas the NMMD model is capable of providing adequate critical loads consistent with experimental data for the size effect of defect.

(2) For brittle materials like PMMA all the parameters in NMMD model can be uniquely determined by standard specimens (e.g., the dog-bone specimen) with and without defect in company. The most important parameter in NMMD model, i.e., the influence radius (also called the characteristic length), is proved to be directly measurable in experiments, and theoretically it is related to the ratio of fracture toughness to tensile strength of the material.

(3) Compared to the modified nonlocalized criteria also involving characteristic length (e.g., line method criterion, quantized fracture mechanics criterion and finite fracture mechanics criterion), the NMMD model is not only able to predict the critical loads in higher accuracy with regarding to experimental results but also able to quantitatively provide crack morphologies from initiation to propagation.

(4) The calibration procedure proposed for NMMD model is suitable for both brittle fracture problems of crack onset at non-singular stress concentration and mixed mode cracking from notch-tip with stress singularity. The results show a remarkable improvement in prediction of fracture toughness and kinking angle in comparison to the celebrated local criteria in linear elastic fracture mechanics, such as the maximum tangential stress criterion and the maximum energy release rate criterion.

These investigations imply that the NMMD model constitutes a promising tool for prediction of failure load and crack evolution in brittle materials. Even so, much more should be done in the future, including extensions of this work to failure analysis of quasi-brittle and heterogeneous materials, e.g., concrete, asphalt and rock mass.