# Research on the Size Effect of Unstable Fracture Toughness by the Modified Maximum Tangential Stress (MMTS) Criterion

Authors: Jian Li; Zhao-Wen Du; Zhong-Ping Guo; De-Chun Ai
Journal: Advances in Civil Engineering
 (2020)

---

## Abstract

This paper investigates the unstable fracture toughness of specimens of different heights using the double-K model for three-point bending tests on notched concrete beams. It is shown that unstable fracture toughness exhibits a significant size effect. The modified maximum tangential stress (MMTS) criterion is used to explain the size effect of unstable fracture toughness. The MMTS criterion considers the higher order terms of the Williams series expansion of the stress field. The results show that the MMTS criterion can reasonably estimate unstable fracture toughness. It is recommended that the minimum height of the specimen be 200 mm when three-point bending tests on notched beams are used to determine unstable fracture toughness.


---

## Body

## 1. Introduction

Fracture toughness is an important property of fracture resistance, especially for quasibrittle materials such as rock and concrete. It is involved in building structure safety, energy mining, and other fields. Therefore, many researchers and organizations have measured fracture toughness with various specimen configurations and methods. For example, the American Society for Testing and Materials (ASTM) [1] proposed ASTM E1820. Three-point bending tests on notched beams were used by the International Union of Laboratories and Experts in Construction Materials (RILEM) [2] to determine the fracture energy of concrete. The International Society for Rock Mechanics (ISRM) recommended four specimens to determine the model I fracture toughness of rock materials (chevron bend (CB) and short rod (SR) specimens [3], cracked chevron-notched Brazilian disc (CCNBD) specimens [4], and semicircular bend (SCB) [5]). As a material property, fracture toughness should be related only to the material itself, regardless of the method used. However, the opposite is true. There are significant differences in fracture toughness between different sample sizes and shapes. Tutluoglu and Keles [6] used CCNBD and SCB specimens to measure the fracture toughness of rocks and found them to be different. Aliha et al. and Ayatollahi et al. [7–9] used different shapes of specimens to measure the fracture toughness of rocks. It was found that fracture toughness was significantly affected by the geometry and loading conditions of the specimens.This is because a fracture process zone (FPZ) exists at the crack tip, and due to its influence, the linear elastic fracture mechanics (LEFM) are no longer suitable for quasibrittle materials [10]. Therefore, many nonlinear fracture mechanics theories have been proposed to explain the fracture problem, for example, the fictitious crack model, the two-parameter model (TPM), the effective crack model, the size effect model (SEL) [11–15], and the double-K model [16, 17]. Among these theoretical methods, one of the advantages of the double-K model is that the test process is simple. There is no need for repeated loading and unloading during the test. Only the load versus crack mouth opening-displacement (P-CMOD) curve needs to be measured to calculate the critical effective crack length and the unstable fracture toughness. The unstable fracture toughness Kun, considering the influence of the fracture process zone, can better represent the fracture resistance property of the material. However, the size effect of Kun was also observed by experiments. Xu and Zhang [18] pointed out that a relatively stable fracture toughness value can be measured only when the height of the test specimens is greater than 200 mm. The values measured for specimens with a height less than 200 mm are significantly lower. Kumar and Barai [19] studied double-K parameters by numerical simulation. They concluded that the fracture parameters obtained from the double-K fracture models were influenced by the specimen size. Ruiz et al. [20] researched the effect of size on the double-K parameters in concrete fracture through experiments and numerical methods. They concluded that Kun showed a moderate increase with size.However, the effect of size on the double-K parameters is not sufficient, especially for small size test pieces. In addition to FPZ, the high-order terms of the crack tip stress field (in Williams’s series expansions) are also important factors influencing the stress intensity factor. For example, Aliha et al. and Ayatollahi et al. [7–9] considered that fracture toughness is affected by the size and geometry of the specimen, mainly due to the influence of high-order terms of the crack tip stress field, and the experimental results were analyzed by the modified maximum tangential stress (MMTS) criterion. They showed that the MMTS criterion can well estimate the apparent fracture toughness.It should be noted that the actual fracture process of the crack is very complicated. Not only will the I model fracture but also the II model or even the III model may occur. For example, Golewski [21–23] pointed out that when the crack develops in a complex stress state, II model fracture, or even III model fracture or their combination, such as I-II model and I-III model, may occur. In addition, due to the high heterogeneity of the internal structure of rock and concrete, the random distribution of microcracks, and the uncertainty of the load, these make the fracture process of quasibrittle materials very complicated [24, 25].The purpose of this paper is to study the size effect of the unstable fracture toughnessKun (one of the double-K parameters) using the MMTS criterion. The P-CMOD curves required to calculate unstable fracture toughness are derived from three-point bending tests on notched beams of Hoover et al. [26, 27] in 2012. They showed that the measured unstable fracture toughness (Kun) values were significantly correlated with the specimen size. The size effect of unstable fracture toughness can be explained by the MMTS criterion.


## 2. Theory

### 2.1. Stress Intensity Factor and CMOD

Since a span-to-depth ratio (S/H) of 2.17 was used in the three-point bending tests on notched beams instead of the usual ratio of 4, the formulas of stress intensity factor and crack mouth opening displacement (CMOD) cannot be directly found in The Stress Analysis of Cracks Handbook [28]. Guinea et al. [29] provide general expressions for the stress intensity factor and CMOD, but they claimed that these expressions are applicable to span-to-depth ratios larger than 2.5. It is therefore necessary to derive the expressions for stress intensity factor and CMOD when S/H = 2.17.The stress intensity factor can be obtained directly in Abaqus finite element software using the contour integral method [30]. Therefore, Abaqus 2016 was employed to simulate three-point bending tests on notched beams with S/H = 2.17 and different initial crack ratios.The finite element models of the notched beam were simulated with plane stress conditions, as shown in Figure1. The detailed parameters are as follows (see Table 1): geometric parameters of the model are length (L) = 2.4 m, span (S) = 2.17, thickness (B) = 0.04, depth (H) = 1, and initial crack ratio (α0) = 0.2, 0.3, 0.4, 0.5, and 0.6; material properties are elastic modulus (E)  =  1 and Poisson’s ratio (v)  =0.25. Taking α = 0.2 as an example, the total number of nodes: 8601, the total number of elements: 8417, and the type of elements: CPS4, as shown in Figure 1. The results of the final simulation are shown in Table 1.Figure 1
Model of three-point bending tests on notched beams in Abaqus (L, length; H, depth; S, span; a, notch length), where the span-depth ratio (S/H) is 2.17.Table 1
Parameters for simulating three-point bending test of notched beam in Abaqus and results of simulation.
L (m)S (m)B (m)H (m)a0 (m)α0E (Pa)vP (Pa)CMOD (m)K (Pam)Y(α)V(α)2.42.170.0410.20.210.254.58E − 074.19371E − 052.70E − 051.621.412.42.170.0410.30.310.254.09E − 076.53230E − 053.14E − 051.721.632.42.170.0410.40.410.253.35E − 078.82119E − 053.32E − 051.932.022.42.170.0410.50.510.252.35E − 071.05873E − 043.16E − 052.342.772.42.170.0410.60.610.251.59E − 071.28538E − 043.05E − 053.054.15B, thickness; E, elastic modulus; v, Poisson’s ratio; CMOD, crack mouth opening displacement.The stress intensity factor and CMOD expressions for notched beams are as follows:(1)K=3PS2BH2Yαa,where Y(α) is the normalized stress intensity factor, and it is a function of the relative crack length ratio α = a/H, depending on the geometry of the specimen and the loading conditions, and(2)CMOD=6PSaEBH2Vα,where V(α) is a dimensionless shape function depending on the span-to-depth ratio and relative crack length.According to equations (1) and (2), the finite element simulation results can be used to calculate the corresponding Y(α) and V(α) values, and they are listed in Table 2.Table 2
Y(α) and V(α) values obtained by the finite element method compared with values obtained by formulas of Guinea et al. and relative errors.
α0Y(α)V(α)Y(α) (Guinea et al.)V(α) (Guinea et al.)V(α) (equation (5))Relative errorY(α) (%)Relative errorV(α) (%)Relative errorV(α) (equation (5)) (%)0.21.62161.44340.83921.37801.44341.242.104.750.31.71921.59321.13241.60081.59321.252.060.470.41.92872.01891.52042.02882.01891.940.270.490.52.33892.77542.09852.78912.77541.480.780.490.63.05364.14453.05764.21704.14451.541.651.72BothY(α) and V(α) can be fitted to a polynomial for relative crack length α, as follows:(3)Yα=6.1283α4+6.2444α3−3.399α2+1.0911α+1.4796,R=1,(4)Vα=37.677α4−20.794α3+6.1234α2+0.7097α+1.1261,R=1,where R is the correlation coefficient.In order to verify the validity of the method, theY(α) and V(α) values were compared with the general expression values according to Guinea et al. [29] and are drawn in Figures 2 and 3. The corresponding values of S/H = 4 are also given in the figures.Figure 2
Y(α) versus relative crack length α (square red points indicate values obtained from formulas of Tada et al. for S/H = 4).Figure 3
V(α) versus relative crack length α (square red points indicate values obtained from formulas of Tada et al. for S/H = 4).It can be seen that the relative errors of theY(α) and V(α) values calculated by the finite element method and the general expression value according to Guinea et al. are at most 2.10%, which is highly accurate.SinceV(α) (equation (4)) is a polynomial of a higher degree, its inverse function is not easy to solve, so it can be further simplified as equation (5). The relative error of the V(α) value is also analyzed and listed in Table 2. In addition to α = 0.2, the relative error of other values is not more than 2%, which can meet requirements:(5)CMOD=PEB15.9978tan2π2α+2.0698.According to equation (5), the relative length of the crack expression can be rewritten as(6)α=2πarctanEBCMOD15.9978P−0.1294.It can be seen that a span-to-depth ratio of 4 corresponds to largerY(α) and V(α) values than S/H = 2.17 in Figures 2 and 3. The geometry has a significant effect on Y(α) and V(α). The stress intensity factor and CMOD expressions obtained by the finite element method are very consistent with the results of Guinea et al., who provide general expressions. It is noted that equations (3) and (4) are only valid for S/H = 2.17 and relative crack length ratio α in the range of 0.2–0.6.


### 2.2. MMTS Criterion

For the stress field at the tip of the model I crack (Figure4), stress σyy (in the Y direction) in the orthogonal coordinate system can be expressed as follows, according to Williams’s series expansions:(7)σy=∑n=1∞n2Anrn/2−1−n2−2+−1ncosn2−1θ+n2−1cosn2−3θ,where r=x2+y2 and θ=arctany/x (as illustrated in Figure 4), n is the order of the term, and An represents the coefficients of the terms in Williams's series expansion.Figure 4
Stress components at crack tip shown in the rectangular coordinate system.Since the coefficientsAn are affected by the geometry of the specimen and the loading conditions, they can be dimensionless for convenience according to the following equation:(8)A1=KIc2π,An∗=a0.5n−1σNAn,where KIc is the exact value of the stress intensity factor and σN is the nominal stress. An∗ is the dimensionless coefficient of the terms in Williams’s series expansion. For the three-point bending test of notched beam, nominal stress σN can be defined by(9)σN=3PS2BH2.According to the maximum tangential stress (MTS) criterion [31], only a singular term of the stress field (A1) is considered. In recent years, many researchers have found that the high-order terms of the stress field have a significant effect on the stress intensity factor for quasibrittle materials such as rock and concrete. The modified maximum tangential stress (MMTS) criterion is based on an improvement of the MTS criterion, taking into account not only the singular term of the crack tip stress field but also the high-order terms. For example, Aliha et al. and Ayatollahi et al. [7–9] found that the fracture parameters estimated by the MMTS criterion agree well with the experimental results.For model I crack,θ = 0, σyy at the extension line of the crack tip is(10)σyyθ=0=A1r1+3A3A1r+5A5A1r2+….According to equation (3), the above equation can be further written as(11)σyy=KIc2πr1+3A3∗A1∗ra+5A5∗A1∗ra2+….Before fracture failure of quasibrittle materials, there are inelastic zones (fracture process zones) at the tip of the original crack, as shown in Figure5. Due to the FPZ, linear elastic fracture mechanics are no longer applicable. If the FPZ is equivalent to a stress-free crack Δa, which is equivalent to moving the tip of the original crack O to the right O′, linear elastic fracture mechanics continue to apply [32, 33].Figure 5
Stress distributions for fictitious crack model (Hillerborg et al.) and equivalent elastic stress field.The double-K model proposes a convenient method to determine the critical effective crack lengthac = a0 + Δa (Δa is the critical effective crack growth) according to the P-CMOD loading curve. Specifically, the peak load Pmax and critical CMODC are first obtained according to the measured P-CMOD loading curve. Then, the critical relative crack length ratio αc can be reversed according to equation (2) of CMOD. Finally, the unstable fracture toughness Kun can be calculated according to the LEFM (equation (1)).If the stress intensity factor is calculated using the MMTS criterion, the equivalent plastic zone lengthΔr needs to be estimated. Irwin [34] proposed a method for estimating Δr; since the FPZ is equivalent to Δa, the influence of stress relaxation is no longer considered.For plane strain, Δr is calculated as follows:(12)Δr=12πKICft21−2v2,where KIc is the fracture toughness, ft is the tensile strength, and v is Poisson’s ratio.Critical stressσyyc is taken to be tensile strength ft according to the fictitious crack model [12]. Plugging σyyc=ft into equation (11), it can be written as(13)KIc=2πΔrft1+3A3∗/A1∗Δr/a+5A5∗/A1∗Δr/a2+….Ignoring the high-order terms, the apparent fracture toughnessKif is(14)Kif=2πΔrft,and then, Kif/KIc is derived as(15)KifKIc=1+3A3∗A1∗Δrac+5A5∗A1∗Δrac2+…,where Δr is the plastic zone length and ac = a0 + Δa is the critical effective crack length.Therefore, the ratioKif/KIc is affected by the high-order terms and Δr/ac.The ratioKif/KIc can be calculated if the high-order terms (A3∗/A1∗, A5∗/A1∗, …), and critical effective crack length ac are determined.Ayatollahi and Nejati [35] proposed a method to determine higher order terms by the finite element method, called the finite element overdeterministic (FEOD) method. This method is simple and reliable, and the dimensionless coefficients of any higher order terms (A1∗, A3∗, …) can be determined if sufficient; accurate displacement field data are available. The calculation process and detailed steps can be found in [35].According to the results in Table3, the critical effective crack length ratio in this test is about 0.4. Therefore, the coefficients of the higher order terms in Williams’s series expansion are calculated by the FEOD method for the three-point bending tests on notched beams with relative crack length ratio α = 0.4, and the highest term is A19. The required displacement field data can be obtained with the finite element simulation results used to derive the stress intensity factor in Section 2.1. The dimensionless coefficients of the higher order terms are shown in Table 4.Table 3
Results calculated from the double-K model.
LabelS (m)B (m)H (m)a0/Hac/HΔa/HPmax (N)Kun (Pa)Aa1.09E + 003.88E − 025.00E − 013.29E − 013.7548E − 014.66E − 021.08E + 041.4615E + 06Ba4.69E − 013.89E − 022.15E − 013.36E − 014.1117E − 017.54E − 026.40E + 031.4596E + 06Ca2.02E − 014.13E − 029.30E − 023.45E − 014.2080E − 017.53E − 023.58E + 031.1967E + 06Da8.65E − 024.04E − 024.00E − 023.65E − 014.4377E − 017.90E − 021.80E + 031.0001E + 06Table 4
Dimensionless coefficient ratioA3∗/A1∗–A19∗/A1∗ for relative crack length ratio α = 0.4.
Ratio of high-order termsValueA3∗/A1∗–0.32644A5∗/A1∗–0.22786A7∗/A1∗0.228546A9∗/A1∗–0.39727A11∗/A1∗0.916394A13∗/A1∗–1.86374A15∗/A1∗3.640217A17∗/A1∗–9.39833A19∗/A1∗23.50735It is worth noting that ratioAn∗/A1∗ begins to alternate between positive and negative, starting with A5∗/A1∗.


## 2.1. Stress Intensity Factor and CMOD

Since a span-to-depth ratio (S/H) of 2.17 was used in the three-point bending tests on notched beams instead of the usual ratio of 4, the formulas of stress intensity factor and crack mouth opening displacement (CMOD) cannot be directly found in The Stress Analysis of Cracks Handbook [28]. Guinea et al. [29] provide general expressions for the stress intensity factor and CMOD, but they claimed that these expressions are applicable to span-to-depth ratios larger than 2.5. It is therefore necessary to derive the expressions for stress intensity factor and CMOD when S/H = 2.17.The stress intensity factor can be obtained directly in Abaqus finite element software using the contour integral method [30]. Therefore, Abaqus 2016 was employed to simulate three-point bending tests on notched beams with S/H = 2.17 and different initial crack ratios.The finite element models of the notched beam were simulated with plane stress conditions, as shown in Figure1. The detailed parameters are as follows (see Table 1): geometric parameters of the model are length (L) = 2.4 m, span (S) = 2.17, thickness (B) = 0.04, depth (H) = 1, and initial crack ratio (α0) = 0.2, 0.3, 0.4, 0.5, and 0.6; material properties are elastic modulus (E)  =  1 and Poisson’s ratio (v)  =0.25. Taking α = 0.2 as an example, the total number of nodes: 8601, the total number of elements: 8417, and the type of elements: CPS4, as shown in Figure 1. The results of the final simulation are shown in Table 1.Figure 1
Model of three-point bending tests on notched beams in Abaqus (L, length; H, depth; S, span; a, notch length), where the span-depth ratio (S/H) is 2.17.Table 1
Parameters for simulating three-point bending test of notched beam in Abaqus and results of simulation.
L (m)S (m)B (m)H (m)a0 (m)α0E (Pa)vP (Pa)CMOD (m)K (Pam)Y(α)V(α)2.42.170.0410.20.210.254.58E − 074.19371E − 052.70E − 051.621.412.42.170.0410.30.310.254.09E − 076.53230E − 053.14E − 051.721.632.42.170.0410.40.410.253.35E − 078.82119E − 053.32E − 051.932.022.42.170.0410.50.510.252.35E − 071.05873E − 043.16E − 052.342.772.42.170.0410.60.610.251.59E − 071.28538E − 043.05E − 053.054.15B, thickness; E, elastic modulus; v, Poisson’s ratio; CMOD, crack mouth opening displacement.The stress intensity factor and CMOD expressions for notched beams are as follows:(1)K=3PS2BH2Yαa,where Y(α) is the normalized stress intensity factor, and it is a function of the relative crack length ratio α = a/H, depending on the geometry of the specimen and the loading conditions, and(2)CMOD=6PSaEBH2Vα,where V(α) is a dimensionless shape function depending on the span-to-depth ratio and relative crack length.According to equations (1) and (2), the finite element simulation results can be used to calculate the corresponding Y(α) and V(α) values, and they are listed in Table 2.Table 2
Y(α) and V(α) values obtained by the finite element method compared with values obtained by formulas of Guinea et al. and relative errors.
α0Y(α)V(α)Y(α) (Guinea et al.)V(α) (Guinea et al.)V(α) (equation (5))Relative errorY(α) (%)Relative errorV(α) (%)Relative errorV(α) (equation (5)) (%)0.21.62161.44340.83921.37801.44341.242.104.750.31.71921.59321.13241.60081.59321.252.060.470.41.92872.01891.52042.02882.01891.940.270.490.52.33892.77542.09852.78912.77541.480.780.490.63.05364.14453.05764.21704.14451.541.651.72BothY(α) and V(α) can be fitted to a polynomial for relative crack length α, as follows:(3)Yα=6.1283α4+6.2444α3−3.399α2+1.0911α+1.4796,R=1,(4)Vα=37.677α4−20.794α3+6.1234α2+0.7097α+1.1261,R=1,where R is the correlation coefficient.In order to verify the validity of the method, theY(α) and V(α) values were compared with the general expression values according to Guinea et al. [29] and are drawn in Figures 2 and 3. The corresponding values of S/H = 4 are also given in the figures.Figure 2
Y(α) versus relative crack length α (square red points indicate values obtained from formulas of Tada et al. for S/H = 4).Figure 3
V(α) versus relative crack length α (square red points indicate values obtained from formulas of Tada et al. for S/H = 4).It can be seen that the relative errors of theY(α) and V(α) values calculated by the finite element method and the general expression value according to Guinea et al. are at most 2.10%, which is highly accurate.SinceV(α) (equation (4)) is a polynomial of a higher degree, its inverse function is not easy to solve, so it can be further simplified as equation (5). The relative error of the V(α) value is also analyzed and listed in Table 2. In addition to α = 0.2, the relative error of other values is not more than 2%, which can meet requirements:(5)CMOD=PEB15.9978tan2π2α+2.0698.According to equation (5), the relative length of the crack expression can be rewritten as(6)α=2πarctanEBCMOD15.9978P−0.1294.It can be seen that a span-to-depth ratio of 4 corresponds to largerY(α) and V(α) values than S/H = 2.17 in Figures 2 and 3. The geometry has a significant effect on Y(α) and V(α). The stress intensity factor and CMOD expressions obtained by the finite element method are very consistent with the results of Guinea et al., who provide general expressions. It is noted that equations (3) and (4) are only valid for S/H = 2.17 and relative crack length ratio α in the range of 0.2–0.6.


## 2.2. MMTS Criterion

For the stress field at the tip of the model I crack (Figure4), stress σyy (in the Y direction) in the orthogonal coordinate system can be expressed as follows, according to Williams’s series expansions:(7)σy=∑n=1∞n2Anrn/2−1−n2−2+−1ncosn2−1θ+n2−1cosn2−3θ,where r=x2+y2 and θ=arctany/x (as illustrated in Figure 4), n is the order of the term, and An represents the coefficients of the terms in Williams's series expansion.Figure 4
Stress components at crack tip shown in the rectangular coordinate system.Since the coefficientsAn are affected by the geometry of the specimen and the loading conditions, they can be dimensionless for convenience according to the following equation:(8)A1=KIc2π,An∗=a0.5n−1σNAn,where KIc is the exact value of the stress intensity factor and σN is the nominal stress. An∗ is the dimensionless coefficient of the terms in Williams’s series expansion. For the three-point bending test of notched beam, nominal stress σN can be defined by(9)σN=3PS2BH2.According to the maximum tangential stress (MTS) criterion [31], only a singular term of the stress field (A1) is considered. In recent years, many researchers have found that the high-order terms of the stress field have a significant effect on the stress intensity factor for quasibrittle materials such as rock and concrete. The modified maximum tangential stress (MMTS) criterion is based on an improvement of the MTS criterion, taking into account not only the singular term of the crack tip stress field but also the high-order terms. For example, Aliha et al. and Ayatollahi et al. [7–9] found that the fracture parameters estimated by the MMTS criterion agree well with the experimental results.For model I crack,θ = 0, σyy at the extension line of the crack tip is(10)σyyθ=0=A1r1+3A3A1r+5A5A1r2+….According to equation (3), the above equation can be further written as(11)σyy=KIc2πr1+3A3∗A1∗ra+5A5∗A1∗ra2+….Before fracture failure of quasibrittle materials, there are inelastic zones (fracture process zones) at the tip of the original crack, as shown in Figure5. Due to the FPZ, linear elastic fracture mechanics are no longer applicable. If the FPZ is equivalent to a stress-free crack Δa, which is equivalent to moving the tip of the original crack O to the right O′, linear elastic fracture mechanics continue to apply [32, 33].Figure 5
Stress distributions for fictitious crack model (Hillerborg et al.) and equivalent elastic stress field.The double-K model proposes a convenient method to determine the critical effective crack lengthac = a0 + Δa (Δa is the critical effective crack growth) according to the P-CMOD loading curve. Specifically, the peak load Pmax and critical CMODC are first obtained according to the measured P-CMOD loading curve. Then, the critical relative crack length ratio αc can be reversed according to equation (2) of CMOD. Finally, the unstable fracture toughness Kun can be calculated according to the LEFM (equation (1)).If the stress intensity factor is calculated using the MMTS criterion, the equivalent plastic zone lengthΔr needs to be estimated. Irwin [34] proposed a method for estimating Δr; since the FPZ is equivalent to Δa, the influence of stress relaxation is no longer considered.For plane strain, Δr is calculated as follows:(12)Δr=12πKICft21−2v2,where KIc is the fracture toughness, ft is the tensile strength, and v is Poisson’s ratio.Critical stressσyyc is taken to be tensile strength ft according to the fictitious crack model [12]. Plugging σyyc=ft into equation (11), it can be written as(13)KIc=2πΔrft1+3A3∗/A1∗Δr/a+5A5∗/A1∗Δr/a2+….Ignoring the high-order terms, the apparent fracture toughnessKif is(14)Kif=2πΔrft,and then, Kif/KIc is derived as(15)KifKIc=1+3A3∗A1∗Δrac+5A5∗A1∗Δrac2+…,where Δr is the plastic zone length and ac = a0 + Δa is the critical effective crack length.Therefore, the ratioKif/KIc is affected by the high-order terms and Δr/ac.The ratioKif/KIc can be calculated if the high-order terms (A3∗/A1∗, A5∗/A1∗, …), and critical effective crack length ac are determined.Ayatollahi and Nejati [35] proposed a method to determine higher order terms by the finite element method, called the finite element overdeterministic (FEOD) method. This method is simple and reliable, and the dimensionless coefficients of any higher order terms (A1∗, A3∗, …) can be determined if sufficient; accurate displacement field data are available. The calculation process and detailed steps can be found in [35].According to the results in Table3, the critical effective crack length ratio in this test is about 0.4. Therefore, the coefficients of the higher order terms in Williams’s series expansion are calculated by the FEOD method for the three-point bending tests on notched beams with relative crack length ratio α = 0.4, and the highest term is A19. The required displacement field data can be obtained with the finite element simulation results used to derive the stress intensity factor in Section 2.1. The dimensionless coefficients of the higher order terms are shown in Table 4.Table 3
Results calculated from the double-K model.
LabelS (m)B (m)H (m)a0/Hac/HΔa/HPmax (N)Kun (Pa)Aa1.09E + 003.88E − 025.00E − 013.29E − 013.7548E − 014.66E − 021.08E + 041.4615E + 06Ba4.69E − 013.89E − 022.15E − 013.36E − 014.1117E − 017.54E − 026.40E + 031.4596E + 06Ca2.02E − 014.13E − 029.30E − 023.45E − 014.2080E − 017.53E − 023.58E + 031.1967E + 06Da8.65E − 024.04E − 024.00E − 023.65E − 014.4377E − 017.90E − 021.80E + 031.0001E + 06Table 4
Dimensionless coefficient ratioA3∗/A1∗–A19∗/A1∗ for relative crack length ratio α = 0.4.
Ratio of high-order termsValueA3∗/A1∗–0.32644A5∗/A1∗–0.22786A7∗/A1∗0.228546A9∗/A1∗–0.39727A11∗/A1∗0.916394A13∗/A1∗–1.86374A15∗/A1∗3.640217A17∗/A1∗–9.39833A19∗/A1∗23.50735It is worth noting that ratioAn∗/A1∗ begins to alternate between positive and negative, starting with A5∗/A1∗.


## 3. Experimental Data

Hoover et al. [26, 27] conducted a comprehensive fracture test on a batch of concrete specimens in 2012, including 128 beams in 4 sizes and similar geometries. Many concrete property parameters were obtained with a small coefficient of variation, and the results were very reliable. All the test data are publicly posted on Professor Bazant’s personal web page, including detailed loading curves (http://www.civil.northwestern.edu/people/bazant/). In this paper, the results of the three-point bending tests on notched beams with an initial crack length ratio of α = 0.3 were selected to calculate the unstable fracture toughness Kun, a total of 28 notched beams, as shown in Table 5.Table 5
Geometric parameters of specimens for three-point bending tests on notched beams [26].
LabelDepthH (mm)WidthB (mm)SpanS (mm)a0/HAa0150038.1333310890.3Aa0250037.6333310880.3Aa0350041.7777810880.3Aa0450037.088891088.50.3Aa0550038.8888910890.3Aa0650039.5444410880.3Ba0121540.244444680.3Ba0221538.38889468.50.3Ba0321539.166674700.3Ba0421539.11111468.50.3Ba0521538.84690.3Ba0621537.911114690.3Ca019341.808201.50.3Ca029341.474201.50.3Ca039341.176200.830.3Ca049340.12201.40.3Ca059342.028202.950.3Ca069341.118202.420.3Ca079341.548201.2950.3Ca089341.4201.840.3Da014040.29286.820.3Da024040.89884.850.3Da034040.77687.010.3Da04404086.030.3Da054040.07487.140.3Da064040.2187.070.3Da074039.88286.60.3Da084040.79486.430.3Wendner et al. [36] collated the raw data of Hoover et al., including the mean response curves of beams, and the results were published in the database at http://www.baunat.boku.ac.at/cd-labor/downloads/versuchsdaten.For convenience, the test data of the database are used directly, and the nominal stress-strain curves are transformed into a load-crack opening curve (P-CMOD), as shown in Figure6. These data are the key to analyzing the calculation of unstable fracture toughness based on the double-K model. Since the coefficient of variation of these data is small, the results are reliable.Figure 6
Average response curves of load versus crack mouth opening-displacement (P-CMOD) curves by three-point bending tests on notched beams withα = 0.3.The authors of this paper are very grateful to Hoover et al. and Wendner et al. for publishing the detailed test results.


## 4. Results and Discussion

### 4.1. Unstable Fracture Toughness

The maximum loadPmax and critical CMODc are directly obtained from the P-CMOD curves, and the critical effective crack length ac and unstable fracture toughness Kun can be calculated according to equations (1)–(6). The results are shown in Table 3. The initial crack length of the concrete specimen will increase slightly due to the influence of the storage and handling process. Xu et al. proposed a method for correcting the initial crack a0, which can be calculated by equation (6). This method was also used to recalculate the initial crack a0 of the test piece, and the results are shown in Table 3. It is noted that only the unstable fracture toughness Kun is calculated without the initial fracture toughness Kini; that is, only one parameter of the double-K model is calculated.In order to obtain the law that the critical effective crack lengthac increases with the height of the specimen, the critical effective crack growth ac is fitted according to the data of Table 3 and is plotted in Figure 7:(16)ac=0.368H+0.0052,R2=0.9986.Figure 7
Critical effective crack length ac versus heightH.


### 4.2. Bazant SEL

For comparison purposes, the unstable fracture toughnessKun values were analyzed using Bazant’s size effect law (SEL) [15]:(17)KUN=KIcSEL1+H0/H,where KIcSEL is the fracture toughness when the specimen height is infinite according to Bazant’s size effect law and H0 is an empirical coefficient.Equation (17) can be rewritten as(18)1KUN2=H0HKIcSEL2+1KIcSEL2.This is a linear function about 1/H. The results of the unstable fracture toughness Kun obtained in Table 4 can be plotted in Figure 8. The equation is as follows:(19)1KUN2=2.2619E−14H+4.1008E−13,R²=9.7542E−01.Figure 8
Linear regression for calculatingKIcSEL and H0.Then, the values ofKIcSEL and H0 can be obtained as follows:(19)KIcSEL=1.56MPa,H0=0.0552mm.The above values are substituted into equation (17), and the Kun values for H in 0.4–1.0 are plotted in Figure 9. KIcSEL is very close to the value KIc = 1.53 MPa calculated by the initial fracture energy [36].Figure 9
Variation in values ofKif/KIc versus height H of specimens (An∗ is the dimensionless coefficients of the terms in Williams’s series expansion).


### 4.3. MMTS SEL

According to the results of Wendner et al. [36], the parameters of concrete used in this test are as follows: Poisson’s ratio v = 0.172, KIc = 1.53 MPa, and tensile strength ft = 5 MPa. Plugging them into equation (12), ∆r is 6.3934E − 03. Substituting ∆r into equation (16), Δr/ac is computed as follows and drawn in Figure 10:(20)Δrac=6.3934E−030.368H+0.0052.Figure 10
Variation in ratioΔr/ac versus height H of specimen.As can be seen from Figure10, as the height of the test piece increases, Δr/ac decreases continuously. Especially when H is less than 200 mm, Δr/ac decreases rapidly. The ratio Δr/ac tends to be stable for H greater than 200 mm.Thus, the values ofKif/KIc can be calculated, substituting equation (20) into equation (15). The results are shown in Table 6 and Figure 9. Columns in Table 6 represent the values of Kif/KIc calculated with respect to the highest degree of high-order terms (except the first two columns). For example, the column labeled A7∗/A1∗ represents that when the value of Kif/KIc is calculated according to equation (15), and the highest order term is A7.Table 6
Values ofKif/KIc calculated by modified maximum tangential stress (MMTS) criterion.
H (m)Δr/acA3∗/A1∗A5∗/A1∗A7∗/A1∗A9∗/A1∗A11∗/A1∗A13∗/A1∗A15∗/A1∗A17∗/A1∗A19∗/A1∗0.040.3300.6760.5520.6100.5670.6070.5750.5990.5760.5970.0930.1670.8360.8050.8120.8090.8110.8100.8100.8100.8100.150.1090.8930.8800.8820.8810.8810.8810.8810.8810.8810.2150.0780.9240.9170.9170.9170.9170.9170.9170.9170.9170.30.0570.9440.9410.9410.9410.9410.9410.9410.9410.9410.350.0490.9520.9490.9490.9490.9490.9490.9490.9490.9490.40.0430.9580.9560.9560.9560.9560.9560.9560.9560.9560.450.0390.9620.9610.9610.9610.9610.9610.9610.9610.9610.50.0350.9660.9650.9650.9650.9650.9650.9650.9650.9650.60.0290.9710.9710.9710.9710.9710.9710.9710.9710.9710.70.0250.9750.9750.9750.9750.9750.9750.9750.9750.9750.80.0220.9780.9780.9780.9780.9780.9780.9780.9780.9780.90.0200.9810.9800.9800.9800.9800.9800.9800.9800.98010.0180.9830.9820.9820.9820.9820.9820.9820.9820.982


### 4.4. Analysis and Comparison of Two Size Effect Laws

(1)
As shown in Figure9, the unstable fracture toughness (Kun) exhibits a significant size effect; that is, the values of Kun increase with the heights (H) of the specimens. But when H > 200 mm, the size effect is no longer obvious and the values of Kun did not increase significantly. This is consistent with the results observed in the experiment by Xu and Zhang [18]. This indicates that larger specimens should be selected when the double-K model is adopted.(2)
The size effect law obtained according to the modification maximum tangential stress (MMTS) criterion is basically consistent with Bazant’s size effect law, and both can reasonably explain the test results of the unstable fracture toughness. However, the values ofKun estimated by the MMTS criterion are slightly larger than the SEL for the H < 1 m. This difference decreases as the height of the specimen increase. When H > 1 m, both are consistent.(3)
According to expression (15), the value of Kif/KIc varies with different highest order terms (An), especially when the height of specimen is less than 200 mm. The most obvious is that when H = 40 mm, the value of Kif/KIc fluctuates significantly. The value of Kif/KIc is the maximum when the highest order term takes only the A3 term. The minimum value is when the highest order item is A5. This difference gradually decreases with the increase in the height of specimen. The Kif/KIc values fluctuate between the maximum and minimum values when more high-order terms are taken. This should be related to the positive and negative alternation of the An∗/A1∗ ratio. The fluctuations gradually decrease with the increase in the height of specimen. The values of Kif/KIc are stable when the H is greater than 200 mm. This is because Δr/ac decreases continuously with the increase in H; the corresponding Δr/acn decreases rapidly so that the influence of high-order terms on Kif/KIc values becomes smaller and smaller. It can be seen that when H > 200 mm, the highest order item only takes the A3 item, and the values of Kif/KIc are almost the same as the higher order item.(4)
ForH > 200 mm, the values of Kif/KIc are greater than 0.9, which can meet the needs of engineers, according to the size effect law obtained by MMTS criterion. Therefore, when three-point bending tests on notched beams are used to determine the unstable fracture toughness, the minimum height of the specimen is suggested to be 200 mm. It is noted that the size effect is due to the effect of Δr/ac, according to expression (15). Since Δr directly depends on (KIc/ft)2, if the ratio of (KIc/ft) of the material is large, Δr will also increase rapidly. Thus, the minimum size required will also increase.


## 4.1. Unstable Fracture Toughness

The maximum loadPmax and critical CMODc are directly obtained from the P-CMOD curves, and the critical effective crack length ac and unstable fracture toughness Kun can be calculated according to equations (1)–(6). The results are shown in Table 3. The initial crack length of the concrete specimen will increase slightly due to the influence of the storage and handling process. Xu et al. proposed a method for correcting the initial crack a0, which can be calculated by equation (6). This method was also used to recalculate the initial crack a0 of the test piece, and the results are shown in Table 3. It is noted that only the unstable fracture toughness Kun is calculated without the initial fracture toughness Kini; that is, only one parameter of the double-K model is calculated.In order to obtain the law that the critical effective crack lengthac increases with the height of the specimen, the critical effective crack growth ac is fitted according to the data of Table 3 and is plotted in Figure 7:(16)ac=0.368H+0.0052,R2=0.9986.Figure 7
Critical effective crack length ac versus heightH.


## 4.2. Bazant SEL

For comparison purposes, the unstable fracture toughnessKun values were analyzed using Bazant’s size effect law (SEL) [15]:(17)KUN=KIcSEL1+H0/H,where KIcSEL is the fracture toughness when the specimen height is infinite according to Bazant’s size effect law and H0 is an empirical coefficient.Equation (17) can be rewritten as(18)1KUN2=H0HKIcSEL2+1KIcSEL2.This is a linear function about 1/H. The results of the unstable fracture toughness Kun obtained in Table 4 can be plotted in Figure 8. The equation is as follows:(19)1KUN2=2.2619E−14H+4.1008E−13,R²=9.7542E−01.Figure 8
Linear regression for calculatingKIcSEL and H0.Then, the values ofKIcSEL and H0 can be obtained as follows:(19)KIcSEL=1.56MPa,H0=0.0552mm.The above values are substituted into equation (17), and the Kun values for H in 0.4–1.0 are plotted in Figure 9. KIcSEL is very close to the value KIc = 1.53 MPa calculated by the initial fracture energy [36].Figure 9
Variation in values ofKif/KIc versus height H of specimens (An∗ is the dimensionless coefficients of the terms in Williams’s series expansion).


## 4.3. MMTS SEL

According to the results of Wendner et al. [36], the parameters of concrete used in this test are as follows: Poisson’s ratio v = 0.172, KIc = 1.53 MPa, and tensile strength ft = 5 MPa. Plugging them into equation (12), ∆r is 6.3934E − 03. Substituting ∆r into equation (16), Δr/ac is computed as follows and drawn in Figure 10:(20)Δrac=6.3934E−030.368H+0.0052.Figure 10
Variation in ratioΔr/ac versus height H of specimen.As can be seen from Figure10, as the height of the test piece increases, Δr/ac decreases continuously. Especially when H is less than 200 mm, Δr/ac decreases rapidly. The ratio Δr/ac tends to be stable for H greater than 200 mm.Thus, the values ofKif/KIc can be calculated, substituting equation (20) into equation (15). The results are shown in Table 6 and Figure 9. Columns in Table 6 represent the values of Kif/KIc calculated with respect to the highest degree of high-order terms (except the first two columns). For example, the column labeled A7∗/A1∗ represents that when the value of Kif/KIc is calculated according to equation (15), and the highest order term is A7.Table 6
Values ofKif/KIc calculated by modified maximum tangential stress (MMTS) criterion.
H (m)Δr/acA3∗/A1∗A5∗/A1∗A7∗/A1∗A9∗/A1∗A11∗/A1∗A13∗/A1∗A15∗/A1∗A17∗/A1∗A19∗/A1∗0.040.3300.6760.5520.6100.5670.6070.5750.5990.5760.5970.0930.1670.8360.8050.8120.8090.8110.8100.8100.8100.8100.150.1090.8930.8800.8820.8810.8810.8810.8810.8810.8810.2150.0780.9240.9170.9170.9170.9170.9170.9170.9170.9170.30.0570.9440.9410.9410.9410.9410.9410.9410.9410.9410.350.0490.9520.9490.9490.9490.9490.9490.9490.9490.9490.40.0430.9580.9560.9560.9560.9560.9560.9560.9560.9560.450.0390.9620.9610.9610.9610.9610.9610.9610.9610.9610.50.0350.9660.9650.9650.9650.9650.9650.9650.9650.9650.60.0290.9710.9710.9710.9710.9710.9710.9710.9710.9710.70.0250.9750.9750.9750.9750.9750.9750.9750.9750.9750.80.0220.9780.9780.9780.9780.9780.9780.9780.9780.9780.90.0200.9810.9800.9800.9800.9800.9800.9800.9800.98010.0180.9830.9820.9820.9820.9820.9820.9820.9820.982


## 4.4. Analysis and Comparison of Two Size Effect Laws

(1)
As shown in Figure9, the unstable fracture toughness (Kun) exhibits a significant size effect; that is, the values of Kun increase with the heights (H) of the specimens. But when H > 200 mm, the size effect is no longer obvious and the values of Kun did not increase significantly. This is consistent with the results observed in the experiment by Xu and Zhang [18]. This indicates that larger specimens should be selected when the double-K model is adopted.(2)
The size effect law obtained according to the modification maximum tangential stress (MMTS) criterion is basically consistent with Bazant’s size effect law, and both can reasonably explain the test results of the unstable fracture toughness. However, the values ofKun estimated by the MMTS criterion are slightly larger than the SEL for the H < 1 m. This difference decreases as the height of the specimen increase. When H > 1 m, both are consistent.(3)
According to expression (15), the value of Kif/KIc varies with different highest order terms (An), especially when the height of specimen is less than 200 mm. The most obvious is that when H = 40 mm, the value of Kif/KIc fluctuates significantly. The value of Kif/KIc is the maximum when the highest order term takes only the A3 term. The minimum value is when the highest order item is A5. This difference gradually decreases with the increase in the height of specimen. The Kif/KIc values fluctuate between the maximum and minimum values when more high-order terms are taken. This should be related to the positive and negative alternation of the An∗/A1∗ ratio. The fluctuations gradually decrease with the increase in the height of specimen. The values of Kif/KIc are stable when the H is greater than 200 mm. This is because Δr/ac decreases continuously with the increase in H; the corresponding Δr/acn decreases rapidly so that the influence of high-order terms on Kif/KIc values becomes smaller and smaller. It can be seen that when H > 200 mm, the highest order item only takes the A3 item, and the values of Kif/KIc are almost the same as the higher order item.(4)
ForH > 200 mm, the values of Kif/KIc are greater than 0.9, which can meet the needs of engineers, according to the size effect law obtained by MMTS criterion. Therefore, when three-point bending tests on notched beams are used to determine the unstable fracture toughness, the minimum height of the specimen is suggested to be 200 mm. It is noted that the size effect is due to the effect of Δr/ac, according to expression (15). Since Δr directly depends on (KIc/ft)2, if the ratio of (KIc/ft) of the material is large, Δr will also increase rapidly. Thus, the minimum size required will also increase.


## 5. Conclusions

In this paper, three-point bending tests of notched beams with span-depth ratioS/H = 2.17 are simulated by the finite element method. The stress intensity factor, CMOD expression, and its inverse function are derived in detail. Then, unstable fracture toughness Kun is calculated according to the double-K model using the experimental data published by Hoover et al. [26, 27] and Wendner et al. [36]. Finally, the size effect of unstable fracture toughness is investigated by the MMTS criterion and compared with Bazant’s size effect law. It is summarized as follows:(1)
The normalized stress intensity factorY(α) and dimensionless shape function V(α) (equations (3) and (4)) are derived for notched beams with S/H = 2.17. It is convenient and efficient to use Abaqus finite element software, which can solve the stress intensity factor of any geometry and loading conditions of test specimens while ensuring accuracy.(2)
As the heightH of the specimen increases, the effective crack length ac increases, but the relative effective crack length ratio ac/H decreases. The ratio Δr/ac also decreases. Δr/ac is the inverse proportional function of H (equation (20)).(3)
The size effect of unstable fracture toughness can be reasonably explained by the MMTS criterion. Three-point bending tests on notched beams are used to determine unstable fracture toughness. The minimum height of the test piece is recommended to be 200 mm. Of course, the minimum height should be related to theKIc/ft value of the material. If KIc/ft is large, a piece with greater height should be used.(4)
The higher order terms of the stress field at the crack tip (Williams’s series expansion) are important factors affecting unstable fracture toughness. In addition to the A3 term, the higher order term (A5, A7, …) also plays a role, especially for small specimens. This is because more high-order terms can more accurately describe the crack tip stress field, which is consistent with Wei et al. [37]. The influence of the higher order term (A5, A7, …) on unstable fracture toughness gradually decreases as the height of the specimen decreases. Considering only the A3 term, unstable fracture toughness will be overestimated for small specimens. However, for large specimens, just A3 can be considered.



---

*Source: 3986367-2020-09-10.xml*