• <tr id="yyy80"></tr>
  • <sup id="yyy80"></sup>
  • <tfoot id="yyy80"><noscript id="yyy80"></noscript></tfoot>
  • 99热精品在线国产_美女午夜性视频免费_国产精品国产高清国产av_av欧美777_自拍偷自拍亚洲精品老妇_亚洲熟女精品中文字幕_www日本黄色视频网_国产精品野战在线观看 ?

    Influence of layer orientation and interlayer bonding force on the mechanical behavior of shale under Brazilian test conditions

    2018-04-18 02:56:00JianmingHeLekanOlatayoAfolagboye
    Acta Mechanica Sinica 2018年2期

    Jianming He·Lekan Olatayo Afolagboye

    1 Introduction

    Studies have shown that clastic sedimentary deposits such as shale and clay are anisotropic in nature[1].Anisotropy can be either(1)inherent anisotropy that is a physical characteristic inherent in the material or(2)induced anisotropy that occurred due to the strain associated with the applied stresses.According to Ref.[2],“Shale is a clastic sedimentary rock mainly composed of silt-size and clay-size particles.Most shales are laminated and display fissility;the rock has a tendency to split along relatively smooth and flat surfaces parallel to the bedding”.Sequence of laminations occur as a result of variations in depositional environment that leads to differences grain size,clay percentage,microfossil content,organic material content,or mineral content[3].The presence of these laminar planes or laminations can have a serious effect on the mechanical behavior and strength of shale rock.This is because the lamination causes the shale rock to be anisotropic.The properties of any anisotropic rock are directionally dependent,which is different from isotropic rocks that have similar characteristics in all directions.Rock engineering activities such as hydraulic fracturing and fracture propagation,stability of underground and surface excavations,interpretation of micro-seismic monitoring,borehole stability and deformation,fluid flow and contaminant transportation are affected by rock anisotropy[4].Hence,there is proper need to understand the anisotropic behaviors of shale rock in the design of rock engineering works.

    Generally,the mechanical properties of anisotropic rocks are typically determined by testing the sample of anisotropic rock cut at varying orientation with respect to the planes of anisotropy.The testing method can be uniaxial compression,triaxial compression,Brazilian test,hollow cylinder test,resonant bar method,or ultrasonic pulse method[4].As remarked by Ref.[4],“the preparation of the test specimens,their dimensions,and the minimum number of test specimens necessary to determine the properties of anisotropic rocks depend on the method used”.Of all the mechanical properties that can be derived from the above-mentioned tests,uniaxial tensile strength is a vital parameter since rocks are naturally much weaker in tension than in compression.The Brazilian tensile test is the generally accepted laboratory method to determine the tensile strength of rock materials[5–7].The method of specimen preparation in the Brazilian tensile test is not strict compared with other methods,and the analysis of the test is not difficult.The International Society for Rock Mechanics[8]had suggested procedures for determining indirect tensile strength using this test.In the Brazilian test,a disc specimen is compressed with diametrically opposite and symmetric line loads.The justification for the test is based on the experimental fact that most rocks in biaxial stress fields fail in tension at their uniaxial tensile strength when one principal stress is tensile and the other finite principal stress is compressive. The layer orientation with respect to the axis of loading plays very important roles in the mechanical response of shale under tensile stress state.

    Many researchers have studied the tensile strength of anisotropic rocks using laboratory testing and numerical simulation methods.The analytical solution to measure the tensile strength of anisotropic rocks was developed by Ref.[9]based on the theoretical relation proposed by Ref.[10]between the stress and strain within a disc of anisotropic material under diametrical loading.The influences of layer orientation and layer thickness on indirect tensile strength and fracture pattern of anisotropic rocks such as sandstone and slate under diametrical loading conditions has been investigated by conducting Brazilian tests[11–15].The general conclusion from these works reveal that failure of anisotropic rocks are complicated as the tensile strength is related to the orientation of the plane of anisotropy and the failure pattern are complex.In their study,Lisjak et al.[16]used combined finite element method/discrete element method(FEM/DEM)to study the anisotropic mechanical behaviors of Opalinus Clay under Brazilian and compression tests.In the model developed by the authors,a transversely isotropic elastic constitutive law was used to describe the elastic response,and the rock fracturing was captured using discrete element method(DEM)algorithms and non-linear fracture mechanics principles.Work by Duan and Kwok[17]studied the behaviors of inherently anisotropic rocks under Brazilian test condition based on DEM simulations.In this study,DEM models were generated to represent inherently anisotropic rocks based on the bonded particle model(BPM)and smooth-joint model provided by particle flow code in two dimensions(PFC2D).Furthermore,Lanaro et.al.[18]used boundary element method code in two dimensions(FRACOD 2D),that implements the displacement discontinuity method(DDM),crack initiation,and propagation algorithms,to model Brazilian tensile tests.

    This work investigates the anisotropic behavior of transversely isotropic shale rocks under the Brazilian test method using experimental and numerical approaches,as well as the structure and composition of the laminations through digital microscopy and micro-computed tomography(micro-CT).Brazilian tests on the shale were carried out in the laboratory according to the suggested method by Ref.[8].PFC2Ddeveloped by Itasca was employed to simulate the mechanical behavior of shale under Brazilian tensile condition.Particle flow models of the shale were built according to the experimental investigation of the anisotropy characteristics.The paper presents the influence of orientation of the laminations and interlayer bonding force of these fine layers on the mechanical behavior of the shale.

    2 Anisotropy characteristics of the shale

    Shale samples used in this study are taken from an outcrop of fossil-bearing shale of the Longmaxi Formation in the Silurian group in ShiZhu County,Chongqin,China.The shale samples are very fine-grained,medium to dark grey,hard,and transversely isotropic.The samples were cored from a block with no macro-cracking.The anisotropy of the shale depends on the relative weakness of the laminar plane and the orientation of the plane with respect to the applied stress[19].The shale has dense laminations of weak planes,normally occurring as parallel planar structures(parallel laminations).Thin rock slices were made in order to show lamination structure and the shale slices were magnified using the digital microscope for the measurement of the thickness of the shale laminations.Figure 1 shows a magnified view of the laminations in the shale.The laminations have a thickness of approximately 1.0 mm.Miniature cubic size of the shale samples(2 mm×2 mm)were scanned using micro-CT and the scanned image of the shale is shown in Fig.2.The mineral content of the laminations was investigated by tracing the boundaries of the different minerals in the micro-view of the shale and are listed in Table 1.Laminations in the shale consist of small differences in the type of sediment that occur in different grain size and mineral content.The distribution of the constituent minerals shows the difference in composition of the shale laminations as indicated by red line.The presence of relatively weak laminar plane can result in a remarkable anisotropic effect on the mechanical behavior and strength of the shale.

    Fig.1 Magnified view of the shale laminations using digital microscope

    Fig.2(Color online)Scanned image of the shale sample using micro-CT

    3 Brazilian testing and results

    The shale samples taken from the outcrop were cut into discs with a diameterDof 50 mm and the thickness of approximately 25 mm.The two ends of the cores were ground in order to meet the standard specification.The axis of the cores was set parallel to the strike of the laminar planes.The orientation angle of the lamination is defined as the angle between the laminar plane and the horizontal direction.Brazilian test was carried out in the laboratory according to the specification of Ref.[8].Two steel loading jaws were designed so as to contact the disc-shaped rock sample at diametrically opposed surfaces over an arc of contact of approximate 10°at failure(Fig.3).The triaxial rock mechanical test system TAW-2000 was used for the Brazilian tests.The rock disc was loaded continuously at a constant loading rate of200 N/s until failure.The orientation angle of the laminations of the specimen was changed from 0°to 90°at an interval of 15°.Five specimens were tested for each orientation angle since one sample may not be adequate to represent the general failure behavior corresponding to a lamination orientation.The tensile strength(Table 2)for each specimen was calculated as

    Fig.3 Brazilian test set-up with sample and loading jaws

    wherePis the splitting load when the rock disc fails;Dis the diameter;Bis the thickness of the rock disc.

    Variation in the Brazilian tensile strength(BTS)with the inclination angle of the lamination is plotted in Fig.4.The graph shows a downward trend that is the failure strength decrease with increase in orientation of the lamination.How-ever,the tensile strength increases a little when the layer orientation increased from 75°to 90°.In addition,the average failure strength of the samples with lamination angle less than 30°(6.46 MPa)was larger compared with samples with lamination angle greater than 30°(3.46 MPa).The trend of decreasing failure strength with increase in orientation angle is similar to the experimental results of Refs.[4,12,20]on transversely isotropic sandstone.Tavallali and Vervoort[12]explained this phenomenon as follows:“when the layers are horizontal or semi-horizontal the fracture is mainly through the stronger material,while by increasing the inclination angle the fracturing processes make use of the layers,of which one could expect that they have weaker mechanical properties.In addition,it is possible that in the latter cases samples fail also in shear and not purely in tension.”Similar to the experimental results of Ref.[4],the difference between the maximum and minimum value of strength was significant(4.16 MPa).Tavallali and Vervoort[12],however,observed a relatively small difference between the maximum and minimum value of strength.

    Table 2 Tensile strength of the shale specimens

    Table 1 Minerals present in the laminations and the corresponding content

    Different modes of fracture propagation(Table 3)led to the failure of different shale specimens,and more than one failure mode can be observed in specimens with equal lamination orientation.The following types of fracture propagation were observed after specimen failure:

    Fig.4 Tensile strength of the laboratory results in comparison with the results of simulation

    (1)Fracture propagation that deviated from the loading direction of the rock disks termed as curved fracture;

    (2)Propagation partly along and partly through the laminar plane termed broken-linear fracture;

    (3)Fractures that propagated along the laminar plane and termed as layer-activated fracture;

    (4)Fracture that passed roughly along the loading direction of the rock disks termed central-linear fracture.

    Curved fracture propagated in tension mode along the rock matrix.Layer-activated or broken-linear fracture propagated in mixed mode that is shearing along the laminar plane and tension in the rock matrix.Central-linear/curved fracture mainly propagated in tensile mode along the laminar plane when the orientation was 0°and in shear mode when the orientation was 90°.Some other forms of fractures werealso observed apart from the main failure fracture especially in samples with low lamination orientations.The fracture pattern is transformed from curved to broken-linear,then layer-activated,and finally central-linear fracture as the layer orientation increases.Curved fractures were commonly observed when the inclination angle was between 0°and 30°.Layer-activated and broken linear fractures were observed when the inclination angle was between 30°and 60°.Central linear fractures were observed when the inclination angle was between 75°and 90°.This indicates that the laminar structures in the shale play an important role with the increase in layer orientation.The rock matrix is replaced gradually by the laminar planes as the main controlling factor for the tensile failure process of shale.

    Table 3 Some typical fracturing patterns of the shale specimens with different layer orientation

    4 Numerical modeling

    DEM has been widely used in civil,mining,petroleum,and geo-environmental engineering applications as it is recognized as a tool capable of simulating the progressive failure of rock and micro-mechanisms underlying the deformation behavior of rocks.It has been used to study the role of grain interlocking on the strength of a rock[21],effect of porosity on the deformability and strength[22],the impacts of fluid viscosity,particle size distribution,and the brittleness of shale on fracture propagation[23].DEM has a greater advantage over other continuum methods because it can simulate crack initiation and propagation from microscale to macro-scale without applying complex constitutive laws[24].

    In this study,PFC2Ddeveloped by Itasca was used for the numerical modelling of the mechanical behavior of the shale under Brazilian test conditions.The two types of bonded particle model(parallel bond model and contact bond model)found in the program were used.In PFC2Dmodels,rocks are represented as an assembly of discs bonded at their point of contacts[25].The overall constitutive behavior of shale is simulated as an assembly of particles bonded to one another at their contact points,and the breakage of individual bonds represents damage.The stress will be redistributed when a bond breaks and this may cause further bond breakage.The induced micro-cracks can coalesce and later form macro fractures.The numerical modelling is based on the laboratory results of the anisotropy characteristics under Brazilian test conditions.The 2D numerical models used in this study were established in strict accordance with the scanning results of the different laminations by micro-CT.The layer orientation,interlayer bonding force,and mineral composition difference between layers were considered in building the numerical model for studying the mechanical behavior of shale as influenced by layer orientation and interlayer bonding force under the Brazilian tensile conditions.

    The specimens used for the numerical simulation of the indirect tensile tests is a disc shape with a diameter of50 mm and unit thickness.The total number of particle is 42,102 and the particle size follows a uniform distribution withRmin=0.15 mm andRmax/Rmin=2.The particle distribution of the specimens with different inclination angles were the same.The diameter of the rock specimen was 50 times that of the laminations and the spacing between each laminae plane was about 1.0 mm.The parallel-bond model was applied to rock matrix and laminar plane,defined as an interface between two neighboring laminations,was assigned with the contact bond along which sliding and separation can occur.A contact bond is defined by the following two parameters:normal contact bond strength and shear contact bond strength.Different bond strength values were given to the rock matrix and the lamination layers to simulate the strength difference between the layers.Different microparameters were assigned to different zones,corresponding to layer and matrix,in the specimen to study the different orientation angles of the lamination planes.Parameters for the parallel bond used in the rock matrix and contactbond used in the laminar planes are listed in the Tables 4 and 5,respectively.The material strengths are picked from a Gaussian(normal)distribution for both shale matrix and laminar planes,which are specified in terms of mean and a standard deviation.A computer program based on the finite element method was employed to estimate the lamination stiffness under the condition of uniaxial compression.According to the scanned image of the shale laminations in Fig.2,different minerals,pores,and micro-cracks were considered in the numerical model and the weakness due to cementation that exist between the mineral grains were ignored in the simulation.The uniaxial compression test was simulated in the program on the neighboring laminations respectively and the different lamination stiffness was obtained.The ratio of stiffness(1.114)between the neighboring laminations was used instead of the absolute value in the subsequent modeling of the shale using PFC2D.In this numerical simulation,the anisotropy angles considered were 0°,15°,30°,45°,60°,75°,90°.The interlayer bonding force was adjusted through modifying the shear and tensile strength of the contact bonds based on the simulation,which can compare with the laboratory results.The influence of interlayer bonding force was studied through modifying the shear and tensile strength of the contact bonds to the integral multiple of the initial value from×1 to×8.

    Table 4 Micro-parameters of the parallel bond to define the rock matrix

    Micro-cracks that form between the bonded particles at peak stress level during the BTS test were counted and used to determine the mechanism that controlled the failure of each specimen.In PFC,micro cracks that can form when(1)the normal strength of contact bond(CN)has been exceeded;(2)the shear strength of contact bond(CS)has been exceeded;(3)the normal strength of parallel bond(PN)has been exceeded;(4)the shear strength of parallel bond(PS)has been exceeded.The number of micro-tensile and micro-shear cracks in the bonds of the parallel bond model and contact bond model were used determine the mechanism by which the models fail.The bond breaks in shearing and normal directions are plotted in red and white color respectively while the particles are plotted in black color.Some typical simulation results of the micro cracks distribution with different layer orientation are shown in Table 6.

    Table 5 Micro parameters of the contact bond that define the laminar plane

    4.1 Comparison of the tensile results

    The variation in numerical BTS with the inclination angle of the lamination is plotted in Fig.4.As shown in Fig.4,laboratory result and simulation results show a good agreement.Apart from when the orientation of the laminar plane is 0°,the BTS of the numerical simulation is lower than the average BTS obtained from the laboratory test in other orientation angles.The numerical simulation also shows that the numerical BTS were very sensitive to the variation in the layer orientation.In addition,similar fracture patterns were observed in both the laboratory and simulation results(Table 6).Some of the typical fracture patterns observed in the numerical simulation tests are shown in row×1 in Table 6.The curved fracture pattern is observed when the layer orientation is 0°.Overall,the anisotropic behavior of the model specimens becomes more pronounced with the increase in layer orientation from 0°to 90°.The fracture pattern changes from curved to the broken-linear,then to layer-activated,and finally central-linear fracture.This indicates that the failure through rock matrix is gradually replaced by the laminar planes failure as the orientation angle increases.In addition,the specimens with smaller layer orientation angles always failure with zigzagged fractures.

    The effects of layer orientation on the micro-cracks formed during the loading process in the numerical simulation are shown in row×1,Table 6.The distribution of the micro-cracks in specimens with layer orientation less than 30°were relatively symmetrical and uniform.Different deformability of the neighboring laminations result in the shearing stress accumulation in the laminar planes,thus the micro-cracks were formed due to the break of the contact bond.With increase in layer orientation,the distribution of micro-cracks became asymmetrically and non-uniformly.The ratio of micro-cracks due CS to micro-cracks due to CN showed a decreasing trend.When the layer orientation exceeded 45°,the ratio dropped to less than 1(Fig.5).This shows that the shearing breaks has transformed to the normal breaks in the laminar planes.The ratio of micro-cracks due PN to micro-cracks due to PS was more than 5 in all layer orientations(Fig.6).This indicates that bond breaks in the normal direction due to the tensile stress dominate the micro-cracks formed in the rock matrix.Most of the ratios of the breaks of PN and PS are more than 5.0,and as shown in Fig.6,this indicates that bond breaks in the normal direction due to the tensile stress dominate the micro-cracks formed in the rock matrix.

    Table 6 Some typical micro-crack distributions of the simulation with different layer orientation

    Fig.5 Variation of the ratio of the breaks of CS/CN against layer orientation

    Fig.6 Variation of the ratio of the breaks of PN/PS against layer orientation

    4.2 Influence of the interlayer bonding force

    Special attention was paid to the influence of interlayer bonding force on the fracturing of shale under Brazilian tensile conditions.Because the rock specimens were taken from the outcrop,the bonding conditions of the lamination may vary due to heterogeneity introduced by weathering.Considering the bonding conditions of the rock disks,study on the influence of interlayer bonding force was achieved by using numerical simulation.As shown in Fig.7,the tensile strength of the different interlayer bonding forces against the layer orientations changed significantly.The layer orientation and interlayer bonding force of the shale play a very important role in the anisotropic mechanical behavior of shale.The amplitude of strength variation is much smaller in cases with layer orientation less than 30°(3.05 MPa)than cases with layer orientation more than 30°(6.44 MPa).The numerical tensile strength shows a decreasing trend in general when the interlayer bonding strength was set smaller(?×4).

    Table 6(row ×3 to×7)shows the influence of interlayer bonding force on the micro-cracks formed during the loading process.The total percentage of micro-cracks due to CN and CS in laminar planes decreased as the interlayer bonding force increased in all the layer orientation(Fig.8).In contrast,the total percentage of micro-cracks due to PN and PS increased as the interlayer bonding force increased(Fig.9).The number of micro-cracks decreased sharply with the increasing of the interlayer bonding strength from×1 to×7(Table 6).This indicates that the increase in interlayer bonding strength leads to reduction in the anisotropic behavior of the shale.The bonding force distributions in the specimens corresponding to the micro-cracks simulation results(Table 7)also indicates that the anisotropic mechanical response is higher in the case of lower interlayer bonding strength(×1 and×3)compared to case of higher interlayer bonding strength(×5 and×7).

    Fig.7 Tensile strength against layer orientation of the simulation results with increasing inter layer bonding strength

    Fig.8(Color online)Variation in the percentage of micro-cracks of CN+CS in the total

    Fig.9(Color online)Variation in the percentage of micro-cracks of PN+PS in the total

    5 Conclusions

    Based on experimental studies and numerical simulations,the mechanical behavior of inherently anisotropic shale rocks under Brazilian test conditions is investigated in this study.The effects of the weak lamination planes and interlayer bonding force of these layers on the failure strength and fracture patterns are examined.The laboratory failure strength decreased with increase in layer orientation angle.Four types of fracture patterns were observed:curved fracture,broken-linear fracture,layer-activated fracture,and central linear fracture.The fracture pattern is transformed from curved to broken-linear,then to layer-activated,and finally to central-linear fracture as the layer orientation increases.The observed fracture propagated in tensile and/or shear mode.The numerical simulations showed a good agreement with the experimental results in the failure strength and fracture patterns.The numerical tensile strength shows a decreasing trend in general when the interlayer bonding strength were set smaller(?×4).The distribution of the micro-cracks in specimens with layer orientation less than 30°were relatively symmetrical and uniform.The quantity of micro-cracks decreased sharply with increase in the interlayer bonding strength and this directly led to reduction in the anisotropic behavior in the mechanical properties of the shale.The failure strength,fracture pattern,micro-cracks distribution,and contact force in the specimens show obvious anisotropic characteristics in both of the laboratory and simulation results.The layer orientation and interlayer bonding force of the shale thus play a very important role in the anisotropic behavior of shale.

    Table 7 Some typical bonding force distribution of the simulations with different layer orientation

    AcknowledgementsThe project was supported by the National Natural Science Foundation of China(Grants 41572310,41272351,and 41227901)and the Strategic Priority Research Program of the Chinese Academy of Sciences(Grants XDB10030301 and XDB10030304).Authors acknowledge their sincere thanks to Prof.Jijin Yang and Prof.Lihui Li for their helpful microscopic information on the shale.

    1.Mitchell,J.K.,Soga,K.:Fundamentals of Soil Behavior.Wiley,New York(2005)

    2.Terzagi,K.,Peck,R.B.,Gholamreza,M.:Soil Mechanics in Engineering Practice.Wiley,New York(1996)

    3.Boggs,S.J.:Principles of Sedimentology and Stratigraphy.Merrill Publishing Company,Westerville(1987)

    4.Chen,C.-S.,Pan,E.,Amadei,B.:Determination of deformability and tensile strength of anisotropic rock using Brazilian tests.Int.J.Rock Mech.Min.Sci.35,43–61(1998)

    5.Andreev,G.E.:A review of the Brazilian test for rock tensile strength determination.Part II:contact conditions.Min.Sci.Technol.13,457–465(1991)

    6.Andreev,G.E.:A review of the Brazilian test for rock tensile strength determination.PartI:calculation formula.Min.Sci.Technol.13,445–456(1991)

    7.Perras,M.A.,Diederichs,M.S.:A review of the tensile strength of rock:concepts and testing.Geotech.Geol.Eng.32,525–546(2014)

    8.ISRM:Suggested methods for determining tensile strength of rock materials.Int.J.Rock Mech.Min.Sci.15,99–103(1978)

    9.Amadei,B.,Rogers,J.D.,Goodman,R.E.:Elastic Constants and Tensile Strength of Anisotropic Rocks.In:5th ISRM Congress,ISRM-5CONGRESS-1983-030.Melbourne,Australia(1983)

    10.Lekhnitskii,S.G.:Anisotropic Plates.Gordon and Breach,New York(1968)

    11.Chen,C.,Pan,E.,Amadei,B.:Determination of deformability and tensile strength of anisotropic rock using Brazilian tests.Int.J.Rock Mech.Min.Sci.Geomech.Abstr.35,43–61(1998)

    12.Tavallali,A.,Vervoort,A.:Failure of transversely isotropic rock material:effect of layer orientation and material properties.In:The 6th International Symposium on Ground Support in Mining and Civil engineering Construction,317–328.Cape Town,South Africa(2008)

    13.Vervoort,A.,Tavallali,A.:Effect of layerorientation on the failure of layered sandstone under Brazilian test conditions.Int.J.Rock Mech.Min.Sci.47,313–322(2010)

    14.Tavallali,A.,Vervoort,A.:Failure of layered sandstone under Brazilian test conditions:effect of micro-scale parameters on macro-scale behaviour.Rock Mech.Rock Eng.43,641–653(2010)

    15.Tan,X.,Konietzky,H.,Frühwirt,T.,et al.:Brazilian tests on transversely isotropic rocks:laboratory testing and numerical simulations.Rock Mech.Rock Eng.48,1341–1351(2015)

    16.Lisjak,A.,Tatone,B.S.A.,Grasselli,G.,et al.:Numerical modelling of the anisotropic mechanical behaviour of opalinus clay at the laboratory-scale using FEM/DEM.Rock Mech.Rock Eng.47,187–206(2014)

    17.Duan,K.,Kwok,C.Y.:Discrete element modeling of anisotropic rock under Brazilian test conditions.Int.J.Rock Mech.Min.Sci.78,46–56(2015)

    18.Lanaro,F.,Sato,T.,Stephansson,O.:Microcrack modelling of Brazilian tensile tests with the boundary element method.Int.J.Rock Mech.Min.Sci.46,450–461(2009)

    19.Mark,D.Z.:Reservoir Geomechanics.Cambridge University Press,Cambridge(2007)

    20.Tavallali,A.,Vervoort,A.:Behaviour of layered sandstone under Brazilian test conditions:layer orientation and shape effects.J.Rock Mech.Geotech.Eng.5,366–377(2013)

    21.Scholtès,L.,Donzé,F.-V.:A DEM model for soft and hard rocks:role of grain interlocking on strength.J.Mech.Phys.Solids 61,352–369(2013)

    22.Weng,M.C.,Li,H.H.:Relationship between the deformation characteristicsand microscopic propertiesofsandstone explored by the bonded-particle model.Int.J.Rock Mech.Min.Sci.56,34–43(2012)

    23.Shimizu,H.,Murata,S.,Ishida,T.:The distinct element analysis for hydraulic fracturing in hard rock considering fluid viscosity and particle size distribution.Int.J.Rock Mech.Min.Sci.48,712–727(2011)

    24.Cundall,P.A.:A discontinuous future for numerical modelling in geomechanics?Proc.Inst.Civ.Eng.Geotech.Eng.149,41–47(2001)

    25.Potyondy,D.O.,Cundall,P.A.:A bonded-particle model for rock.Int.J.Rock Mech.Min.Sci.41,1329–1364(2004)

    日韩av在线大香蕉| 每晚都被弄得嗷嗷叫到高潮| 免费在线观看日本一区| 久9热在线精品视频| 在线观看日韩欧美| 在线看三级毛片| 国产单亲对白刺激| 成人三级做爰电影| 久久久久久久久免费视频了| 精品国内亚洲2022精品成人| 国产亚洲精品一区二区www| 超碰成人久久| 女同久久另类99精品国产91| 成人av一区二区三区在线看| 久久欧美精品欧美久久欧美| 午夜福利高清视频| 国产精品永久免费网站| 久久婷婷成人综合色麻豆| 国产高清videossex| 精品一区二区三区四区五区乱码| 久久婷婷成人综合色麻豆| 欧美乱色亚洲激情| 国产成人av激情在线播放| 亚洲国产精品成人综合色| 哪里可以看免费的av片| 美女午夜性视频免费| 神马国产精品三级电影在线观看 | 精品无人区乱码1区二区| 一区二区三区精品91| 99国产极品粉嫩在线观看| 欧美绝顶高潮抽搐喷水| 国产精品永久免费网站| 国产成人一区二区三区免费视频网站| 非洲黑人性xxxx精品又粗又长| 欧美日韩乱码在线| 国产精品日韩av在线免费观看| 亚洲三区欧美一区| 欧美黄色淫秽网站| 午夜老司机福利片| 久久久久精品国产欧美久久久| 人人妻人人看人人澡| www日本在线高清视频| 亚洲专区字幕在线| 香蕉久久夜色| 亚洲自偷自拍图片 自拍| 免费一级毛片在线播放高清视频| 亚洲午夜理论影院| 动漫黄色视频在线观看| 波多野结衣高清无吗| 一区福利在线观看| 国产黄色小视频在线观看| 色精品久久人妻99蜜桃| 99久久综合精品五月天人人| netflix在线观看网站| 国产1区2区3区精品| 亚洲国产精品合色在线| 亚洲一区二区三区色噜噜| 精品国产美女av久久久久小说| 国产三级黄色录像| 搡老岳熟女国产| 天堂动漫精品| 精品久久久久久久末码| 精品久久久久久久末码| 99久久综合精品五月天人人| 国内精品久久久久精免费| 老熟妇乱子伦视频在线观看| 我的亚洲天堂| 亚洲一卡2卡3卡4卡5卡精品中文| 亚洲专区中文字幕在线| 淫妇啪啪啪对白视频| 三级毛片av免费| 国产在线观看jvid| 制服人妻中文乱码| 亚洲国产欧美网| 国产真实乱freesex| 精品少妇一区二区三区视频日本电影| 精品国产乱码久久久久久男人| 久久久久久九九精品二区国产 | а√天堂www在线а√下载| 波多野结衣高清作品| 欧美黑人巨大hd| 成人国语在线视频| 丁香六月欧美| 黄色成人免费大全| 成人av一区二区三区在线看| 久久国产精品影院| 亚洲成国产人片在线观看| 99热这里只有精品一区 | 亚洲色图av天堂| 成人av一区二区三区在线看| 大型黄色视频在线免费观看| 国产亚洲精品第一综合不卡| 亚洲精品久久成人aⅴ小说| 18禁国产床啪视频网站| 久久久久久久久久黄片| 国产精品久久视频播放| 成人三级做爰电影| 欧美 亚洲 国产 日韩一| 麻豆久久精品国产亚洲av| 国产av一区二区精品久久| 欧美在线一区亚洲| 亚洲欧美精品综合久久99| 日韩国内少妇激情av| 波多野结衣高清无吗| 亚洲中文av在线| 视频在线观看一区二区三区| 久久草成人影院| 热re99久久国产66热| 亚洲熟女毛片儿| 亚洲国产精品sss在线观看| 一进一出好大好爽视频| 亚洲自偷自拍图片 自拍| 国产一级毛片七仙女欲春2 | 欧美黑人巨大hd| 国内精品久久久久久久电影| 成人18禁在线播放| 亚洲国产日韩欧美精品在线观看 | 国产一区二区激情短视频| 亚洲成av人片免费观看| 成人三级做爰电影| 国产成人精品久久二区二区免费| 亚洲电影在线观看av| 亚洲精品一卡2卡三卡4卡5卡| 在线永久观看黄色视频| 一个人免费在线观看的高清视频| 好男人在线观看高清免费视频 | 变态另类丝袜制服| 男人的好看免费观看在线视频 | 免费看日本二区| 亚洲成a人片在线一区二区| 国产欧美日韩一区二区精品| 国产精品 欧美亚洲| 国产伦在线观看视频一区| 久久婷婷成人综合色麻豆| 日韩成人在线观看一区二区三区| 日韩欧美一区二区三区在线观看| 久久久久久大精品| 久久久久九九精品影院| 亚洲国产中文字幕在线视频| 搡老熟女国产l中国老女人| 12—13女人毛片做爰片一| 午夜福利欧美成人| 国产精品自产拍在线观看55亚洲| 18禁裸乳无遮挡免费网站照片 | 一级a爱视频在线免费观看| 精品久久久久久,| 中文字幕人妻丝袜一区二区| bbb黄色大片| 在线十欧美十亚洲十日本专区| 国产又色又爽无遮挡免费看| 免费在线观看黄色视频的| 亚洲专区字幕在线| 亚洲av成人av| www.精华液| 麻豆久久精品国产亚洲av| av在线天堂中文字幕| 午夜福利成人在线免费观看| 动漫黄色视频在线观看| 亚洲精品中文字幕一二三四区| 精品久久久久久久久久久久久 | 麻豆av在线久日| 怎么达到女性高潮| 久久精品夜夜夜夜夜久久蜜豆 | xxxwww97欧美| 一二三四在线观看免费中文在| 亚洲国产欧洲综合997久久, | 国产一卡二卡三卡精品| 欧美人与性动交α欧美精品济南到| 色播在线永久视频| 男女下面进入的视频免费午夜 | 成年版毛片免费区| 久久香蕉国产精品| 亚洲男人的天堂狠狠| 亚洲va日本ⅴa欧美va伊人久久| 黄网站色视频无遮挡免费观看| 国产爱豆传媒在线观看 | 久久国产精品人妻蜜桃| 久久久久久久久免费视频了| 日本成人三级电影网站| 国产又色又爽无遮挡免费看| 精品免费久久久久久久清纯| x7x7x7水蜜桃| 色播在线永久视频| 18禁观看日本| 国产午夜精品久久久久久| 成人手机av| 精品久久久久久久久久久久久 | 在线免费观看的www视频| 亚洲五月婷婷丁香| 欧美大码av| 夜夜看夜夜爽夜夜摸| 国产私拍福利视频在线观看| 91av网站免费观看| 国内精品久久久久精免费| 在线视频色国产色| 别揉我奶头~嗯~啊~动态视频| 国产高清视频在线播放一区| 亚洲男人天堂网一区| 免费在线观看黄色视频的| 曰老女人黄片| 90打野战视频偷拍视频| 正在播放国产对白刺激| 午夜福利18| 视频在线观看一区二区三区| 99久久99久久久精品蜜桃| 18禁黄网站禁片免费观看直播| 久久久久国内视频| 成人永久免费在线观看视频| 长腿黑丝高跟| 国产精品永久免费网站| 搡老岳熟女国产| 99在线视频只有这里精品首页| 欧美精品啪啪一区二区三区| 久久久精品欧美日韩精品| 伦理电影免费视频| 中文字幕人成人乱码亚洲影| 日本一区二区免费在线视频| 国产午夜福利久久久久久| 丝袜美腿诱惑在线| 国产在线精品亚洲第一网站| 日韩精品免费视频一区二区三区| 欧美黄色淫秽网站| 成人国产一区最新在线观看| 亚洲五月婷婷丁香| 亚洲欧美精品综合一区二区三区| 国内毛片毛片毛片毛片毛片| 日日夜夜操网爽| 精品国产一区二区三区四区第35| e午夜精品久久久久久久| 久久精品91蜜桃| 少妇裸体淫交视频免费看高清 | 中文字幕精品免费在线观看视频| 久久国产精品男人的天堂亚洲| 黄片大片在线免费观看| 曰老女人黄片| 男女之事视频高清在线观看| 中文字幕av电影在线播放| 国产又色又爽无遮挡免费看| 成人国语在线视频| 欧美又色又爽又黄视频| 又大又爽又粗| 亚洲最大成人中文| 精品国产乱码久久久久久男人| 精品久久久久久久久久免费视频| 亚洲片人在线观看| 免费在线观看成人毛片| 亚洲成人久久性| 免费人成视频x8x8入口观看| 国产成人精品久久二区二区免费| 黄色视频不卡| 50天的宝宝边吃奶边哭怎么回事| 久久精品国产亚洲av香蕉五月| 国产视频内射| ponron亚洲| 亚洲五月天丁香| 美国免费a级毛片| 久久久水蜜桃国产精品网| 午夜成年电影在线免费观看| 日韩 欧美 亚洲 中文字幕| 后天国语完整版免费观看| 757午夜福利合集在线观看| 哪里可以看免费的av片| 黄色a级毛片大全视频| 国产成人系列免费观看| 最新在线观看一区二区三区| 久久精品人妻少妇| av福利片在线| 悠悠久久av| av视频在线观看入口| 真人做人爱边吃奶动态| 国产午夜精品久久久久久| 国产高清videossex| 99国产极品粉嫩在线观看| 精品欧美国产一区二区三| 亚洲av五月六月丁香网| 国产成人精品无人区| 国产1区2区3区精品| 国产精品免费一区二区三区在线| 日本撒尿小便嘘嘘汇集6| 两性午夜刺激爽爽歪歪视频在线观看 | 色精品久久人妻99蜜桃| 亚洲最大成人中文| 51午夜福利影视在线观看| 免费一级毛片在线播放高清视频| 99在线人妻在线中文字幕| 国产精品永久免费网站| 婷婷六月久久综合丁香| 精品欧美一区二区三区在线| 精品人妻1区二区| 亚洲色图 男人天堂 中文字幕| 免费在线观看视频国产中文字幕亚洲| 天堂√8在线中文| √禁漫天堂资源中文www| 老熟妇仑乱视频hdxx| 亚洲无线在线观看| 免费看日本二区| 国产精品98久久久久久宅男小说| 免费在线观看完整版高清| 午夜福利视频1000在线观看| 国产亚洲精品av在线| 国产成人欧美| 亚洲九九香蕉| 色尼玛亚洲综合影院| 久热爱精品视频在线9| 亚洲国产精品999在线| 国产精品精品国产色婷婷| 国内精品久久久久久久电影| 欧美最黄视频在线播放免费| 国产高清视频在线播放一区| 亚洲精品久久成人aⅴ小说| 丝袜美腿诱惑在线| 91麻豆av在线| 黑人操中国人逼视频| 男人舔奶头视频| 国产精品久久视频播放| 国产aⅴ精品一区二区三区波| www.自偷自拍.com| 久久久久久大精品| 久久久久久大精品| 国产精品一区二区精品视频观看| 91字幕亚洲| 久久精品aⅴ一区二区三区四区| 亚洲av电影在线进入| 欧美日韩一级在线毛片| 国产精品久久电影中文字幕| 欧美三级亚洲精品| 亚洲一码二码三码区别大吗| 在线播放国产精品三级| www.自偷自拍.com| av视频在线观看入口| 夜夜躁狠狠躁天天躁| 午夜福利在线在线| 99国产精品一区二区蜜桃av| 免费观看精品视频网站| 美女大奶头视频| 欧美成狂野欧美在线观看| 热re99久久国产66热| 久热这里只有精品99| 久久香蕉激情| 亚洲va日本ⅴa欧美va伊人久久| 88av欧美| 国产成人欧美在线观看| 亚洲国产欧美日韩在线播放| 免费女性裸体啪啪无遮挡网站| 亚洲人成网站高清观看| 女人高潮潮喷娇喘18禁视频| 18禁裸乳无遮挡免费网站照片 | 中文亚洲av片在线观看爽| 在线永久观看黄色视频| aaaaa片日本免费| 精品福利观看| 色精品久久人妻99蜜桃| 亚洲成人免费电影在线观看| 成年人黄色毛片网站| 亚洲专区国产一区二区| 十八禁人妻一区二区| 欧美日韩精品网址| 久久亚洲真实| 欧美日韩精品网址| 日韩高清综合在线| 又紧又爽又黄一区二区| 18美女黄网站色大片免费观看| 黄色成人免费大全| 午夜福利在线观看吧| 欧美日韩精品网址| a在线观看视频网站| 亚洲三区欧美一区| 欧美色视频一区免费| 一区二区三区激情视频| 麻豆国产av国片精品| 亚洲国产欧美日韩在线播放| 亚洲一卡2卡3卡4卡5卡精品中文| 嫩草影视91久久| 久久香蕉激情| 欧美又色又爽又黄视频| 免费看日本二区| 精华霜和精华液先用哪个| 精品日产1卡2卡| 色在线成人网| 校园春色视频在线观看| 久久天堂一区二区三区四区| 九色国产91popny在线| 亚洲精品色激情综合| 夜夜躁狠狠躁天天躁| 亚洲熟妇熟女久久| 中文字幕高清在线视频| 国产精品av久久久久免费| 欧美乱码精品一区二区三区| 久久精品国产亚洲av高清一级| 久久久久国产精品人妻aⅴ院| 国产精品av久久久久免费| 中文字幕精品免费在线观看视频| 久久久久久久精品吃奶| 黑人巨大精品欧美一区二区mp4| 一本久久中文字幕| 亚洲中文字幕一区二区三区有码在线看 | 亚洲片人在线观看| 观看免费一级毛片| 19禁男女啪啪无遮挡网站| 50天的宝宝边吃奶边哭怎么回事| 久久久久久久午夜电影| 成人精品一区二区免费| 国产亚洲精品久久久久5区| 亚洲最大成人中文| 久久久久国内视频| 99久久久亚洲精品蜜臀av| 黄色片一级片一级黄色片| 国产精品久久久久久亚洲av鲁大| 老司机午夜福利在线观看视频| 午夜福利欧美成人| 最新美女视频免费是黄的| 亚洲国产欧洲综合997久久, | 99国产综合亚洲精品| 91九色精品人成在线观看| 亚洲人成77777在线视频| 波多野结衣巨乳人妻| 国产精品日韩av在线免费观看| 91在线观看av| 日日干狠狠操夜夜爽| 99国产精品一区二区三区| 夜夜爽天天搞| 欧美日韩黄片免| 欧美日韩一级在线毛片| 亚洲自拍偷在线| 亚洲精品国产一区二区精华液| 一进一出好大好爽视频| 窝窝影院91人妻| 久久久精品欧美日韩精品| 又黄又粗又硬又大视频| 久久天躁狠狠躁夜夜2o2o| 久久久久久免费高清国产稀缺| 亚洲色图 男人天堂 中文字幕| 男人舔女人的私密视频| 很黄的视频免费| 精品一区二区三区视频在线观看免费| 午夜免费激情av| 一区二区三区高清视频在线| 国产亚洲av嫩草精品影院| 国产精品98久久久久久宅男小说| 一边摸一边抽搐一进一小说| 亚洲精品在线美女| 久久这里只有精品19| 丁香欧美五月| 国产激情久久老熟女| 国产成人精品久久二区二区91| 此物有八面人人有两片| 欧美激情高清一区二区三区| 夜夜躁狠狠躁天天躁| 男人舔女人下体高潮全视频| 在线观看免费午夜福利视频| 美女国产高潮福利片在线看| 国产成年人精品一区二区| 成年女人毛片免费观看观看9| 亚洲精品在线观看二区| 国产精华一区二区三区| 99久久99久久久精品蜜桃| 国产又爽黄色视频| 啪啪无遮挡十八禁网站| 亚洲午夜理论影院| 黄色女人牲交| 久久香蕉激情| 一本综合久久免费| 国产高清激情床上av| 免费人成视频x8x8入口观看| 久99久视频精品免费| 不卡av一区二区三区| 久久久国产精品麻豆| 校园春色视频在线观看| 每晚都被弄得嗷嗷叫到高潮| 午夜亚洲福利在线播放| 精品卡一卡二卡四卡免费| 啦啦啦 在线观看视频| 国产人伦9x9x在线观看| 在线视频色国产色| 精品欧美国产一区二区三| 亚洲av成人一区二区三| 亚洲精品色激情综合| 免费看日本二区| 女人高潮潮喷娇喘18禁视频| 无人区码免费观看不卡| 日韩欧美 国产精品| 久久久久国内视频| 婷婷亚洲欧美| 1024香蕉在线观看| 在线看三级毛片| 亚洲熟女毛片儿| 欧美zozozo另类| 成人永久免费在线观看视频| 长腿黑丝高跟| 国产aⅴ精品一区二区三区波| 亚洲美女黄片视频| 无人区码免费观看不卡| 成人av一区二区三区在线看| 中亚洲国语对白在线视频| 脱女人内裤的视频| 久久狼人影院| 亚洲人成电影免费在线| 白带黄色成豆腐渣| 搡老岳熟女国产| 亚洲五月天丁香| 无限看片的www在线观看| 亚洲精品色激情综合| 久久午夜综合久久蜜桃| 午夜福利欧美成人| 亚洲国产精品久久男人天堂| 色在线成人网| 久久久国产成人精品二区| 日韩欧美国产一区二区入口| 国产精品 国内视频| 91在线观看av| 他把我摸到了高潮在线观看| 精品久久久久久久毛片微露脸| 99在线视频只有这里精品首页| 欧美不卡视频在线免费观看 | 手机成人av网站| 久久久久久大精品| xxx96com| 国产私拍福利视频在线观看| 精品久久蜜臀av无| 精品人妻1区二区| 精品欧美国产一区二区三| 亚洲中文字幕日韩| 国产精品 国内视频| 少妇熟女aⅴ在线视频| 成人免费观看视频高清| 成人精品一区二区免费| 亚洲精品美女久久久久99蜜臀| 久久人妻av系列| 99热只有精品国产| 正在播放国产对白刺激| 国内揄拍国产精品人妻在线 | 黄色毛片三级朝国网站| 两个人免费观看高清视频| 国产亚洲欧美精品永久| 欧美日韩黄片免| 午夜亚洲福利在线播放| 欧美 亚洲 国产 日韩一| 亚洲无线在线观看| 脱女人内裤的视频| 国产精品二区激情视频| 高清在线国产一区| 少妇被粗大的猛进出69影院| 成人18禁在线播放| 男女做爰动态图高潮gif福利片| 成人亚洲精品av一区二区| 国产精品香港三级国产av潘金莲| 一区二区三区精品91| 免费在线观看视频国产中文字幕亚洲| 亚洲无线在线观看| 亚洲av中文字字幕乱码综合 | 一区二区三区高清视频在线| 男人操女人黄网站| 亚洲电影在线观看av| 精品欧美国产一区二区三| 老汉色∧v一级毛片| 国产真实乱freesex| 老汉色∧v一级毛片| 操出白浆在线播放| 国产精品久久久久久亚洲av鲁大| 十分钟在线观看高清视频www| 亚洲欧美日韩无卡精品| 久久草成人影院| 日韩精品免费视频一区二区三区| 午夜成年电影在线免费观看| 黄色女人牲交| 人人澡人人妻人| 亚洲成人国产一区在线观看| xxx96com| 国产黄a三级三级三级人| 国产亚洲精品av在线| 天天添夜夜摸| 欧美乱色亚洲激情| 亚洲欧美一区二区三区黑人| 亚洲国产精品合色在线| 一边摸一边抽搐一进一小说| 侵犯人妻中文字幕一二三四区| 99riav亚洲国产免费| 午夜两性在线视频| 麻豆国产av国片精品| 日韩精品青青久久久久久| 很黄的视频免费| 美女免费视频网站| 成人18禁在线播放| 宅男免费午夜| 精品国产乱子伦一区二区三区| 欧美日韩中文字幕国产精品一区二区三区| 欧美乱码精品一区二区三区| 中文字幕另类日韩欧美亚洲嫩草| 亚洲国产精品合色在线| 欧美成人一区二区免费高清观看 | 老鸭窝网址在线观看| 老司机深夜福利视频在线观看| 午夜福利免费观看在线| 久久国产亚洲av麻豆专区| 女性生殖器流出的白浆| 欧美最黄视频在线播放免费| 男女床上黄色一级片免费看| 欧美一级a爱片免费观看看 | 天天躁狠狠躁夜夜躁狠狠躁| aaaaa片日本免费| 亚洲美女黄片视频| 国产成人啪精品午夜网站| 国产成年人精品一区二区| 在线国产一区二区在线| 搡老妇女老女人老熟妇| 曰老女人黄片| 久久久久九九精品影院| 午夜视频精品福利| 国产97色在线日韩免费| 中文字幕人成人乱码亚洲影| 国产精品永久免费网站| 久久中文看片网| 日日爽夜夜爽网站| 国产精品av久久久久免费| 老司机在亚洲福利影院|