This paper investigates the effect of strain rate on the scaling behavior of dynamic tensile strength of quasibrittle structures. The theoretical framework is anchored by a rate-dependent finite weakest link model. The model involves a rate-dependent length scale, which captures the transition from localized damage to diffused damage with an increasing strain rate. As a result, the model predicts a rate- and size-dependent probability distribution function of the nominal tensile strength. The transitional behavior of the strength distribution directly leads to the rate and size effects on the mean and standard deviation of the tensile strength. The model is verified by a series of stochastic discrete element simulations of dynamic fracture of aluminum nitride specimens. The simulations involve a set of geometrically similar specimens of various sizes subjected to a number of different strain rates. Both random microstructure geometry and fracture properties are considered in these simulations. The simulated damage pattern indicates that an increase in the strain rate results in a more diffusive cracking pattern, which supports the theoretical formulation. The simulated rate and size effects on the mean and standard deviation of the nominal tensile strength agree well with the predictions by the rate-dependent finite weakest link model.

## Introduction

Understanding the influence of loading rate on the structural strength plays an important role in the design of engineering structures against impact loading. Over the years, brittle heterogeneous materials, such as ceramics and composites, have extensively been used in many protective technologies including human and vehicle armors [1,2]. These engineering applications have stimulated substantial advances in computational modeling of dynamic fracture of brittle and quasibrittle materials. Continuum-based numerical modeling of dynamic quasibrittle fracture often relies on some micromechanics-based constitutive models, which aim to capture the essential dynamic failure mechanisms at the meso and microscales. In many of these models [3–5], dynamic fracture mechanics was used to describe the rate-dependent microcrack growth and the macroscopic structural response was obtained through some homogenization schemes. By employing the kinetics of microcrack growth, the macroscopic constitutive models could also take into account the effect of random distribution of preexisting flaws [4]. Recent studies further incorporated the effect of material comminution into the continuum constitutive models [6–8], which led to more realistic simulations of structural failure under high velocity impact.

Besides continuum-based modeling, discrete modeling has also been widely used to simulate the fracture behavior of brittle and quasibrittle materials subjected to high strain-rate loading rates [9–14]. The discrete models aim to capture explicitly the microstructural features of the material, such as size distribution and morphology of the inhomogeneities, and randomly distributed preexisting flaws. These models have been used in stochastic computations for predicting the probability distribution of material strength. Meanwhile, these models have also been adopted to simulate material fragmentation under blast loading [11,15], and it was shown that the convergence rate of the overall calculated energy dissipation can be improved significantly if one considers random microstructural features.

For the purpose of design of engineering structures against impact loading, one key design parameter is the dynamic structural strength. It is well known that many brittle and quasibrittle materials, such as ceramics and concrete, exhibit a rate enhancement of structural strength [3–5,12,16,17]. This rate enhancement has been extensively investigated by using the aforementioned computational models. The rate effect on structural strength can be attributed to the dynamics of interactive microcrack growth. In several recent studies, mesoscale discrete models have been used to simulate the tensile fracture of concrete and ceramic materials subjected to high strain-rate loading [12,18,19]. These simulations indicated that the overall failure mechanism strongly depends on the applied strain rates. At low strain rates, the failure pattern at peak load can be characterized by the formation of a localized damage zone, in which many microcracks nucleate to a macroscopic crack, whereas at high strain rates a diffused damage pattern is observed at peak load.

In the actual design process, another critical consideration is the size effect on the structural strength since many full-scale designs are based on the extrapolation of reduced-scale laboratory test results. Over the past three decades, extensive research has focused on the size dependence of the static strength of quasibrittle structures. The salient feature of quasibrittle structures is that, under quasi-static loading, the size of the fracture process zone, which represents a material length scale, is not negligible compared to the overall structure size [20]. The existence of this length scale gives rise to an intricate size effect on the structural strength. So far, two fundamental types of size effect have been well understood: (1) energetic–statistical size effect (type 1 size effect), which applies to structures without any preexisting notches [21,22], and (2) pure energetic size effect (type 2 size effect), which applies to structures with a large preexisting notch [21,23]. It has been acknowledged that both size effects manifest the underlying size dependence of the overall failure mechanism of the structure, which transitions from a quasi-plastic failure manner to a brittle failure manner as the structure size increases [20,22].

As mentioned earlier, there exists clear evidence that, for a given specimen, the failure mechanism is strongly influenced by the strain rate. Therefore, it is expected that the applied strain rate would have an influence on the scaling behavior of the structural strength. However, most of the existing scaling models were developed only for the quasi-static loading condition, in which the rate effect is absent. Therefore, despite the major advances in computational modeling of dynamic failure of brittle and quasibrittle structures, there is a lack of understanding of rate-dependent scaling of structural strength. As an attempt to address this fundamental problem, this study presents a rate-dependent size effect model of the nominal tensile strength of quasibrittle specimens. By using aluminum nitride (AlN) as a model system, the model is verified by a series of stochastic discrete computations.

## Rate-Dependent Finite Weakest Link Model

Consider a series of geometrically similar rectangular quasibrittle specimens, which is loaded under displacement control with a prescribed strain rate (Fig. 1(a)). We define the nominal tensile strength, *σ*_{N}, as the maximum average stress that the specimen can sustain. The dependence of the mean nominal tensile strength $\sigma \xafN$ on the specimen size *D* arises from two sources: (1) energetic size effect due to the nonlinear softening material behavior, and (2) statistical size effect due to the spatial variability of material strength, fracture energy, and microstructure. Previous studies have shown that the small and intermediate-size regimes of this type of size effect can be captured by some nonlinear deterministic models, such as the cohesive crack model, crack band model, and nonlocal models, whereas the large-size regime can be described by using the theory of Weibull statistics [22].

In a series of recent studies, it was shown that this size effect can alternatively be explained by a finite weakest link model, which is anchored by a statistical description of structural failure [24–26]. In this model, the structure is statistically modeled as a chain of material representative volume elements (RVEs) (Fig. 1(b)). This model implies that the structure will fail under controlled loads once one RVE is fully damaged. In other words, in the present context, the RVE is taken as the smallest material element whose failure triggers the failure of the entire structure, which is different from the RVE definition used in the classical homogenization theory, e.g., Refs. [27] and [28]. It should be emphasized that the weakest link model is a statistical representation of the damage localization mechanism, which is a fundamental failure mechanism of quasibrittle structures. Meanwhile, it is noted that this model is designed solely for the purpose of calculating the statistics of the structural strength, but not for determining the overall load–displacement response of the structure.

where $P1(\sigma N,\u03f5\u02d9)$ = strength distribution of either one RVE or the entire specimen, whichever has a smaller size, $n(\u03f5\u02d9)=A/[l0(\u03f5\u02d9)]2,\u2009A=DL$ = area of the specimen, $\u27e8x\u27e9=max(x,0)$ = Macaulay bracket, and $l0(\u03f5\u02d9)=RVE$ size, which is a function of the strain rate. Equation (1) is written in a two-dimensional setting since it is reasonable to consider that the material damage grows simultaneously throughout the thickness of the specimen.

Several recent studies [12,18] have investigated the dynamic tensile fracture behavior of ceramic and concrete materials by considering specimens similar to that depicted in Fig. 1(a). These simulation results consistently showed that, with an increasing strain rate, the damage pattern becomes more diffusive. Therefore, we expect that the RVE size will increase with the applied strain rate. For a given specimen size, an increase in the applied strain rate could possibly make the RVE larger than the specimen size, i.e., $n(\u03f5\u02d9)<1$. In this case, Eq. (1) reduces to $Pf(\sigma N,\u03f5\u02d9)=P1(\sigma N,\u03f5\u02d9)$. In other words, $P1(\sigma N,\u03f5\u02d9)$ becomes the strength distribution of the entire specimen. This implies that the weakest link representation of structural failure vanishes since there is no damage localization in the specimen.

where *m*_{0} and *s*_{0} are the Weibull modulus and the scale parameter of the Weibull tail, respectively, and *μ*_{G0} and *δ*_{G0} are the mean and standard deviation of the Gaussian core if considered extended to −*∞*; *r*_{f0} is a scaling parameter required to normalize the grafted distribution such that $P1(\u221e)=1,\u2009\sigma gr0$ = grafting stress, and $Pgrm0=(\sigma gr0/s0)m$ = grafting probability. Here, the subscript and superscript “0” denote the case of quasi-static loading. The continuity of the probability density function (pdf) at the grafting point requires that $[dP1/d\sigma N]\sigma gr0+=[dP1/d\sigma N]\sigma gr0\u2212$. Therefore, we need four of these parameters to define the grafted distribution of RVE strength.

*l*

_{0 s}and henceforth we refer to the material element of size

*l*

_{0 s}as the static RVE. We can then define a parameter

*n*

_{b}, which is equal to the equivalent number of static RVEs in the RVE or in the entire specimen, whichever is smaller, under high strain-rate loading, i.e.,

In this study, we propose to determine the distribution function $P1(\sigma N,\u03f5\u02d9)$ for dynamic strength by using a fiber-bundle model (Fig. 2(b)), in which each element (also called “fiber”) represents one static RVE. Since *n*_{b} may not be an integer, we consider the bundle contains $\u2308nb\u2309$ fibers, where $\u2308nb\u2309=$ the least integer that is greater than or equal to *n*_{b}. In this bundle, $\u2308nb\u2309\u22121$ fibers have the same cross-sectional area *A*_{0}, and the remaining fiber has a cross-sectional area of $(nb+1\u2212\u2308nb\u2309)A0$. The fiber-bundle model statistically represents the damage redistribution mechanism of the RVE under high strain-rate loading (or of the entire specimen if $n(\u03f5\u02d9)<1$). Since the autocorrelation length of the random strength field is considered to be smaller than the size of the static RVE, the strength of each fiber can be treated as statistically independent.

The investigation of strength statistics of fiber bundles has a long history. The pioneering work of Daniels [31] dealt with brittle bundles with an equal load sharing rule, and a recursive equation was derived for the corresponding strength statistics. A considerable amount of efforts were subsequently devoted to studying brittle fiber bundles with other types of phenomenological load sharing rules [32,33]. A more realistic approach is to determine the load sharing rule from the mechanical behavior of the fiber. Here, we consider a fiber bundle, which consists of $\u2308nb\u2309$ fibers connected by two rigid plates (Fig. 2(b)). As discussed earlier, each fiber represents a static RVE. For quasibrittle materials, the mechanical behavior of static RVE will exhibit some degree of strain softening. On the other hand, we also note that the static RVE should behave in a quasi-plastic manner since it contains only a few material inhomogeneities [24]. This quasi-plastic failure behavior is also manifested by the dominant Gaussian distribution of RVE strength as indicated by Eq. (2*b*) [24–26]. Therefore, the stress–strain curve of the static RVE must have a gentle softening behavior even under quasi-static loading. As the applied strain rate increases, the static RVE would experience more diffusive microcracking, which alleviates the degree of strain softening [12,18,34,35].

where $\sigma i(i=1,\u2026,\u2308nb\u2309)$ are the random strengths of individual fibers. The strength distribution of each fiber is predominantly Gaussian with a very short power-law tail (Eqs. (2*a*) and (2*b*)). Previous studies have shown that, for such a strength distribution of fibers, the strength distribution of the bundle will also contain a power-law tail, and meanwhile the core of the strength distribution of the bundle is Gaussian [24]. Therefore, it is clear that the previously proposed grafted distribution function (Eqs. (2*a*) and (2*b*) provides the correct functional form for the strength distribution function $P1(\sigma N,\u03f5\u02d9)$ though the values of the statistical parameters must now depend on the applied strain rate.

As mentioned earlier, the grafted strength distribution function can be determined by any four statistical parameters in Eqs. (2*a*) and (2*b*). In this study, we incorporate the strain-rate dependence into the strength distribution function $P1(\sigma N,\u03f5\u02d9)$ through the following four parameters:

- The Weibull modulus
*m*measures the exponent of the power-law tail of the bundle. By using the series expansion method, it was shown that the power-law exponent of the bundle is equal to the sum of the power-law exponents of the tail distributions of the individual fibers [24]. In fact, this property also holds for brittle and softening bundles, where the randomness of the stress–strain curve of the fibers is described by an affine transformation [25,26,36]. Meanwhile, the strain rate may also affect the Weibull modulus of the static RVE itself, which is described by an empirical function $f(\u03f5\u02d9)$ to be determined later. Following these considerations, we can express the Weibull modulus of the distribution function $P1(\sigma N,\u03f5\u02d9)$ as$m(\u03f5\u02d9)=m0f(\u03f5\u02d9)\u2308nb(\u03f5\u02d9)\u2309$(5) - The grafting probability
*P*_{gr}determines the extent of the power-law tail. Previous studies have shown that, as the number of fibers increases, the power-law tail of the strength distribution of the bundle shortens at the rate of $(Pgr,f/\u2308nb\u2309)\u2308nb\u2309$ (*P*_{gr,f}is the grafting probability of strength distribution of an individual fiber) [24]. By further considering the potential influence of the strain rate on the grafting probability of the static RVE, the rate-dependent grafting probability is written aswhere $Pgr0g(\u03f5\u02d9)$ describes the rate-dependent grafting probability of the static RVE. Based on Eq. (6), it is clear that the Weibull tail of the RVE strength distribution is shortened drastically as$Pgr(\u03f5\u02d9)\u2248[Pgr0g(\u03f5\u02d9)/\u2308nb(\u03f5\u02d9)\u2309]\u2308nb(\u03f5\u02d9)\u2309$(6)*n*_{b}increases. This means that when the strain rate is sufficiently high, Eq. (5) becomes practically unimportant. - Parameter
*μ*_{G}is approximately equal to the mean strength of the plastic bundle. The power-law tail of the strength distribution of each fiber is too short to influence the mean behavior. Based on the statistics of the weighted sum of independent Gaussian variables (Eq. (4)), it is clear that*μ*_{G}is equal to the mean strength of individual fibers, which leads to the following expression:where $p(\u03f5\u02d9)$ is the rate enhancement function of the mean strength of the static RVE.$\mu G(\u03f5\u02d9)=\mu G0p(\u03f5\u02d9)$(7) - Parameter
*δ*_{G}measures the standard deviation of the bundle strength since the power-law tail of the distribution function has a minimal influence on the second moment of statistics. The statistics of the weighted sum of independent Gaussian variables gives$\delta G(\u03f5\u02d9)=\delta G0q(\u03f5\u02d9)nb[\u2308nb\u2309\u22121+(nb+1\u2212\u2308nb\u2309)2]1/2$(8)

where function $q(\u03f5\u02d9)$ describes the influence of strain rate on the standard deviation of the strength of the static RVE.

*a*) and (2

*b*)), we obtain the rate-dependent strength distribution function $P1(\sigma N,\u03f5\u02d9)$. Based on the finite weakest link model (Eq. (1)), we can then calculate the mean dynamic strength of the specimen as

By considering specimens of different sizes, we can determine the effect of specimen size on the statistics of dynamic strength. Together, with the aforementioned rate dependence, Eqs. (9) and (10) describe the combined size and rate effects on the first and second moments of the statistics of nominal tensile strength of the specimen.

## Stochastic Discrete Element Computational Model

Motivated by the recent success in discrete modeling of dynamic fracture [12,14,18,37–39], we perform a series of stochastic discrete element computations of dynamic fracture of ceramic specimens to verify the proposed analytical model. In the computational model, the specimen is represented by a set of interconnected discrete rigid bodies. A set of nuclei is first randomly placed in the domain, where the mutual distance of two adjacent nuclei is controlled to be approximately equal to the grain size *l*_{min} of the material. Once these nuclei are created, the Voronoi tessellation is used to discretize the domain (Fig. 3(a)). The location of the *i*th nuclei is denoted by $xi$. Each nucleus has three translational and three rotational degrees-of-freedom. Due to the imposition of the length scale *l*_{min}, each Voronoi body can be considered to represent a material grain. The contact surface (henceforth referred to as facet) between the bodies represents the grain boundary (Fig. 3(b)). The macroscopic nonlinear behavior of the material is determined by the mesoscale constitutive behavior of the facets.

In this study, the formulation of the constitutive behavior of the facet follows the lattice discrete particle model, which was originally developed for concrete materials [40–45]. Since this study focuses primarily on macroscopic tensile failure, we adopt a simplified version of the lattice discrete particle model, which is briefly summarized here. The detailed description of the model can be found in Ref. [46].

*k*, and $Ak(k=i,j)$ = transformation matrices, which can be written as

where $x1c,x2c,x3c=$ coordinates of the facet centroid, and $x1k,x2k,x3k=$ coordinates of the nucleus of body *k*.

*l*= distance between the nuclei of bodies, and $n,l,m$ = unit vectors in the normal and two orthogonal tangential directions, respectively. The corresponding traction vector transmitted over the facet can then be related to the strain vector by using damage mechanics, i.e.,

_{ij}*E*

_{0}= elastic stiffness,

*ω*= damage parameter, α = constant, and

*t*

_{n},

*t*

_{l},

*t*

_{m}= tractions in the directions of $n,l$ and

**m**, respectively. The damage parameter is expressed in terms of the equivalent stress and strain, i.e.,

*f*

_{t}= tensile strength,

*f*

_{s}= shear strength, $\psi =tan\u22121[en/\alpha (em2+el2)],\u2009\eta =fs/ft$, and

*μ*= constant. The history variable

*χ*is expressed as

where $emax=max(en2)+\alpha max(em2+el2)$, in which the maxima is calculated for the entire loading history, and $\psi 0=\u2212arctan(\zeta \alpha )$. For the present numerical simulations, the facets mainly experience tensile–shear damage. Therefore, it is not necessary to cover the full scenario of compressive damage.

*K*represents the initial slope of the nonlinear branch of the constitutive law [42], which is given by

*K*

_{t}and

*K*

_{s}are the values of

*K*for pure tensile and shear loading, respectively.

*K*

_{t}and

*K*

_{s}are expressed as

where $lt=E0Gt/ft2,\u2009ls=E0Gs/fs2$, and *G*_{t}, *G*_{s} = mode I and II fracture energies. It should be pointed out that Eq. (19) involves two material length scales *l*_{t} and *l*_{s} so that the energy expended for the fracture of a unit area of facet is independent of the distance between the nuclei. Once the tractions on the facets are determined, the forces and moments of the nuclei can be calculated based on the principle of virtual work [46].

It is worthwhile to mention that the present constitutive model itself is rate-independent, and therefore the rate effect on the macroscopic behavior of the specimen is largely due to inertia and its influence on the interaction of facet failures. Such a modeling approach has been adopted in several previous discrete simulations of dynamic fracture of ceramic materials [10,12,19,47], which were able to capture the rate-dependent failure behavior reasonably well. Meanwhile, it has also been shown that the rate dependence of the constitutive behavior plays a secondary role in the overall failure behavior under high strain-rate loading [48,49], which is of interest to the present study.

*E*

_{0}in the normal direction, (2) tangential to normal stiffness ratio

*α*, (3) mesoscopic tensile and shear strengths

*f*

_{t},

*f*

_{s}, (4) mesoscopic mode I and II fracture energies

*G*

_{t},

*G*

_{s}, and (5) constants

*μ*and

*ζ*. It is well known that the properties of quasibrittle materials usually exhibit a considerable degree of spatial variability. For studying the failure behavior, the most relevant parameters are strengths and fracture energies [29]. In this study, we characterize the spatial random distribution of mesoscopic strengths and fracture energies by assigning a single random field $h(x)$:

*l*

_{t}and

*l*

_{s}are deterministic constants. The underlying probability distribution function governing the random field $h(x)$ follows the Gauss–Weibull grafted distribution (Eqs. (2

*a*) and (2

*b*)) with a mean value equal to one. The spatial autocorrelation of the field $h(x)$ is characterized by a Gaussian function, i.e.,

where $\delta x=|x\u2212x\u2032|$ and *l*_{a} = autocorrelation length. The generation of the random field $h(x)$ consists of three steps [50]: (1) a standard Gaussian random field $h(x)$ is generated on a regular square grid with a spacing of *l*_{a}/4, (2) the value of $h(x)$ on each facet is determined by using the optimal linear estimation method [51], and (3) the standard Gaussian field generated on the facets is converted to the Gauss–Weibull random field $h(x)$ by using the isoprobabilistic transformation. Figure 4 presents a typical realization of the random field $h(x)$.

The aforementioned computational model is numerically implemented with the implicit Newmark method [52,53]. The equation of motion is solved incrementally in an implicit scheme, which is unconditionally stable. The mass matrix is constructed by using the full inertia tensor of Voronoi bodies. The moment of inertia of these polyhedral bodies is calculated by dividing them into tetrahedra [54]. Though a stable numerical solution scheme is chosen, a small time-step is needed to capture the dynamic failure of the facet. For each simulation case, several trial simulations are performed to determine a desirable time-step, which yields a consistent result.

## Simulations of Dynamic Tensile Failure

of Aluminum Nitride Specimens

The aforementioned stochastic computational model is applied to simulate the dynamic tensile failure of AlN specimens. We consider square specimens of different in-plane sizes *D* = *L* = 50 *μ*m, 100 *μ*m, 200 *μ*m, 400 *μ*m, 800 *μ*m (Fig. 1(a)), whereas the out-of-plane thickness is set constant *b* = 10 *μ*m. The average in-plane grain size *l*_{min} of AlN material is taken as 6 *μ*m [12]. Each specimen is subjected to a constant strain-rate loading, which is applied by imposing a nonuniform velocity field described as $vx1=x1\u03f5\u02d9$ (*x*_{1} denotes the horizontal position of each nucleus). In order to investigate the rate-dependent behavior, we consider a wide range of applied strain rates, i.e., $\u03f5\u02d9=1,1000,2500,5000,7500,104,3\xd7104,5\xd7104,105,2\xd7105/s$.

We choose the following mesoscale material parameters for AlN: *E*_{0} = 530 GPa, *α* = 0.17, $f\xaft=150\u2009MPa,\u2009f\xafs=3f\xaft,\u2009G\xaft=2\u2009Jm\u22122,\u2009G\xafs=16G\xaft,\u2009\mu =0.2$, and *ζ* = 0.95. The elastic parameters (*E*_{0} and *α*) are determined to match the macroscopic elastic response of AlN reported in Ref. [12]. The inelastic mesoscale material parameters are chosen so that the model predicts the quasi-static tensile strength of a 50 *μ*m square AlN specimens being about 120 MPa, which is similar to the published results [12]. It is admitted that the mixed-mode loading scenarios will be needed to better calibrate the coupling between tensile and shear damage. However, such a detailed calibration procedure is not necessary for the present purpose. For the random field $h(x)$, we consider that it has an autocorrelation length *l*_{a} of 24 *μ*m, and the underlying Gauss–Weibull distribution function has a mean value of one, a coefficient of variation of 20%, a Weibull modulus of 30 and a grafting probability of 5 × 10^{−4}.

*σ*

_{a}along

*x*

_{1}-direction (Fig. 1(a)). For the present discrete system, we first calculate the fabric stress tensor of a single Voronoi body

*k*as [55,56]

*V*= volume of the

_{k}*k*th body,

*n*= number of facets of the body,

_{k}*A*≡ surface area of the

_{p}*p*th facet, $tp$ = traction vector on the

*p*th facet expressed in the global coordinate system, $xcp$ = position vector of the centroid of the

*p*th facet, and $xk$ = position vector of the nucleus of the body. The nominal stress

*σ*

_{a}can then be calculated as a volume average of the fabric stress components in the

*x*

_{1}-direction, i.e.,

where *V* = volume of the specimen, and *N* = number of Voronoi bodies in the specimen.

Based on the simulations, we obtain a set of random responses of the relationship between the nominal stress *σ*_{a} and the average strain *ϵ*_{a}, where *ϵ*_{a} = *u*/*D* and *u* = displacement measured at the free edge of the specimen. Figure 5 presents the simulated average *σ*_{a}−*ϵ*_{a} curves for different specimen sizes and strain rates. It can be seen that, for a given specimen size, the average nominal strength $\sigma \xafN$, which is the maximum value of the average nominal stress, increases with the applied strain rate. This rate enhancement is primarily due to the inertia effect. Meanwhile, we also observe that, as the strain rate increases, the postpeak softening behavior becomes less pronounced, which indicates that the specimen exhibits a more ductile behavior. This rate dependence can directly be observed from the simulated damage pattern at peak load. As a demonstration, Fig. 6 shows the damage patterns of specimens of *D* = 400 *μ*m at peak load for different applied strain rates. It is seen that, for a given specimen size, the damage is fairly localized at low strain rates while the specimen exhibits a diffused damage pattern as the strain rate increases.

It should be pointed out that the aforementioned brittle-to-ductile transition is based on the current simulation setup with the prescribed initial velocity field. Recent computational studies of ceramics and concrete specimens with the same setup showed similar rate-dependent failure patterns [12,18]. Several experiments on the rate-dependent behavior of concrete also showed that the tensile stress–strain curve exhibits a more gentle postpeak softening behavior at high strain rates [34,35]. However, there is experimental evidence showing that the compressive stress–strain response of brittle composites and ceramics could become more brittle at higher strain rates [57–59]. The setups of these experiments are different from the present simulations. We may expect that the rate dependence of failure behavior of quasibrittle materials could depend strongly on the specimen geometry and test setup. This is consistent with the well accepted fact that, for quasibrittle structures, brittleness is influenced not only by the material properties but also by the structure geometry and loading configuration [20,22,26].

We can also compare the simulated postpeak responses of specimens of different sizes for a given strain rate. It is seen that, at low strain rates, the slope of the postpeak softening curve becomes much steeper as the specimen size increases. This implies that the specimen experiences a more brittle failure as the specimen size increases. Such a size-dependent failure behavior at low strain rates has been well documented for many quasibrittle materials, such as concrete, rock, composites, and ceramics [20,22]. As the strain rate increases, the size dependence of the postpeak behavior starts to diminish. This is consistent with the present simulation results, which show that, at high strain rates, all the specimens considered in this study exhibit a diffused cracking pattern.

## Comparison With Rate-Dependent Weakest Link Model

The foregoing discussion reveals qualitatively the combined rate and size effects on the failure behavior of the simulated specimens. The quantitative description of these effects can be best presented in terms of the size effect on the mean structural strength. In this section, we discuss the simulated size effect on both the mean nominal strength and the standard deviation of the nominal strength. These simulation results are used to compare with the rate-dependent finite weakest link model.

Figure 7 presents the simulated size effect curves of the mean nominal strength at different applied strain rates. The most salient feature is that the mean size effect diminishes as the strain rate increases. At high strain rates, the entire specimen experiences diffused cracking as the peak load is reached. In this case, there is no characteristic length governing the nominal strength. Since the failure is quasi-plastic, the size effect on the nominal strength must be absent [22,26]. As the strain rate decreases, the damage pattern transitions from diffused cracking to localized cracking. The size of the localized damage zone represents a characteristic length scale, which leads to a nonpower-law form of the mean size effect curve.

We now use the rate-dependent finite weakest link model (Eq. (9)) to fit this set of size effect curves. It is seen from Fig. 7 that the model can match well the simulation results for all strain rates. Based on this fitting, we determine the rate effect on the RVE size and the strength distribution of the static RVE. Figure 8 presents the effect of strain rate on the RVE size. It is clear that the RVE size *l*_{0} increases with the strain rate. It is seen that the RVE size is almost a constant for the relatively low strain rates, and increases mildly for a range of intermediate strain rates 3 × 10^{4}–5 × 10^{4}/s. When the strain rate exceeds 5 × 10^{4}/s, the RVE size increases significantly with the strain rate. In fact, for the high strain-rate regime studied here, the analysis indicates that all the specimens are smaller than the RVE. Therefore, the actual RVE size is not known. At high strain rates, the RVE can be purely a mathematical abstraction if the specimen does not exhibit damage localization. In this case, we may set $l0(\u03f5\u02d9)$ to be a very large number. However, introducing the concept of RVE for the present theoretical framework is essential because it captures the rate dependence of the failure statistics, which transitions from a weakest link model signifying a damage localization mechanism to a fiber-bundle model representing a diffused damage mechanism as the strain rate increases.

It is clear that the concept of RVE is central to the present weakest link model. Within the framework of weakest link statistics, RVE is a mathematical concept. In a recent study [60], it has been shown that the RVE size can be best determined from the optimum fitting of the mean size curve of nominal structural strength, and a relationship between the RVE size and basic material length scales (the Irwin characteristic length and crack band width) was established for the case of quasi-static loading. Based on the results of this recent study, we expect that the RVE size will increase with the average grain size. Meanwhile, the RVE size will also increase with an increase in the mesoscale fracture energies. This implies that the RVE size is dependent on the loading configuration. For example, if the applied loading induces more shear-dominant fractures along the grain boundaries, based on the material properties used in the present simulation, the specimen will exhibit a more ductile failure behavior and therefore the RVE size will increase.

Figure 9 presents the influence of the applied strain rate on the strength statistics of the static RVE. We observe that the Weibull modulus increases with the applied strain rate (i.e., $f(\u03f5\u02d9)>1$ in Eq. (5)) and the grafting probability decreases with the strain rate (i.e., $g(\u03f5\u02d9)<1$ in Eq. (6)). It is noted that, when the strain rate is higher than 7500/s, the power-law tail of the strength distribution of the static RVE is too short to be determined for the range of the specimen sizes considered here. The observed rate dependence of the strength statistics of static RVE can be attributed to the fact that, as the strain rate increases, the static RVE experiences more dense microcracking. In the context of the present stochastic discrete element model, microcracking is represented by the failure of grain boundaries (facets). The increase in the Weibull modulus can be qualitatively explained by a fishnet statistics model, which was recently developed to analyze a fishnet structure for predicting the probability distribution of static strength of nacreous imbricated lamellar materials. Though the specimens dealt in this study are very different from the fishnet structure, the same mathematical framework can be applied to both cases.

where *P*_{s,r} = probability that *r* number of facets have experienced damage after stress *σ*_{N} is applied to the static RVE, and *N* = total number of facets. For the practical use of Eq. (25), we just need to retain several dominant terms (*r* = 0,…, *k*), which contribute to the overall failure probability. The aforementioned rate effect on the microcrack density indicates that *k* increases with the applied strain rate. To calculate the individual probability component *P*_{s,r} of Eq. (25), one would need to determine the load distribution among the facets. A realistic calculation of the load distribution pattern is equivalent to performing the nonlinear stochastic discrete element simulations, which is computationally intensive. Nevertheless, the recent results of the fishnet analysis indicate that if the strength statistics of each facet has a power-law tail, then the strength distribution of the system will also exhibit a power-law tail. The power-law exponent increases in proportion to the dominant number of damaged facets (i.e., *k* in Eq. (25)) [61]. This explains the observed increase in the Weibull modulus with an increasing strain rate.

The decrease in the grafting probability indicates that the strength distribution of the static RVE approaches a Gaussian distribution. This can be physically related to the observed increasing number of damaged facets, which implies that the failure becomes more ductile. Since more facets contribute to the failure of the static RVE, we expect that the standard deviation would decrease considerably. However, it is also noted that the random mesostructure causes a considerable variability of the stress field, and this variability is further enhanced by the inertia effect. Therefore the stress field becomes more random at high strain rates, which leads to a larger standard deviation of the overall nominal strength. This effect counteracts the aforementioned decrease in the standard deviation due to the increasing number of damaged facets. As shown in Fig. 9, the calibration results show that the overall standard deviation of the strength of static RVE decrease mildly with an increasing strain rate.

*σ*

_{N}can be written as

where Γ(*x*) = Eulerian gamma function. Equation (26) indicates that, at low strain rates, the asymptotic size effect on the standard deviation is very weak since *m* is usually a large number.

*n*

_{b}= (

*D*/

*l*

_{0s})

^{2}. Therefore, the statistics of the nominal strength of the specimen is modeled by a plastic bundle, and the scaling of the standard deviation is given by Eq. (8). It can be easily shown that, as

*n*

_{b}increases (i.e.,

*n*

_{b}> 4), Eq. (8) is approximately equivalent to

We note that the size effect depicted by Eq. (27) is much stronger than that predicted by the Weibull distribution. From Fig. 10, we see that Eq. (27) agrees well with the simulated size effect on the standard deviation at high strain rates.

The foregoing analysis elucidates the size and rate dependence of the functional form of strength distribution. At low strain rates, the strength distribution of the specimen transitions from a Gaussian distribution modified by a far-left Weibull tail to a Weibull distribution as the specimen size increases. At high strain rates, for the specimens considered in this study, the strength distribution is essentially Gaussian with a constant mean value but a decreasing standard deviation. This transitional behavior is well captured by the present model.

It is also worthwhile to discuss the transition between the low and high strain rates. There could exist a narrow range of intermediate strain rates, where the increase in the RVE size is not significantly large but the grafting probability of the RVE is considerably reduced (Eq. (6)). Consider a very large specimen containing a sufficient number of RVEs subjected to such intermediate strain-rate loading. Since the power-law tail of RVE strength becomes very short, for the failure risk of practical interest (*P*_{f} ≈ 10^{−6}), the strength distribution of the specimen can be considered as a long chain of Gaussian elements, which would converge to a Gumbel distribution as an intermediate asymptote [62]. Such a Gumbel asymptotic distribution is not seen in the finite weakest link model for static loading since the Weibull tail is not sufficiently short.

Finally, it should be emphasized that the present rate-dependent weakest link model is formulated for specimens shown in Fig. 1(a). Therefore, the model is not intended to predict the dynamic strength of specimens under general impact loading, where the strain rate could vary both spatially and temporally. However, this model has important implications for the stochastic finite element (FE) simulation of dynamic quasibrittle fracture. If we consider that each Gauss point represents a material element, whose size is related to the mesh size, the present model indicates that the input probability distribution of tensile strength of each Gauss point would depend on both the strain rate and mesh size. It has recently been demonstrated that, for quasi-static loading, the influence of the mesh size on the input probability distribution of tensile strength is an essential consideration for mitigating in order to mitigate the mesh sensitivity in stochastic FE simulations of tensile fracture of quasibrittle structures [63].

In this regard, the present rate-dependent weakest link model can be used to determine the probability distribution of the tensile strength of each Gauss point for stochastic FE simulations. Meanwhile, as demonstrated in this study, the weakest link model can be calibrated by a set of mesoscale discrete simulations. In this way, the present model provides a link between the variability of mesoscale structural properties and the probability distribution of the macroscopic tensile strength. This link not only provides physical insights into the rate-dependent tensile strength distribution but also could improve the robustness of the stochastic FE simulation of dynamic quasibrittle fracture.

## Conclusions

This study shows that the strain rate has a profound influence on the scaling of the nominal tensile strength of quasibrittle structures. This rate dependence arises from the effect of the strain rate on the failure behavior of the structure. The discrete element simulations show that the structure exhibits a more diffused cracking pattern as the strain rate increases, which explains the observed diminishing mean size effect at high strain rates. Meanwhile, the rate-dependent failure behavior also leads to different scaling behaviors of the standard deviation of the nominal strength at different strain rates. It is shown that, for a given strain rate, the standard deviation decreases with an increasing specimen size. At high strain rates, this decreasing trend becomes stronger and it exhibits a power-law scaling behavior. The result also indicates that, for a given specimen size, an increasing strain rate leads to a reduction in the standard deviation of the nominal strength.

It is shown that the simulated rate and size effects on the mean and standard deviation of the nominal tensile strength can be well captured by a rate-dependent finite weakest link model. The model employs a rate-dependent length scale, which physically represents the transition from damage localization to diffused damage with an increasing strain rate. The model directly predicts the combined rate and size effects on the probability distribution of the nominal strength. At low strain rates, the strength distribution varies from a predominant Gaussian distribution to a Weibull distribution as the specimen increases, while at high strain rates, the strength distribution is Gaussian for a large range of specimen sizes. The present rate-dependent finite weakest link model provides an analytical tool to model the rate and size effects on the probability distribution of nominal tensile strength, which could be implemented into the stochastic FE simulation of dynamic fracture of quasibrittle structures.

## Funding Data

The Solid Mechanics Program of the U.S. Army Research Office (Grant No. W911NF-15-1-0197).

National Science Foundation (Grant No. NSF/CMMI-775 1361868).

Undergraduate Research Apprenticeship Program of the U.S. Army Research Office.

Czech Science Foundation (Project No. 15-19865Y).