Experimental and Numerical Investigation on the Tensile Fracture of Compacted Clay

: This paper performed flexural test and numerical simulation of clay-beams with different water contents to study the tensile fracture of clay soil and the relevant mechanisms. The crack initiation and propagation process and the accompanied strain localization behaviors were all clearly observed and analyzed. The exponential cohesive zone model was proposed to simulate the crack interface behavior of the cohesive-frictional materials. The experimental results show that the bending capacity of clay-beams decrease with the water content, while those of the crack mouth opening displacement, crack-tip strain and the strain localization range increase. The numerical predictions successfully reproduce the evolving tensile cracks and the strain localization phenomenon of the clay beams with different fracture ductility, which demonstrates the validity of the proposed cohesive zone model in modelling clay fractures.


Introduction
Tensile fracture of soils caused by loading widely exists in many geotechnical applications, such as in earth dams, slopes, landfills, and foundation pits, etc. Many researches showed that instabilities and failures of soil structures were triggered by the tension cracks [Terzaghi (1943); Michalowski (2017) ;Tang, Zhao, Luo et al. (2019)], but the relevant work was very limited due to a dominant compression state of soils being in most situations [Zhang, Li, Yuan et al. (2014)]. These hindered us to further study soils' tensile behavior and to acquire a higher accuracy in numerical modeling soil structures. Recently, many efforts both from experimental and numerical aspects have been conducted to enrich our cognition on soil cracking ; Amarasiri, Costa and Kodikara (2011); Barani, Mosallanejad and Sadrnejad (2015); Xu, Wu, Jin et al. (2018)]. Various test apparatus using direct or indirect measurement has been developed to evaluate the fracture properties of soil. The direct method covers tensile mold direct pull tests [Tamrakar, Mitachi, Toyosawa et al. (2005); Tang, Pei, Wang et al. (2014)], uniaxial tensile test [Ibarra, McKyes and Broughton (2005) ;Zhang, Li, Yuan et al. (2014)] and triaxial tensile test [Parry (1960)]. The tensile strength of the specimen can be easily determined by dividing the ultimate load by the cross-sectional area. However, the difficulty lies in how to fix the specimen at its end so as to avoid a stress concentration, which may lead to the first break of the specimen at the boundary. Due to the disadvantages of the direct method as mentioned above, the indirect method was further introduced, such as Brazilian tensile test [Narain and Rawat (1970)], double punch test [Li, Liu, Jiang et al. (2011)] and flexural beam test [Ajaz and Parry (1975)]. Using these methods, the fracture properties of soil under different influence factors (e.g, dry density, water content, mineral constituent, etc.) have been studied comprehensively. However, these studies are mainly focused on the estimation of the tensile strength, and few attentions have been paid to the cracking-process, especially the evolution of the strain localization behavior. For the numerical simulation of soil-cracking, the method based on the linear elastic fracture mechanics (LEFM) has gained many interests in the past [Morris, Graham and Williams (1992); Hallett, Dexter and Seville (1995)]. Fracture parameters, e.g., fracture toughness and stress-intensity factor, were correlated to the tensile strength or specimen size to demonstrate the validity of LEFM [Hanson, Hardin and Mahboub (1994); Wang, Zhu, Chiu et al. (2007)]. However, as Hallett et al. [Hallett and Newson (2005)] and Prat et al. [Prat, Ledesma, Lakshmikantha et al. (2008)] pointed out, LEFM may be more valid for dry soil with a brittle fracture behavior. Actually, the failure of clay soil generally shows an evident post-peak softening behavior, which contradicts with the LEFM [Tamrakar, Mitachi, Toyosawa et al. (2005); Amarasiri, Costa and Kodikara (2011)]. In such cases, the non-linear deformation range ahead of a crack tip becomes remarkable, as compared to the specimen dimensions, which shouldn't be neglected in the analysis. Therefore, by combining with some advanced FEMs, the cohesive zone model (CZM) [Xu and Needleman (1994); Zeng and Li (2010); Tserpes and Koumpias (2012); Park and Paulino (2012)] and the damage mechanics [Li and Ren (2009) ;Ren, Chen, Li et al. (2011); Ren and Li (2012)] become an alternative solution. The CZM was originally proposed by Barenblatt et al. [Barenblatt (1959); Dugdale (1960)] to describe the plastic-zone beyond a crack tip. The bridging stress in the fracture process zone [Carpinteri and Colombo (1989)] is assumed as a function of the opening displacement. In the past, the CZM was widely used for the cracking-process simulation of metals, concretes and composite materials. Recently, Amarasiri et al. Amarasiri, Costa and Kodikara (2011)] and Barani et al. [Barani, Mosallanejad and Sadrnejad (2015)] made the efforts to introduce the CZM into the mode I fracture analysis of soil. Using the mono-linear or bi-linear CZM, they modeled and analyzed the tensile fracture of the clay-beam with different water contents, additives, and consolidation. However, it's still uncertain whether the fractured index, such as the crack depth or width, is matching with the experimental results. Sanborn et al. [Sanborn and Prevost (2011)] performed a similar work mentioned above, but the tensile properties of soil were not considered in their investigations. Although these researches have shown some preliminary success in the fracture analysis of soil by using the CZMs, more specific works that can give a deeper insight into the fracture behaviors is still necessary. In this work, we performed a further study on the tensile fracture of clay soil and the relevant mechanism has been revealed. The paper is organized as follows: in Section 2, the flexural tests of clay-beams with different water contents were carried out; in Section 3, the cracking-process of soil-beams was observed and further analyzed through the particle image velocimetry (PIV) technique; in Section 4, a modified exponential CZM was proposed for the fracture analysis of soil with frictional behaviors and the detailed implementation procedure was also presented; the experimental and predicted results were given in Section 5, where the comparison has been made to demonstrate the validity of the proposed CZM; finally, the main conclusions are outlined in Section 6.
2 Materials, specimen preparation, and test procedures, etc.

Materials
The typical kaolin soil with a mean particle size (d50) of 3.1×10 -3 mm was used in this paper to prepare soil-beam specimens. It is a mixture of white, delicate and soft clay, and the content of clay is about 38%. Fig. 1 presents the particle-size distribution curve of this used soil. Tab. 1 shows main physical properties associated with specimen preparation.

Specimen preparation
The clay soil was pretreated firstly, including the process of drying, crushing, and sieving. Then it was mixed with a specific weight of pure water and was sealed by a polyethylene film for more than 24 h. When the moisture homogenization has achieved, we put the soil by four equal layers into the split mold. Each layer must be roughened before preliminary compaction. The internal size of the mold is 500 mm (Length)×180 mm (Height)×100 mm (Width). Using a large consolidation instrument, the soil in the mold was compacted to the required height (H=100 mm), thus a rectangular specimen with size of 500 mm×100 mm×100 mm (L×H×W) was obtained (see Fig. 2). After removing the specimens from the mold, a stable condition with 20 °C and 95% humidity was provided to maintain them to relieve possible shrinking. Four specimens with water contents of 32.6%, 34.6%, 35.6%, and 37.0% were tested in this paper to analyze the effects of water contents on the cracking properties of the clay soil. The dry density of all specimens was set as 1.29 g/cm 3 .

Figure 2:
The prepared soil-beam specimen

Test procedures
Note: 1-displacement sensor; 2-soil-beam; 3-measurement area of PIV; 4-Stype force sensor; 5-ball screw elevator; 6-reducer; 7-electric motor; 8-servo controller. As shown in Fig. 3, the soil-beam is placed in the flexural beam apparatus, which is consisting of the reaction bracket, loading system, and measuring system. The reaction bracket contains four steel poles and three steel plates. The roof and floor plates were fixed with the steel poles. The middle loading plate can freely slide along the steel poles.
The loading system is constituted by the electric motor, servo controller, reducer and ball screw elevator. The displacement-controlled loading test was carried out in this work. Using the loading system, the loading plate is forced to vertically move with a specific constant velocity, 0.10 mm/min in this test.
The measuring system consists of five displacement sensors and one S-type force sensor. The displacement sensors fixed on the top surface of specimens are used to measure vertical displacement of the loading point, where the measurement range is 50 mm and the accuracy is 0.001 mm. The force sensor placed above the ball screw elevator is used to measure the mid-span load, whose measurement range is 1 kN and accuracy is 0.3% kN. The data acquisition frequency of all sensors is 0.5 Hz. Before testing, all equipment must be examined and calibrated. The soil parameters (e.g., specimen size, weight, and water content) were measured and recorded. After that, we started the apparatus until soil-beams entirely cracked.

Strain measurement
The strain in potential cracking area of soil-beams was measured by the particle image velocimetry (PIV) technique [White, Take and Bolton (2003)  Based on the specimen photos taken before and after the deformation, PIV can identify and trace the texture features of the specimen and present the displacement distribution. Taking the test patch A with initial location (u1, v1) for example: the texture features in moment t1 (see image 1) are identified and stored, then the texture matching and correlation operation is performed on the searching region shown in Fig. 4 (image 2). According to the results of the degree of match, the current locations of test patch A in moment t2, i.e., u2 and v2, is updated. Using a similar process, the locations of the full analysis area at different loading moment can be determined.
Note that the size of test patches directly decides the measurement accuracy of displacement field and strain. Here test patches with the accuracy of 0.052 mm per pixel were used to capture the displacement/strain fields of soil-beams. By using the PIV technique, the location of soil-beams in the analysis area at different loading moment can be obtained, thus the strain of specimens in x-direction can be calculated as illustrated below in Fig. 5.
The first image of test patch The n th image of test patch With the increase of mid-span deflection, a crack irised out by the dotted line appears at the bottom of the specimen, and the initiation point is closed to the mid-point. Although the crack path is slightly tortuous, the trajectory is nearly a straight-line. These results are in agreed with the elastic beam theory. Moreover, it indicates that soil-beams are uniform enough to guarantee the reliability of the experimental results. Similar results can also be acquired from other specimens.
(a) stage 1 (b) stage 2 (c) stage 3 Figure 6: The cracking-process of the soil-beam with a water content of 32.6% Fig. 7 presents the crack-mouth-opening-displacement (CMOD) vs. deflection curves of the soil-beam with different water contents. At the initial stage, the CMOD is at a lowlevel close to zero and no cracks can be observed (as shown in Fig. 6 stage 1). With the mid-span deflection further increases, the crack initiates along the bottom center of soilbeams (i.e., Fig. 6 stages 2 and 3). The growth rate of the CMOD in all specimens is nearly linear with the loading-point deflection and does not have a correlation with the water content. Additionally, the CMOD of soil-beam with high water content is larger than that of low water content specimen during the crack propagation process. As the water content increases from 32.6% to 37%, the CMOD at the fracture failure point increased by 3 times. Results indicate that the increase of water content will aggravate the tensile fracture of soil-beams. where the gravity of all specimens is eliminated in advance from the reactive force. As the water content increases, the initial bending stiffness and peak load of soil-beam decrease. But the mid-span deflection at both of the crack initiation point and fracture failure point increases due to the non-linear deformation of the high water content specimen. It indicates that the soil-beam with low water content has a higher bending capacity, while that of high water content shows better ability in resisting the differential settlement. Furthermore, when crack initiates, the soil-beam with low water content, i.e., w=32.6%, exhibits a sharp decrease in the mid-span load, which shows a brittle fracture behaviors.

Experimental curves
For soil-beam with the high water content, i.e., w=37%, the load-deflection curve shows a strong non-linear behavior and the fracture type is more close to the quasi-brittle fracture. In general, for the clay-beams, the increase of water content leads to the change of fracture types from brittle to quasi-brittle.

The strain localization behavior
To better analyze the strain localization behavior of soil-beam, both the strain vs. deflection curves and load vs. deflection curves are given in Fig. 9, where the strain is derived from the crack onset-point. In addition, the strain contours in the PIV observation area are also added to Fig. 9 to illustrate the evolving cracks in soil-beam. Based on the strain development states, we divided the strain curve of the four specimens into two sections: the linear section and exponential section. The dotted lines marked in Fig. 9 denotes the dividing line. As seen, the strain-deflection curve before the dividing line is nearly linear and the maximum strain of four specimens are less than 0.5%, which means that the soil-beams undergo a continuous deformation without any cracks. However, when the deflection reaches the dividing line, the trend of the strain curves suffer from a sudden transform. The strain localization phenomenon can be observed at the bottom of these specimens, as shown in the strain contours marked as tA. From the point of tA on, the soil-beams exhibits the so-called weak discontinuity (strain discontinuity). What's more, the onset-points of the strain localization behavior are all located in the ascendant stage of the load-deflection curves, which proves the viewpoint that strain localization in the plane problems can initiate at the hardening stage of the stress-strain curve [Rudnicki and Rice (1975); Lade and Wang (2001)]. With a continued loading, the micro-cracks conceiving in the strain localization region develop into one or more macro-cracks (see strain contours in tB). The soil-beams exhibit a strong discontinuity (displacement discontinuity) until the main crack propagates to the top surface. The onset and evolution of the strain localization is the key factor that leads to the tensile fracture of soils.
(a) w=32.6% (c) w=35.6% (d) w=37.0% Figure 9: the strain-deflection curves of soil-beams at the crack onset-point and the corresponding strain contours From the strain contours at time tA, tB, and tC, we can see that the soil-beams with different water contents suffer from a tensile deformation at the bottom and a compression deformation at the top. By comparing the size of strain localization region, we can draw that the higher the water content, the larger the strain localization regions and the higher strain appear at the crack tip. This is mainly because the specimen with higher water content usually shows a lower tensile strength, which is easier to meet the strain localization condition.

CZM and numerical implementation
As shown in the above-mentioned flexural tests, the soil-beam with high water content shows an evident post-peak softening behavior, which makes the fracture type of soils transform from brittle to quasi-brittle. In these cases, the LEFM theories are not applicable anymore. Hence, an alternative method, i.e., the CZM approach, is introduced to the fracture analysis of such the geotechnique problems. In this work, the exponential CZM, which was originally proposed by Xu et al. [Xu and Needleman (1994)], was introduced and modified to consider the compression-shear coupling effects of the frictional materials and to simulate the tensile fracture of the clay beams.

The cohesive law
For the plane problems, the mixed-mode cohesive law shall involve both the normal and tangential behaviors of the fracture surface. Thus, the resultant displacement δ and the resultant traction t are predefined [Ortiz and Pandolfi (1999) where δn, δs is the normal and tangential opening displacement across the fracture surface; tn and ts are the respective tractions in the normal and tangential directions; β is a weight factor evaluating the traction components, 1 β = in this paper.
If the interface is in tension ( n t ＞0), the free energy potential with an exponential form [Xu and Needleman (1994)] is adopted, where σp denotes the tensile strength of materials; δp is the critical displacement corresponding to the peak strength; ( ) exp 1 e = . The introduction of φ brings out much computationally convenient in describing the cracking-process near the crack tip. Note that the fracture energy should be equal to φ in physics.
When the interface is on loading, the cohesive law is obtained by taking the first derivative of the free energy potential in Eq. (2), When the interface is on unloading and reloading, a cohesive law, which is similar to the form reported in Park et al. [Park and Paulino (2012)], was adopted: where n k is a penalty stiffness used for preventing material self-penetration, which is usually set as the gradient of the Eq.
where c and ϕ are the cohesion and internal friction angle of the interface. For simplicity, they are assumed to be equal to the shear strength index of soils. Based on Eq. (3.2) and Eq. (6), the cohesive law in the tangential direction is written as, The cohesive law upon unloading or reloading is similar to the Eq. (4). Considering the tangent stiffness is an essential part for the assembly of the stiffness matrix, we give it by taking the first derivative of the cohesive formulations, n n n s nn ns c sn ss s s n s Here the elements in Eq. (8) is not provided, but they can be easily derived from Eq. (3)-Eq. (7).

Numerical implementation
The proposed CZM was implemented as a four-node user element (UEL) of ABAQUS.
Since the Newton-Raphson method was chosen for the iteration solution of equations system, the cohesive force vector and the tangent stiffness matrix must be formulated to define the above UEL. They were derived from the weak form of the governing equations with a cohesive crack, Here, N is a standard matrix of shape functions of a one-dimensional element; is the relative displacement vector of nodes. For instance, x 1,4 ∆ denotes the relative displacement between node 1 and node 4 along the x-direction. , sin cos where cos sin sin cos is the rotational matrix; θ is the included angle between the local coordinate and the global coordinate. Substituting Eqs. (10) and (11) into the second item of Eq. (9), the cohesive force vector can be expressed as The tangent stiffness matrix was given by, where d is the displacement vector of nodes; c D is the tangent stiffness matrix of materials defined by Eq. (3)-Eq. (7).

Verification of CZM
In this section, the uniaxial tensile test ] and the direct shear test [Bransby, Davies, Nahas et al. (2008)] were used to achieve the following purposes: (i) verify the validity of the proposed CZM as a UEL of ABAQUS; (ii) analyze the sensitivity of the cohesive parameters to the numerical responses and (iii) demonstrate the necessity of correlating the shear strength to the normal tractions.

The uniaxial tension model
As described in Fig. 11, a homogeneous and isotropic panel was stretched on its two ends with relative displacement of δx. P is the reactive force at the loading point. The plane stress condition was assumed. The stress σx of each element can be determined by the ratio of P and the cross-section area A (σx=P/A). As δx increases, σx in all elements simultaneously approaches to σp, thus the panel initiates fracture in the middle of this model. The cohesive elements defined in Section 4 were pre-inserted in the middle joint. The cohesive parameters defining the cohesive element directly refer to those taken in ]. The tensile strength σp was 20.6 kPa and the critical displacement δc was 0.01 mm. The fracture energy Gc was 0.56 N/m. For the bulk material, a linear elastic deformation with Young's modulus E=200 GPa and Poisson's ratio v=0.2 was assumed. Since the deformation mainly concentrates on the pre-inserted cohesive elements, the stress-displacement curve (σx-δx) obtained from Eq.
(3.1) is regarded as the analytical solution. Fig. 12 shows the comparison of the stressdisplacement (σx-δx) curves obtained analytically and numerically. As seen, the predicted curve is nearly the same as the one obtained by the analytical method, which indicates that the inputted cohesive law can be reproduced at the structural scale. In other words, the proposed CZM in Section 4 has been successfully implemented in ABAQUS. The effects of the cohesive parameters, i.e., fracture energy and tensile strength, on the numerical responses were analyzed in this section. First, we set three values of the fracture energy, 0.8Gc, Gc and 1.2Gc, while the tensile strength is constant (σp). Fig. 13(a) presents the load-displacement (P-δx) curves of the panel with different fracture energies. The integral area under the P-δx curve increases with the fracture energy because it reflects the global released energy in creating a new crack. The peak loads of the P-δx curves are the same because the tensile strength of all models is a constant. Changes of the fracture energy have little effects on the shape of the P-δx curves. The second is the tensile strength which was analyzed by setting three values: 0.8σp, σp, and 1.2σp, and the fracture energy was equal to Gc. As shown in Fig. 13(b), the shape of the P-δx curves displays significant differences. As the tensile strength decreases, the P-δx curve becomes more 'short' and 'fat' and the initial gradient is smaller. The integral area of the P-δx curves is nearly identical. In summary, the tensile strength of materials plays a dominant role in characterizing the geometry of the P-δx curve.

The direct shear test
The direct shear test, which satisfies the plane strain condition, was modeled in this section to verify the necessity of the modification for the proposed CZM. As illustrated in Fig. 14, the bottom edge, the left edge and right edge are constrained while the top edge is subjected to a constant compressive stress, P. Under a displacement-controlled loading (i.e., u1), the top box moves horizontally along the shear sliding surface that is bonded by the cohesive elements. Following the British Standards (BS1377-7), the dimensions of this model are 60 mm×30.4 mm (W×H).

Figure 14: Schematic diagram of the direct shear test
Referring to the experimental data reported by Bransby et al. [Bransby, Davies, Nahas et al. (2008)], the material in the top and bottom box is assumed to be linear elastic. Young's modulus E is proportional to compressive stress P (=20, 53, 100, 200 kPa), i.e., E=0.6P+1.2 (MPa) and Poisson's ratio v is 0.3. Since the crack surface suffers from a compressive deformation, the cohesive law shown in Eqs. (5) and (7) are used to describe the fracture behavior of the shear sliding surface. Through the back analysis (presented in subsection 5.2.2), the cohesive fracture parameters can be easily obtained as follows: kn = 1.13×10 9 kN/m 3 , δp=1.2 mm, c=σp=0 kPa (Fontainebleau sand was used in tests), peak internal friction angle φp=38.66 o and residual internal friction angle φr=32.62 o . Fig. 15 compares the shear stress (τxy) vs. tangential displacement (u1) curves under different compressive stress P, where τxy is averaged from all the cohesive elements in the sliding surface. As shown, the trend of the predicted curves in the liner section or the softening and residual section is nearly identical to the experimental curves. Overall, by correlating the shear strength to the normal stress, the interfacial mechanical behavior of the direct shear tests can be accurately reproduced.
However, the τxy-u1 curves obtained by the unmodified CZM show essential differences as compared with the experimental curves. As shown, the predicted curves under different normal stress (P) are nearly straight line. The shear stress is always at a lowlevel close to zero and does not change with the tangential displacement. It cannot reveal the compressive-shear relationship of the frictional materials. Therefore, the unmodified CZM cannot be used to simulate the fracture failure of geotechnique problems with a compressive-shear cracking behavior. This paper, mainly focuses on the tensile fracture of compacted clay. The complex problems with tensile-compressive-shear coupled, such as the instability of slopes, will be analyzed in the subsequent publications by inserting the modified CZM into some advanced FEMs.

The soil-beam model
In this section, the flexural test results in Section 3 were employed to validate the effectiveness of the proposed CZM in simulating the tensile fracture of soils. The physical properties of soil-beams at time of testing were listed in Tab. 2. Referring to the flexural tests, the soil-beam specimen was reasonably simplified (see Fig. 16(a)). A pair of simple supports was added to the lower surface, while a vertical displacement boundary δy was imposed at the mid-span. We assume the soil-beam meets the plane strain condition. To facilitate meshing, the model was divided into two parts (I and II), which was bonded by a cohesive crack in the mid-span (see Fig. 16(b)). Part I/II consists of 2500 structured Q4 elements. The cohesive crack is constituted by 50 cohesive elements with an average size of 3.5 mm.   [Amarasiri, Costa and Kodikara (2011)].
Firstly, an iterative fitting aiming for the initial slope of the load-deflection curves was carried out. When the predicted curve is matched with the experimental curve at the initial section, the desired Young's modulus is acquired (see Tab. 3). Fig. 17(a) describes the relationship between the Young's modulus and water content. With the increases of water content, the elastic modulus of the soil-beam decreases. A remarkable negative correlation between the Young's modulus and the water content can be established. A three-step procedure was used to estimate the cohesive parameters.
Step 1: The global release energy (GF) for creating a unit length crack was determined by dividing the integral area (Sd) of the load-deflection curve by the interface area (H×W) of the soil-beam, i.e., GF=Sd/(H×W). The results of the GF for soil-beams with different water contents were listed in Tab. 2.
Step 2: The fracture energy (Gc) for the proposed CZM was set as an approximation of the GF (although Gc may be smaller than GF because it does not contain the plastic work and the test errors), thus we adjusted the tensile strength σp to make the predicted peakload of the load-deflection curves match the test results.
Step 3: Keeping the fitting values of the tensile strength in Step 2 as a constant, then we adjusted Gc appropriately so that the post-peak softening section of the predicted loaddeflection curve is consistent with the one measured experimentally. By an iterative process of Step 1 to Step 3, the values of cohesive parameters for soilbeams with different water contents can be obtained (see Tab. 3). Fig. 17(b) shows the relationship between tensile strength and water content. The tensile strength of soilbeams decreases as the water content increases. There is an approximate linearity relation between them.  Fig. 18 compares the load-deflection curves obtained numerically and experimentally. By using the fitting parameters, the predicted curve in the hardening section is nearly identical to the experimental curve. In the tail-section, the predicted curve displays a slight difference with the experimental one due to the gravity of the soil-beam, whereas the overall trend of the post-peak softening section is very similar. Meanwhile, the proposed CZM shows good performance in reproducing the transformation of the fracture behavior from the brittle fracture to the quasi-brittle fracture. The obtained loaddeflection curves agree fairly well with the experimental results.
(a) w=32.6% (b) w=34.6% (c) w=35.6% (d) w=37.0% Figure 18: Comparison of the load-deflection curves obtained experimentally and numerically Fig. 19 shows the CMOD-deflection curves obtained numerically and experimentally. As shown, the soil-beam undergoes a linear-elastic or nonlinear-harden deformation in the initial stage. The maximum tension stress is still less than the tensile strength, so no cracks initiate in the soil-beams. Hence, the predicted CMOD is closed to zero, which is consistent with the test results. When the stress exceeds the tensile strength, a quasi-linear relationship between the CMOD and mid-span deflection can be established. The slope of the predicted CMOD curves in the cracking stage is also similar to the test results. Figure 19: Comparison of the CMOD-deflection curves obtained experimentally and numerically Fig. 20 provides the horizontal stress contours of the soil-beam, where the water content is 34.6%. As seen, the horizontal tension stress before cracking (δy=0.5 mm) mainly concentrates on the mid-point bottom, and the tension zone is very small. With the deflection increases (δy=1.0 mm), the tensile zone moves outward and the tensile stress is gradually increased. A crack initiates in the middle of the soil-beam due to the tensile stress attains the soil's tensile strength. Meanwhile, the top surface of this model shows an evident compression deformation. When the deflection further increases from 1.0 mm to 2.1 mm, the stress concentration phenomenon appears at the crack tip. The crack length and CMOD are increased. In summary, by inserting the proposed CZM into the FE method, the crack initiation and propagation process and the stress evolution in the soilbeam can be reproduced well. Fig. 21 compares the stress contours of soil-beams with different water contents, where the crack length is the same. As seen, the water content has essential effects on the stress contour of the soil-beam. The lower the water content, the smaller of the stress concentration zone and the higher of the crack-tip stress display in the soil-beam. With the increase of water content, the stress concentration region near the crack tip shall extend to the bottom of soil-beam. Meanwhile, the soil-beam with high water content shows larger stress than that of the low water content. It indicates that soil-beam with high water content has a lower flexural bearing capacity. The points drawn from the stress contours are consistent with the test results as shown in Fig. 9.

Conclusions
The flexural tests and numerical simulation of clay-beams with different water contents were conducted in this paper to study the tensile fracture of clay soil and relevant mechanisms. The crack initiation and propagation process, the experimental response curves and the companied strain localization phenomenon were observed and analyzed in detail. Based on the modified exponential CZM as proposed in this paper, the fracture characteristics (e.g., crack propagation and strain localizations) of clay-beams with different water contents were simulated accurately. The main conclusions can be drawn as follows: (1) The CMOD of soil-beams after cracking is proportional to the mid-span deflection, and its growth rate has no correlation with the water content. With the increase of water content from 32.6% to 37%, the corresponding CMOD at failure point increased by 3 times.
(2) The soil-beam with low water content has a higher bending capacity, while that with high water content shows better ability in resisting the differential settlement. The increase of water content leads to the change of fracture type of the clay-beams from brittle to quasi-brittle. (3) The occurrence of the strain localization phenomenon is the manifestation of material discontinuity in soil-beams. The fracture responses of soil-beams after cracking are mainly determined by the mechanical behavior at the cracking surface. Both The strain localization zone and the strain near the crack tip increase with the water content.
(4) The proposed exponential cohesive crack model shows high accuracy in simulating the tensile fracture of clay-beams. The crack propagation process and the accompanied strain localization behaviors of the clay-beams with different water content were all reproduced well by the proposed numerical strategy.