Numerical Analysis Method of Shear Properties of Infilled Joints under Constant Normal Stiffness Condition

School of Mines, Key Laboratory of Deep Coal Resource Mining, Ministry of Education of China, China University of Mining and Technology, Xuzhou 221116, China State Key Laboratory of Mining Disaster Prevention and Control Co-founded by Shandong Province and the Ministry of Science and Technology, Shandong University of Science and Technology, Qingdao 266590, China Work Safety Key Lab on Prevention and Control of Gas and Roof Disasters for Southern Goal Mines, Hunan University of Science and Technology, Xiangtan 411201, China


Introduction
With the combined action of tectonics, groundwater movement, and weathering, some joints are filled with soft geologic materials, such as clay and sand.ese soft fillings significantly weaken the shear strength of jointed rock and make the rock present more obvious heterogeneity and anisotropy.Infilled joints generally become the weakest position under shear which will easily trigger the whole instability of engineering rock [1][2][3][4].erefore, the study of shear properties of infilled joints is of significance in both stability evaluation and reinforcement of engineering rock.Many direct shear tests have been conducted to study the influence of different factors on the shear failure mode of infilled joints.
ese influenced factors mainly included material type, surface morphology, initial normal stress, infilled ratio, fillings saturation degree, and overconsolidation degree [1,[5][6][7].Previous test results showed that soft fillings, such as Green mudstone, serpentinite, and clay, will significantly reduce the shear strength of infilled joints [8], while the high-strength fillings, such as mortar and concrete, enhance the shear strength [6].It was found that different morphology of infilled joints may bring about different shear failure mechanism.Furthermore, the large joint undulating degree will prevent the roll of filling grain.It is easy to understand that the shear failure of infilled joints should overcome the sliding friction of fillings.However, when the joint surface is flat providing little resistance to the movement of filling grain, the shear failure only needs to overcome the rolling friction [9,10].Laboratory results from both Ladanyi and Archambault's [11] and Indraratna's [12] direct shear tests showed that the increase of joint undulation was able to enhance the shear strength of infilled joints.On the contrary, normal stress controls the shear failure mode of infilled joints in a certain degree.With the small normal stress, the sliding of infilled joints only occurs in soft fillings and along the surface of protruding block in the sliding face.However, the large normal stress-induced sliding will cut the protruding block which means the yield of both soft llings and rock [11].Meanwhile, the in lled ratio, that is, the ratio of lling thickness, t, and the asperity height, a, is a signi cant factor controlling the shear properties of in lled joints.A large number of direct shear test results [1,3,[12][13][14][15][16] showed that the in lled ratio, t/a, has a critical value.If the in lled ratio is less than the critical value, the shear failure strength of in lled joints is typically higher than the shear strength of llings but lower than that of asperity.If the in lled ratio is greater than the critical value, the shear behaviors are totally controlled by llings that means the shear strength of the in lled joints will be equal to that of the llings.Xu and Ren [17] also obtained similar conclusions from the numerical direct shear test.Indraratna [12] pointed out that the critical value of the in lled ratio under CNS condition was less than that under constant normal load/stress (CNL) condition while small initial normal stress corresponded to a small critical value of the in lled ratio.
Most of previous direct shear tests were conducted under CNL condition [12,[18][19][20][21].However, for underground rock mass, the normal stress applied on rock joints will increase due to both dilatancy e ect and restraining e ect of surrounding rock.
us, the laboratory results under CNL condition cannot accurately re ect the shear failure mode of joints in underground rock mass.Furthermore, current research achievements are mainly based on laboratory direct shear tests.As there are many factors in uencing the shear properties of in lled joints as mentioned above, studying the e ect of all these factors require a large number of direct shear tests which is inconvenient.erefore, numerical simulation, as a repeatable and highly e cient research method, is of great advantages.e paper will provide a numerical simulation method to implement the direct shear test under CNS condition.e capability of the numerical simulation method for evaluating e ects of di erent factors including undulating angles, in lled ratio, and normal stress levels on shear properties of in lled joints will be studied and discussed.

Modeling of In lled Joint.
is paper developed a concept model for in lled joints (shown in Figure 1 2 Advances in Civil Engineering the joint surface is idealized as sawtooth shape.e model consists of two rock blocks lying top and bottom sides, respectively, with soft llings in the middle.An interface is added between the top rock and the lling to simulate the shear slip phenomenon. is study used a wedge-shaped element to build a sawtooth-shaped joint in the middle part of the model while the block element was applied to build the top and bottom parts of the model for improving the calculation e ciency.In the numerical model developed by FLAC3D, as shown in Figure 1

Mechanical Parameters.
e Mohr-Coulomb failure criterion is adopted as the constitutive model for both rock and lling.e bulk modulus and shear modulus used in FLAC3D can be calculated from the elastic modulus and Poisson's ratio based on the generalized Hook's Law.
e normal sti ness and shear sti ness of interface refer to the recommendation value in the FLAC3D manual [22], namely, ten times the equivalent sti ness of the sti est neighboring zone, as shown in the following equation: where k n is the normal sti ness, k s is the shear sti ness, K is the bulk modulus, G is the shear modulus, and Δz min is the smallest width of an adjoining zone in the normal direction.en, the normal displacement of the top block, Δδ n , generated in this cycle is monitored.If the constant normal sti ness, k n , is set to be 4.0 MPa/m, the increment of normal stress, Δσ n , should be Δδ n × k n in this cycle.e normal stress applied on the top boundary in the next shearing cycle will be adjusted to be σ ′ n σ n + Δσ n .e above shearing cycle should be executed 5500 times so that the shearing displacement can reach 0.22 m (larger than the half of bottom length of sawtooth).In addition, the CNL condition can be achieved by setting the constant normal sti ness as 0.

Boundary
In order to validate the reliability of this numerical model, we rotate this model to imitate Indraratna's shear strength experiment of in lled rock joints [14] and compare it with the experimental results (Figure 2).Obviously, comparison shows a good agreement between the numerical model and experimental results.
erefore, this model is considered capable of being further used to analyze the shear mechanical properties of in lled joints under CNS condition.

Variables Monitoring.
e main purpose of the direct shear test is to obtain the changing tendency of three variables (shear stress, normal stress, and normal displacement) along with shearing displacement.In shearing process, the three variables were calculated by the FISH code and extracted by the HISTORY command.e shear stress is the ratio of the maximum x-axis unbalance force of the bottom block to the shearing area of the in lled joint.

Numerical Results and Discussion
In this section, the numerical model of the direct shear test is used to study the in uence of di erent factors, including di erent loading conditions (CNS and CNL), joint undulation, and in lled ratio, on shear properties of rock joints.
en, the numerical results combined with the previous studies [12,[15][16][17]23] are discussed.e mechanical parameters of rock and llings are listed in Table 1.

Shear Properties under CNS and CNL Conditions.
In this test, the in lled ratio, t/a, is 0. e undulating angle is 26.6 °.
e initial normal stresses were set to be 0.1 MPa, 0.3 MPa, 0.5 MPa, 0.7 MPa, 0.9 MPa, 1.1 MPa, and 1.3 MPa.Direct shear tests are conducted under CNS and CNL conditions, respectively, as shown in Figure 3. e shear stress-shear displacement curve, normal displacement-shear displacement curve, and normal stress-shear displacement curve obtained from numerical simulation are shown in Figure 4.
Figure 4(a) illustrates the shear stress evolution of a clean (without llings) joint under CNS and CNL conditions.It is found that if the direct shear test is conducted under CNS  condition or CNL condition, there are two in ection points in the stress-displacement curve.e whole curve can be divided into three parts by these two points which represent elastic part, yield part, and failure part, respectively.In the elastic part, shear stress increases linearly and grows to the rst in ection point at a very small shear displacement.It should be noted that the rst in ection point under di erent loading conditions (CNS and CNL) shows the same shear stress value.In the yield part lying between the two in ection points, the shear stress stays stable under CNL condition.However, it continues increasing linearly with displacement under CNS condition.With a certain shear displacement, the shear stress reaches the second in ection point whose value is regarded as the shear strength of joints.en, the curve comes into the failure part in which the shear stress reduces gradually under both types of conditions.It can be seen that, at any level of initial normal stress, the shear strength under CNS condition is greater than that under CNL.However, the corresponding shear displacement under CNS condition is smaller.In addition, the shear strength of joint shows an upward trend with the increase of initial normal stress, while the corresponding shear displacement reduces gradually under both loading conditions.
Figure 4(b) displays the dilation curves of joint under CNS and CNL conditions.It can be seen that the dilation curve shows a linearly increasing trend at the initial shearing stress stage.It is worth mentioning that the dilatancy angle is equal to the angle of sawtooth on the bottom.is means that there is no deformation failure occurring in the sawtooth, and the top block is sliding along the contour line of sawtooth from the root to the top of sawtooth.With the shear proceeding, the contact area between two sawteeth reduces gradually, causing a decrease in shear strength of sawtooth.Once the strength drops to a certain value, the top block will cut o the sawtooth rather than sliding along the contour line which makes the dilation curve appear an in ection point.en the normal displacement rate gradually decreases, approaching to a certain value eventually.At the same initial normal stress, the cuto position on the sawtooth is more close to the tooth root for CNS condition which means the dilatancy e ect is restrained due to the increase of normal stress.
Figure 4(c) shows the normal stress evolution under two types of conditions.It can be seen that the normal stress under CNL always keeps stable during the shear process.In contrast, the normal stress under CNS experiences a linear increasing trend at the beginning and then decelerates slowly to approach a certain value.
In conclusion, the dilatancy e ect of the joint will be restrained under CNS condition due to the increase of normal stress, causing the squeezing action between two sawteeth to grow gradually.erefore, the sheared position of sawtooth is more close to the tooth root under CNS condition which is good for the sawtooth to exert its antishear ability.Indraratna [12] obtained the similar conclusions through conducting laboratory direct shear tests under CNS and CNL conditions.e di erence is that shear stress curves obtained by Indraratna present obvious strain softening characteristics which were not observed in the numerical results.
is is because that the Mohr-Coulomb criterion applied in the numerical model is not able to re ect the postpeak mechanical characteristics.
In addition, previous studies [24][25][26] showed that in lled joints demonstrate three di erent failure modes based on discrete element numerical simulation which are failure of the interface between rock and llings, sliding along the interface, and the cuto failure of the sawtooth as shown in Figure 5.It can be found that the shear failure mode of in lled joints is in uenced by mechanical properties of rock, llings, and interface.Meanwhile, based on the investigation of the failure mode under di erent conditions (e.g., normal stress, undulating angle, and in lled ratio), this   Advances in Civil Engineering study reveals that shear failure of in lled joints can be divided into three phases which are elastic phase, yield phase, and failure phase while each phase corresponds to a speci c failure mode as shown in the literatures.
erefore, the numerical model in this study shows a good agreement with previous studies [24][25][26].

In uence of Undulating Angle on Shear Properties of
In lled Joint

Shear Properties within Di erent Undulating Angles.
In this section, the in lled ratio of the joint is 0.5.Undulating angles are 16.7 °and 26.6 °.Initial normal stresses are 0.1 MPa, 0.3 MPa, 0.5 MPa, 0.7 MPa, 0.9 MPa, 1.1 MPa, and 1.3 MPa.Direct shear tests are conducted under CNS condition.e obtained shear stress-shear displacement curve, normal displacement-shear displacement curve, and normal stressshear displacement curve are shown in Figure 6.
Figure 6(a) shows shear stress curves of in lled joints within di erent undulating angles.It can be seen if the normal stress is small (e.g., 0.1 MPa, 0.3 MPa, 0.5 MPa, 0.7 MPa, and 0.9 MPa), the two in ection points and three parts illustrated in the previous section can be observed in the shear stress-shear displacement curve.At the same initial normal stress, the shear strength within an undulating angle of 16.7 °is obviously smaller than that of 26.6 °.However, the corresponding displacement shows a contrary relationship which is the displacement with the undulating angle of 16.7 Advances in Civil Engineering failure mode will change from sliding failure to shear failure of sawtooth due to the increase in the undulating angle at a relatively high normal stress condition.
e reason causing the above results seem to be that the small undulating angle of sawtooth has a relatively small resistance to the shear motion of the in lled joint, so the joint is more likely to slide along the contour line.In contrast, the large undulating angle contributes to a large resistance to the shear motion resulting the crush of sawtooth is more likely to occur in the shear process with the large undulating angle.

Shear Strength Curves within Di erent Undulating
Angles. Figure 7 shows the shear strength curve of in lled joints with di erent undulating angles.
It can be seen that the shear strength curve of 16.7 °can be tted in a linear relationship perfectly while the curve of 26.6 °also can be tted in a linear relationship at a low normal stress level.However, once the normal stress increases to a certain value, the shear strength curve with the large undulating angle presents nonlinear distribution.
is is caused by the change of the failure mode that is changing from sliding to cuto failure.In addition, the slope of the strength curve with 26.6 °in low normal stress level is obviously larger than that with 16.7 °which means that the increase of the undulating angle can e ectively raise the friction angle of in lled joints.Furthermore, an additional investigation of in uence of the interface (between rock and llings) friction angle on shear strength of in lled joints was also conducted based on the numerical model.It can be seen in Figure 7  1.8 Fitting curve within Advances in Civil Engineering of in lled joints under this condition (undulating angle: 16.7 °and interface friction angle: 20 °) also increased linearly which is the same as that when the friction angle is 25 °and undulating angle is 16.7 °.
Patton [27] proposed a mathematical model, as shown in the following equation, to predict the shear strength of the sawtooth joint: where τ is the shear strength, σ is the normal stress, φ is the friction angle of interface, and i is the undulating angle.As shown in (2), the equivalent friction angle of in lled joints can be seen as the sum of the interface friction angle and the undulating angle.Put the corresponding friction angle (25 °and 20 °) and undulating angle (16.7 °and 26.6 °) into (2), respectively, and draw three kinds of curves in the coordinate system.It can be seen from Figure 7 that the tted curves are in accordance with the Patton shear strength curve, which indicates that the numerical results have a very good agreement with the shear strength model proposed by Patton [27].

In uence of In lled Ratio on Shear Properties of
In lled Joint

Shear Mechanical Curves within Di erent In lled
Ratios.In the contrast tests of this section, the undulating angles of the model is 26.6 °, and the in lled ratio of joint is 0, 0.25, 0.5, 0.75, 1.0, 1.25, 1.5, and 2.0, respectively.Direct shear tests are conducted under CNS condition, and the initial normal stresses are 0.3 MPa, 0.5 MPa, and 0.7 MPa. e obtained shear stress-shear displacement curve, normal displacement-shear displacement curve, and normal stress-shear displacement curve are shown as Figure 8. Figure 8(a) shows the shear stress-shear displacement curves with the initial normal stress of 0.3 MPa, 0.5 MPa, and 0.7 MPa.It can be seen that if the in lled ratio remains stable at a certain value, the shear strength rises in a certain degree with the increase of initial normal stress.On the contrary, if the initial normal stress remains stable at a certain value, the failure mode of in lled joints will change with the increase of the in lled ratio.
It can be seen from Figures 8(a) and 8(b) that when the normal stress is 0.3 MPa, the shear stress-shear displacement for di erent in lled ratios can be divided into three parts as mentioned previously.e shear displacement in the yield part decreases with the increase of the in lled ratio.e reason of this result is that the proportion of lling in a sawtooth height will rise due to the increase of the in lled ratio.Correspondingly, the failure position will move to the sawtooth root gradually which will reduce the shear displacement in the yield part.
As shown in Figure 8(a), as for the initial normal stress is 0.5 MPa, the shear stress curve shows the three stages as discussed above for the in lled ratio less than 1.0.On the contrary, the shear stress curve only has two parts for the in lled ratio larger than 1.0.It can be seen that once the in lled ratio rises to or exceeds 1.0, the corresponding shear strength reduces signi cantly.If the in lled ratio is less than 1.0, the top block will slide along the contour line of sawtooth at the beginning.en, once the force between sawteeth reaches the shear strength of sawtooth containing that of llings, the sawtooth will be cut o .In contrast, when the in lled ratio is larger than 1.0, the thickness of llings  Advances in Civil Engineering exceeds the height of sawtooth.In this situation, the top block will cut o the sawtooth directly instead of sliding along the contour line of sawtooth.is results in a significant reduction of in lled joints shear strength.erefore, the yield section cannot be observed in the stressdisplacement curve when the in lled ratio is larger than 1.0.A physical sketch as shown in Figure 9 was also drawn to re ect three di erent failure modes of in lled joints based on the shear stress-shear displacement curve under the in lled ratio of 0.25 and normal stress of 0.5 MPa.It can be seen that failure of the interface between rock and lling happens in the rst stage which can be seen as an elastic stage.Furthermore, rock slides along the interface which is the yielding stage.en, as the nal failure stage, the sawtooth will be cut o . is is because that the shear strength of in lled joint is smaller than that of llings at the beginning of shear motion.
It should be noted that when the in lled ratio equals to 1.0, which means the thickness of llings is the same as the height of sawtooth, there is no yield section in the stress-displacement curve.However, the shear strength of the in lled joint is larger than that of llings, as shown in Figure 8(a).is is because the compression of llings during shear process will induce the antishear ability of rock.Compared with the stressdisplacement curve under 0.5 MPa, the curve under 0.7 MPa only shows an increase in shear strength, and the other characteristics have few di erences.
As shown in Figure 8(b), with 0.3 MPa of initial normal stress, the top block slides along the contour line of sawtooth for each in lled ratio.For 0.5 MPa and 0.7 MPa of initial normal stress, the relative thick lling has a large compression deformation during shearing process.Under CNS condition, the normal stress applied on the in lled joint decreases due to shear compression, and correspondingly, the shear strength also shows a signi cant reduction as shown in Figure 8(c).

In uence of In lled Ratio on the Shear Stress Peak.
Figure 10 illustrates the relationship between the in lled ratio and peak of shear stress.
It can be seen from Figure 10 that if the initial normal stress keeps stable, the peak of shear stress of the in lled joint decreases gradually with the increase of the in lled ratio.When the in lled ratio is less than 1.0, the rock height within a sawtooth height decreases gradually with the increase of the in lled ratio which will bring a reduction in shear strength of the in lled joint correspondingly.When the in lled ratio is larger than 1.0, the shear strength of the in lled joint depends on the mechanical properties of lling, so the change of peak of shear stress is not obvious when the in lled ratio is larger than 1.0.Meanwhile, it can be seen that when the in lled ratio reaches 1.0, the shear stress peak decreases to the minimum value for the initial normal stress of 0.3 MPa.However, for the higher initial normal stress such as 0.5 MPa and 0.7 MPa, the critical value of the in lled ratio with which the shear stress peak reduces to minimum value rises to 1.25.
In conclusion, numerical results of the paper have a good agreement with the experimental results obtained by Indraratna, Papaliangas et al., and Goodman and the numerical results of Xu and Ren [12,[15][16][17].Numerical results in this paper show that there is a certain critical value for the in lled ratio.If the in lled ratio is less than the value, the shear strength of the in lled joint is controlled by both rock and llings and will decrease gradually with the growth of the in lled ratio.When the in lled ratio reaches the critical value, the shear strength of in lled joints is mainly in uenced by the mechanical properties of llings and only has little changes with the increase of the in lled ratio.Meanwhile, the results also indicate that this critical value of the in lled ratio is not necessarily to be 1.0 which will increase along with initial normal stress.is is presumably owing to the exertion of rock shear behaviors due to the compression deformation of llings.Based on a series of laboratory tests, Indraratna [12] also acquired the similar conclusion that the critical value is between 1.4 and 1.8, and the value of 1.4 corresponds to the larger initial normal stress.As for the Advances in Civil Engineering larger initial normal stress, the shear dilatancy of the infilled joint will be weakened gradually and even changes to be compression.Furthermore, it will eventually reduce the shear strength of the infilled joint.

Discussion
By comparing the numerical results with the previous experimental results, it can be found that the numerical model is able to conduct the direct shear test under CNS condition.Numerical results are able to reflect the corresponding shear behaviors of infilled joints under CNS condition.However, the numerical model still has several drawbacks.First of all, the infilled joints modeled in this paper are idealized as sawtooth shape for simplifying the research.erefore, the simulation of effects of irregular surface on the shear behaviors of natural joints is expected to be studied in future.Secondly, interface was established between the infilling and the upper rock block, while no interface was considered on the lower side of the infilling.
is is because firstly, lab measurements indicate that shear failure of rock joints typically occurs at one side of the fillings where the cohesion is smaller [28][29][30][31].Secondly, establishing interface elements on both sides of the fillings usually will lead to many errors during calculating and is computationally expensive.erefore, it is not necessary to build interface on both sides of the fillings in the model.Finally, the Mohr-Coulomb criterion used in the numerical model is not able to reflect the postpeak mechanical characteristics of rock.
erefore, the future study will focus on the model establishment of natural irregular joint, mechanical parameters determination, and representative constitutive models.Based on these studies, the numerical direct shear test is expected to reflect the shear behaviors of the infilled joint more veritably.

Conclusions
In this paper, the direct shear test under CNS condition for infilled joints was implemented in the numerical method.Based on the contrast tests, the paper studied the effects of different factors including normal constraints, undulating angle, and infilled ratio on the shear behaviors of infilled joints.A few conclusions are highlighted.Numerical simulation results reveal that the shear failure mode of infilled joints is influenced by initial normal stress, undulating angle, and infilled ratio.If the shear strength of interface between rock and filling is smaller than that of fillings, infilled joints will generate sliding failure along the interface.Otherwise, then the sawtooth will be cut off.e change of the shear failure model, that is, transforming from sliding failure to sawtooth cutoff failure, will induce the nonlinear distribution of the shear strength curve of infilled joints.is indicates that the numerical model shows a good agreement with experimental results available in literatures.Compared with traditional rock shear experiment, numerical simulation is also able to quantify the effect of multiple factors on shear properties of infilled joints such as initial normal stress and infilled ratio more effectively.erefore, this study proposed and verified an effective numerical analysis method capable of studying the effects of normal stress, undulating angles, and infilled ratio on the shear behavior of infilled rock joints.

Figure 1 :
Figure 1: (a) Concept model and (b) numerical model of the in lled joint.
(b), the length × width × height (without llings) of the model is 2.0 m × 0.2 m × 2.1 m. e height and bottom length of a sawtooth are 0.1 m and 0.4 m, respectively.

Figure 2 :
Figure 2: Comparison of the stress-strain curve obtained by the numerical model to that of Indraratna's experiment.

Figure 6 :
Figure 6: Shearing characteristics curve of the in lled joint with di erent undulating angles.(a) Shear stress-shear displacement curve.(b) Normal displacement-shear displacement curve.(c) Normal stress-shear displacement curve.

i 1 =Figure 7 :Figure 8 :
Figure 7: Shearing strength curve of the in lled joint with di erent undulating angles and friction angles.

Figure 8 :
Figure 8: Shearing characteristics curve of the in lled joint with di erent in lled ratios t/a.(a) Shear stress-shear displacement curve.(b) Normal displacement-shear displacement curve.(c) Normal stress-shear displacement curve.

Figure 9 :Figure 10 :
Figure 9: ree di erent shear failure stages of in lled joints corresponding to di erent failure modes.
Conditions. e bottom boundary and right boundary are xed in the z-axis and x-axis directions, respectively.Both the front boundary and back boundary are xed in the y-axis direction.A certain normal stress and CNS condition are applied on the top boundary of the model.ena very small shearing rate (4 × 10 −7 m/step) is applied on all the nodes in the top block.In this model, CNS condition is implemented through a servoprogram developed by the FISH code.e executing processes of the servoprogram are illustrated as follows.e whole shearing process is divided into many shearing cycles.In each cycle, the test runs 100 steps under the constant normal stress, σ n .

Table 1 :
Mechanical parameters of in lled joints.