Discrete Element Simulation of Orthogonal Machining of Soda-Lime Glass

Transparent, brittle materials, like glass, are used in various applications due to their advantages of mechanical and optical properties. However, their hard and brittle nature causes significant challenges to researchers when they design and test a new machining process. In order to optimize this time-consuming process, discrete element method (DEM) is applied to simulate the cutting process of soda-lime glass in this study. The first step is to create a synthetic material that behaves like soda-lime glass. Then, the macroproperties are calibrated by adjusting the microparameters of the DEM model to match the mechanical properties of the real soda-lime glass. Finally, orthogonal machining simulations are conducted and model validation are conducted by comparing the predicted cutting forces with those from the orthogonal cutting experiments.


Introduction
Glass materials are widely used in aerospace, art, cookware and pharmaceutics due to their superior properties such as high strength, hardness, wear resistance, large optical transmission range, and biocompatibility [1].In modern architecture and civil engineering applications, glass composite, such as laminated glass and glass-fiber-reinforced-polymer, are playing important roles [2,3].However, due to their high hardness and extreme brittleness, these materials are generally difficult to machine [4].Glass is classified into several types depending on chemical composition.Soda lime glass is one of these types which is prepared by melting and mixing of silica with alumina, dolomite, lime, and soda at a temperature of up to 1675 • C. The soda lime glass material is used for a wide range of products, such as mirrors, photo masks, data storage disks, microscopic slides, touch screen filters, printed circuit substrates, photographic plates, optical windows, and microfluidic devices [5].There are several machining processes used to machine soda lime glass including electrochemical discharge machining, precision milling, microwave drilling, precision lathe machining, abrasive jet machining, micro-grinding, and laser machining.Each of them has its unique advantages and disadvantages.Electrochemical discharge machining can achieve high geometrical accuracy and a high machining rate, but micro-cracks and thermal damage are noted on machined surfaces [6].In a precision machining process, the feed rate significantly affects the cutting forces and surface finish.It is also observed that the higher critical uncut chip thickness at higher cutting speed affect the selection of the feed rate [7].Drilling of micro-holes can be achieved with microwave drilling.However, high heat is generated due to microwave energy.This caused cracking and deformation defects due to melting in the drilling zone [8].One example of precision lathe machining is to use an ultra-precision lathe with a single crystal diamond tool.With this machining method, surface roughness is highly dependent on tool wear, which is caused by diffusion, mechanical friction, thermo-chemical action, and abrasive wear [9].Due to a lack of heat affected zones, the depth of damage to the surface by abrasive jet machining is very little.However, the material removal rate is too low.For example, the typical material removal rate for cutting glass using abrasive jet machining is 16.4 mm 3 /min [10].Although a critical depth of 2 to 5 nm and surface roughness of 78 to 980 nm are achieved by a micro-grinding process on soda lime glass material [11,12], the low material removal rate limits the possibility of adopting micro grinding in mass production.Therefore, a hybrid machining methodology, which combines laser micromachining and mechanical machining, is introduced to reduce the cutting force, subsurface damage, and tool wear [1].Hybrid machining, which combines different machining actions or phases, is able to improve the machining performance.It can take the combined advantages and avoid or reduce some adverse effects the constituent processes produce when they are individually applied [10].However, due to the high brittleness of glass, the process of random crack initiation and propagation is difficult to study through theoretical analysis and experimental observation.Hence, numerical methods have been used more and more to simulate glass behavior under dynamic loading in recent years [2].
The finite element method (FEM) and extended finite element method (XFEM) have been widely used for numerical simulation of glass cracking and progressive damage in structural and civil engineering applications.For example, in order to optimize the pedestrian protection design for laminated windshields and reduce the severity of head injuries resulting from pedestrian-vehicle crashes, the headform impact response for the laminated windshield is simulated by FEM [13].With the help of XFEM simulation, researchers are also able to perform a mechanical delamination analysis of the three-dimensional laminated glass model [14].In addition, FEM has been used in simulating cutting processes in the past few decades, especially in terms of different material properties and cutting conditions on material separation, cutting force, chip formation and heat generation [15].However, when element distortion becomes large or nodes separate, the failure grids will be deleted and taken out of the cutting zone.It also cannot express the cracks randomly and intuitively because the separation line is usually predefined.Thus, simulating the initiation of cracks in the cutting process is not suitable for FEM.When cracks initiate during machining, the material model for the workpiece would be essentially changed from continuous to discontinuous medium.Discrete element method (DEM) which was introduced by Cundall perfectly solved this type of problems [16].Due to its meshless feature, large bulks can be separated into small units and it is convenient to deal with large deformation and fracture in cutting processes.Moreover, since the DEM model can describe non-linear behavior, it has great advantages to model the microstructures of brittle materials and handle the complex particle contact physical process with coupled shear and bulk deformation effects.In recent years, DEM has been widely used for simulating cutting processes of brittle materials, such as rock, ceramics, and carbon fiber-reinforced polymer [17][18][19].In addition to these approaches, a hybrid numerical method is used lately to model the fracture of brittle materials with specific reference to glass [20], in which FEM and DEM are combined satisfactorily in the simulation of a glass beam subjected to a low-velocity hard body impact.
In this paper, the orthogonal cutting process of soda-lime glass is simulated with a DEM model.The purpose is to find out how the DEM model can accurately simulate the cutting process at different cutting conditions.The first step is to create a synthetic material that behaves like soda-lime glass.Then the macro-properties are calibrated by adjusting the micro-parameters of the DEM model to match the mechanical properties of the real soda-lime glass.Finally, orthogonal machining simulations are conducted and model validation are conducted by comparing the predicted cutting forces with those from the orthogonal cutting experiments.Through this study, the DEM simulation of the orthogonal cutting process of soda-lime glass is identified and validated.

DEM Model of Soda-Lime Glass
Particle flow code in two dimensions (PFC2D) [21] is used to create a synthetic material that behaves like soda-lime glass.PFC2D is a commercial DEM simulation software which is widely used for different engineering applications.PFC2D models the interactions and movements of circular particles by distinct element method, which allows finite displacements and rotations of discrete bodies, including complete detachment, and recognizes new contacts automatically as the calculation progresses.PFC2D can model a solid by bonding every particle to its neighbors and compacting the particles within four walls.The resulting assembly is able to not only possess elastic properties, but also maintains plastic properties which can generate fractures when bonds between these particles break.Although the standard bonded-particle model (BPM) has been successfully applied to solve a wide range of problems, there are still some intrinsic problems that cannot be avoided using the standard BPM, such as an unrealistically low unconfined compressive strength to tensile strength (UCS/UTS) ratio, an excessively low internal friction angle, and a linear envelope [22].In order to avoid those problems in the simulating process, it is decided to apply the flat-joint model to model soda-lime glass.

Fundamentals of Flat-Joint Model
Flat-joint is an interface of ball-ball or ball-face contact (Figure 1) which can be deformable, partially damaged, and breakable.The breakage of each bonded element contributes partial damage to the interface and is denoted as a crack.The flat-joint model provides the macroscopic behavior of a finite-size, linear elastic, and either bonded or frictional interfaces that may sustain partial damage.A bonded element has a linear elastic behavior until the strength limit is exceeded and bond breaks, making the element unbonded.An unbonded element has linear elastic and frictional behavior with slip accommodated by imposing a Coulomb limit on the shear force.
progresses.PFC2D can model a solid by bonding every particle to its neighbors and compacting the particles within four walls.The resulting assembly is able to not only possess elastic properties, but also maintains plastic properties which can generate fractures when bonds between these particles break.Although the standard bonded-particle model (BPM) has been successfully applied to solve a wide range of problems, there are still some intrinsic problems that cannot be avoided using the standard BPM, such as an unrealistically low unconfined compressive strength to tensile strength (UCS/UTS) ratio, an excessively low internal friction angle, and a linear envelope [22].In order to avoid those problems in the simulating process, it is decided to apply the flat-joint model to model soda-lime glass.

Fundamentals of Flat-Joint Model
Flat-joint is an interface of ball-ball or ball-face contact (Figure 1) which can be deformable, partially damaged, and breakable.The breakage of each bonded element contributes partial damage to the interface and is denoted as a crack.The flat-joint model provides the macroscopic behavior of a finite-size, linear elastic, and either bonded or frictional interfaces that may sustain partial damage.A bonded element has a linear elastic behavior until the strength limit is exceeded and bond breaks, making the element unbonded.An unbonded element has linear elastic and frictional behavior with slip accommodated by imposing a Coulomb limit on the shear force.After a flat-joint contact is installed, the force and moment in each element are initialized to zero.Then they will be updates by generalized force-displacement law which demonstrate internal force to the relative motion at the contact.Each element carries a force and moment (F (e) and M (e) ) that act at the element centroid (Figure 2), and the formulation is expressed in terms of a set of internal rheological parameters (Figure 3).After a flat-joint contact is installed, the force and moment in each element are initialized to zero.Then they will be updates by generalized force-displacement law which demonstrate internal force to the relative motion at the contact.Each element carries a force and moment (F (e) and M (e) ) that act at the element centroid (Figure 2), and the formulation is expressed in terms of a set of internal rheological parameters (Figure 3).The element forces and moments produce a statically equivalent force and moment at the center of the interface by: where r (e) = x (e) − xc is the relative position, x (e) is the centroid location of element (e) and xc is the center of the interface.The element normal stress and shear stress are given by: where Fn and Fs are the normal and shear forces acting at the element, A is the area of the element.There are two main types of elements for a flat-joint contact, which are bonded elements and unbonded elements.The bonded element's strength envelope is as follows (Figure 4).
is the tensile strength and the elements will break with a tensile crack when > .The shear strength is which follows the Coulomb criterion given by: where is the bond cohesion, is the normal stress, and is the local friction angle.The shear stress will lock to when < , or the bond will break in shear and the residual friction takes effects:  The element forces and moments produce a statically equivalent force and moment at the center of the interface by: where r (e) = x (e) − xc is the relative position, x (e) is the centroid location of element (e) and xc is the center of the interface.The element normal stress and shear stress are given by: where Fn and Fs are the normal and shear forces acting at the element, A is the area of the element.There are two main types of elements for a flat-joint contact, which are bonded elements and unbonded elements.The bonded element's strength envelope is as follows (Figure 4).
is the tensile strength and the elements will break with a tensile crack when > .The shear strength is which follows the Coulomb criterion given by: where is the bond cohesion, is the normal stress, and is the local friction angle.The shear stress will lock to when < , or the bond will break in shear and the residual friction takes effects: The element forces and moments produce a statically equivalent force and moment at the center of the interface by: where r (e) = x (e) − x c is the relative position, x (e) is the centroid location of element (e) and x c is the center of the interface.The element normal stress and shear stress are given by: where F n and F s are the normal and shear forces acting at the element, A is the area of the element.There are two main types of elements for a flat-joint contact, which are bonded elements and unbonded elements.The bonded element's strength envelope is as follows (Figure 4).σ b is the tensile strength and the elements will break with a tensile crack when σ max > σ b .The shear strength is τ c which follows the Coulomb criterion given by: where c b is the bond cohesion, σ is the normal stress, and φ b is the local friction angle.The shear stress will lock to τ max when τ max < τ c , or the bond will break in shear and the residual friction takes effects: where τ r is the residual friction strength and φ r = tan −1 (µ) is the residual friction angle.This type of element can be broken with a tensile or shear crack.Additionally, the unbonded elements follow such a force-displacement law.The normal stress is given by: where k n is the normal stiffness and g is the interface gap.The shear strength follows the Coulomb sliding criterion: where φ r is the residual friction angle.When the shear stress τ max < τ r , it remains as τ max , or the shear stress will change and enforces sliding at the interface.Within this type of element, there is no tensile or shear crack.
where is the residual friction strength and = ( ) is the residual friction angle.This type of element can be broken with a tensile or shear crack.Additionally, the unbonded elements follow such a force-displacement law.The normal stress is given by: = 0, unbonded and ̅ ≥ 0 − ̅ , o t h e r w i s e where kn is the normal stiffness and is the interface gap.The shear strength follows the Coulomb sliding criterion: where is the residual friction angle.When the shear stress < , it remains as , or the shear stress will change and enforces sliding at the interface.Within this type of element, there is no tensile or shear crack.

Creation of the Bonded Particle Specimen
The first, and most crucial, step of a successful numerical simulation is to create the specimen that behaves like the real material.There are three general steps to generate the specimen in PFC2D: (1) creating a dense particle assembly; (2) installing flat-joint bonds at contacts; and (3) removing confinement for the specimen.The specimen produced by this procedure is a well-connected particle assembly with low locked-in forces.There are two sets of parameters which need to be defined when creating a synthetic specimen (Table 1).One includes the geometric and physical parameters, which are the specimen height, specimen width, particle size, and particle density.The other contains constitutive parameters to define the micro properties of the particles and their contact and bond, which contain the particle contact stiffness, particle stiffness ratio, particle friction coefficient, particle damping coefficient, bond shear strength, bond normal strength, and friction angle.
A synthetic specimen is created following the procedure and considerations mentioned above, as shown in Figure 5.The micro-parameters needed for the specimen generation are estimated values at this point and will be adjusted to match the actual material properties in the next section.The specimen is 2 mm in height and 1 mm in width.It contains 7406 circular disks of unit thickness and the radii of the disks range from 15 to 20 μm.The particle density is set to be 2400 kg m −3 according

Creation of the Bonded Particle Specimen
The first, and most crucial, step of a successful numerical simulation is to create the specimen that behaves like the real material.There are three general steps to generate the specimen in PFC2D: (1) creating a dense particle assembly; (2) installing flat-joint bonds at contacts; and (3) removing confinement for the specimen.The specimen produced by this procedure is a well-connected particle assembly with low locked-in forces.There are two sets of parameters which need to be defined when creating a synthetic specimen (Table 1).One includes the geometric and physical parameters, which are the specimen height, specimen width, particle size, and particle density.The other contains constitutive parameters to define the micro properties of the particles and their contact and bond, which contain the particle contact stiffness, particle stiffness ratio, particle friction coefficient, particle damping coefficient, bond shear strength, bond normal strength, and friction angle.
A synthetic specimen is created following the procedure and considerations mentioned above, as shown in Figure 5.The micro-parameters needed for the specimen generation are estimated values at this point and will be adjusted to match the actual material properties in the next section.The specimen is 2 mm in height and 1 mm in width.It contains 7406 circular disks of unit thickness and the radii of the disks range from 15 to 20 µm.The particle density is set to be 2400 kg m −3 according to soda-lime glass properties [23].The network shown in Figure 6 indicates the flat-joint bonds, which are the straight lines connecting the particles together.
to soda-lime glass properties [23].The network shown in Figure 6 indicates the flat-joint bonds, which are the straight lines connecting the particles together.to soda-lime glass properties [23].The network shown in Figure 6 indicates the flat-joint bonds, which are the straight lines connecting the particles together.Table 1.Macro-and micro-parameters of the synthetic specimen.Table 1.Macro-and micro-parameters of the synthetic specimen.

Calibration of the DEM Model (Mechanical Properties)
It is known that the macro-properties of the DEM model are directly related to the micro-parameters of the flat-joint contact model.However, there is no definite relationship between macro-properties of the sample and the micro-parameters of the flat-joint model.Thus, a series of tests need to be conducted to calibrate the model and make sure that the simulation model has the equivalent properties with those of soda-lime glass.The macro-properties to be matched are elastic modulus (E), Poisson's ratio (ν), ultimate tensile strength (σ t ), and compressive strength (σ c ).

Uniaxial Tensile Test
The same specimen created in the previous section is used for the uniaxial tensile test.A row of particles at the top and bottom of the sample are identified as "grips".As shown in Figure 7a, the top grip is moved upward and the bottom grip is moved downward with the same speed of 0.02 mm/min to apply tensile loading.The model is loaded until the axial stress drops to 70% of the peak stress.The fracturing is tracked by monitoring the broken bonds and the cracks are plotted during the test.The cracks can be generated by both shear and tensile force, which are also shown as green and red lines in Figure 7a.Additionally, there are two important quantities being recorded during the test.One is the stress, which is the average stress in the measurement region of the sample [24].The other is the strain, which is the ratio of total deformation to the initial height of the sample.Since PFC2D is able to track the positions of each particle during the simulation, the displacement of the particles can be calculated and retrieved by manipulating the FISH functions which is an embedded programming language of PFC2D.As shown in Figure 7b, the stress-strain curve can be plotted.Elastic modulus and tensile strength can be identified through the plot.

Calibration of the DEM Model (Mechanical Properties)
It is known that the macro-properties of the DEM model are directly related to the microparameters of the flat-joint contact model.However, there is no definite relationship between macroproperties of the sample and the micro-parameters of the flat-joint model.Thus, a series of tests need to be conducted to calibrate the model and make sure that the simulation model has the equivalent properties with those of soda-lime glass.The macro-properties to be matched are elastic modulus (E), Poisson's ratio (ν), ultimate tensile strength (σt), and compressive strength (σc).

Uniaxial Tensile Test
The same specimen created in the previous section is used for the uniaxial tensile test.A row of particles at the top and bottom of the sample are identified as "grips".As shown in Figure 7a, the top grip is moved upward and the bottom grip is moved downward with the same speed of 0.02 mm/min to apply tensile loading.The model is loaded until the axial stress drops to 70% of the peak stress.The fracturing is tracked by monitoring the broken bonds and the cracks are plotted during the test.The cracks can be generated by both shear and tensile force, which are also shown as green and red lines in Figure 7a.Additionally, there are two important quantities being recorded during the test.One is the stress, which is the average stress in the measurement region of the sample [24].The other is the strain, which is the ratio of total deformation to the initial height of the sample.Since PFC2D is able to track the positions of each particle during the simulation, the displacement of the particles can be calculated and retrieved by manipulating the FISH functions which is an embedded programming language of PFC2D.As shown in Figure 7b, the stress-strain curve can be plotted.Elastic modulus and tensile strength can be identified through the plot.

Uniaxial Compressive Test
The same specimen is used for the unconfined compressive test.Instead of "grips", two rigid walls are created at the top and bottom surface of the sample as shown in Figure 8a.These two walls are considered as rigid supports, which cannot be deformed or broken during the process.During the compressive test, the top wall is moved downward and bottom wall moved upward with the same speed of 0.02 mm/min to apply compressive loading.The model is loaded until the axial stress

Uniaxial Compressive Test
The same specimen is used for the unconfined compressive test.Instead of "grips", two rigid walls are created at the top and bottom surface of the sample as shown in Figure 8a.These two walls are considered as rigid supports, which cannot be deformed or broken during the process.During the compressive test, the top wall is moved downward and bottom wall moved upward with the same speed of 0.02 mm/min to apply compressive loading.The model is loaded until the axial stress drops to 70% of the peak stress.With the same method as the uniaxial tensile test, the stress-strain curve can be plotted, as shown in Figure 8b.
drops to 70% of the peak stress.With the same method as the uniaxial tensile test, the stress-strain curve can be plotted, as shown in Figure 8b.Since the strain in transverse direction can be calculated as well, we can determine the Poisson's ratio, which is the ratio of transverse strain to axial strain.The normal and most straightforward method for measuring Poisson's ratio is conducting a normal mechanical tensile test [25].Thus, the Poisson's ratio is calibrated through the uniaxial tensile test.Through a series of adjustments, the DEM model is calibrated to match the soda-lime glass properties as listed in Table 2.The next step is to simulate the cutting processes of the specimen.

Cutting Simulation
The orthogonal cutting with different cutting conditions is simulated and each simulation is repeated three times using specimens generated with different random particle arrangements to reduce the influence of particle arrangement.The cutting conditions are shown in Table 3.The model geometry and boundary conditions for the cutting simulation are shown in Figure 9.The glass workpiece is 1 mm in length and 1 mm in height, and it contains 2085 particles.The particles on the left side and bottom side of the workpiece are fixed.The cutting tool is modeled as a rigid body using wall elements which has the same rake angle (−15°) and same clearance angle (15°) during the cutting simulation.The depths of cut are 0.1 mm and 0.2 mm with the fixed cutting speed of 4.23 mm/s.
The glass specimens have two different surface conditions.One is laser-treated with grooves on the top surface induced by femtosecond laser micromachining [1].The groove depths are 0.1 mm and 0.2 mm according to the depth of cut.The groove distance and width are 0.2 mm and 0.06 mm (Figure 10).The other type of specimen is the untreated sample which does not have microgrooves.Since the strain in transverse direction can be calculated as well, we can determine the Poisson's ratio, which is the ratio of transverse strain to axial strain.The normal and most straightforward method for measuring Poisson's ratio is conducting a normal mechanical tensile test [25].Thus, the Poisson's ratio is calibrated through the uniaxial tensile test.Through a series of adjustments, the DEM model is calibrated to match the soda-lime glass properties as listed in Table 2.The next step is to simulate the cutting processes of the specimen.

Cutting Simulation
The orthogonal cutting with different cutting conditions is simulated and each simulation is repeated three times using specimens generated with different random particle arrangements to reduce the influence of particle arrangement.The cutting conditions are shown in Table 3.The model geometry and boundary conditions for the cutting simulation are shown in Figure 9.The glass workpiece is 1 mm in length and 1 mm in height, and it contains 2085 particles.The particles on the left side and bottom side of the workpiece are fixed.The cutting tool is modeled as a rigid body using wall elements which has the same rake angle (−15 • ) and same clearance angle (15 • ) during the cutting simulation.The depths of cut are 0.1 mm and 0.2 mm with the fixed cutting speed of 4.23 mm/s.
The glass specimens have two different surface conditions.One is laser-treated with grooves on the top surface induced by femtosecond laser micromachining [1].The groove depths are 0.1 mm and 0.2 mm according to the depth of cut.The groove distance and width are 0.2 mm and 0.06 mm (Figure 10).The other type of specimen is the untreated sample which does not have microgrooves.Thus, there are a total of four different conditions to simulate for the cutting process: (1) 0.1 mm depth of cut without microgrooves; (2) 0.1 mm depth of cut with microgrooves; (3) 0.2 mm depth of cut without microgrooves; and (4) 0.2 mm depth of cut with microgrooves.The main (horizontal direction) and thrust (vertical direction) cutting forces are recorded during the simulation.Figures 11-14 show the results for each cutting condition: (1) A 0.1 mm depth of cut without microgrooves:   Thus, there are a total of four different conditions to simulate for the cutting process: (1) 0.1 mm depth of cut without microgrooves; (2) 0.1 mm depth of cut with microgrooves; (3) 0.2 mm depth of cut without microgrooves; and (4) 0.2 mm depth of cut with microgrooves.The main (horizontal direction) and thrust (vertical direction) cutting forces are recorded during the simulation.Figures 11-14 show the results for each cutting condition: (1) A 0.1 mm depth of cut without microgrooves:  Thus, there are a total of four different conditions to simulate for the cutting process: (1) 0.1 mm depth of cut without microgrooves; (2) 0.1 mm depth of cut with microgrooves; (3) 0.2 mm depth of cut without microgrooves; and (4) 0.2 mm depth of cut with microgrooves.The main (horizontal direction) and thrust (vertical direction) cutting forces are recorded during the simulation.Figures 11-14 show the results for each cutting condition: (1) A 0.1 mm depth of cut without microgrooves: (2) A 0.1 mm depth of cut with microgrooves: (3) A 0.2 mm depth of cut without microgrooves:   (2) A 0.1 mm depth of cut with microgrooves: (3) A 0.2 mm depth of cut without microgrooves:  (2) A 0.1 mm depth of cut with microgrooves: (3) A 0.2 mm depth of cut without microgrooves:   In all the chip formation images, the green, short lines denote broken bonds between the particles which are caused by shear failure; the blue, short lines are also broken bonds, but they are caused by tensile failure.It is obvious that chip formation is due to initiation and propagation of brittle cracks in the material and varied size of chip segments are the evidence for the brittle material removal process.As the depth of cut increases from 0.1 mm (e.g.Fig. 11) to 0.2 mm (e.g.Fig. 13), chips become larger and subsurface cracks become deeper in general.Additionally, there are distinct differences between the regular and grooved workpieces under the same cutting conditions.First, the broken bonds are less for the grooved workpiece, which is particularly true for the 0.2 mm depth of cut when comparing Figures 13 and 14.Second, the cracks left on the machined part are shorter and shallower for the grooved workpiece.The main benefit of installing microgrooves on the workpiece is reducing the subsurface cracks dramatically.The workpiece is very strong without the microgrooves and needs a large force to initiate a crack.With microgrooves on the surface, the workpiece is significantly weakened.Since cracks initiate mainly at the bottom of the microgrooves, which is the weakest part of the workpiece, less force is needed to initiate a crack and subsurface cracks are almost nonexistent.
The cutting forces in both the horizontal and vertical directions are recorded during the simulation process.It can be seen that numerous force peaks exist for both the regular and grooved workpieces, which are caused by the initiation and propagation of cracks due to bond breakage between particles.This behavior is typical for brittle material removal processes, characterized by random peaks and valleys which correspond to force build-up followed by sudden fracture occurrence.It is found that the force magnitude for the grooved workpiece is much smaller than for the regular one because the grooves reduce the material strength.It is also reasonable that when the depth of cut is reduced, the cutting force is reduced correspondingly.The average forces are listed in Table 4 for comparison.In all the chip formation images, the green, short lines denote broken bonds between the particles which are caused by shear failure; the blue, short lines are also broken bonds, but they are caused by tensile failure.It is obvious that chip formation is due to initiation and propagation of brittle cracks in the material and varied size of chip segments are the evidence for the brittle material removal process.As the depth of cut increases from 0.1 mm (e.g., Figure 11) to 0.2 mm (e.g., Figure 13), chips become larger and subsurface cracks become deeper in general.Additionally, there are distinct differences between the regular and grooved workpieces under the same cutting conditions.First, the broken bonds are less for the grooved workpiece, which is particularly true for the 0.2 mm depth of cut when comparing Figures 13 and 14.Second, the cracks left on the machined part are shorter and shallower for the grooved workpiece.The main benefit of installing microgrooves on the workpiece is reducing the subsurface cracks dramatically.The workpiece is very strong without the microgrooves and needs a large force to initiate a crack.With microgrooves on the surface, the workpiece is significantly weakened.Since cracks initiate mainly at the bottom of the microgrooves, which is the weakest part of the workpiece, less force is needed to initiate a crack and subsurface cracks are almost nonexistent.
The cutting forces in both the horizontal and vertical directions are recorded during the simulation process.It can be seen that numerous force peaks exist for both the regular and grooved workpieces, which are caused by the initiation and propagation of cracks due to bond breakage between particles.This behavior is typical for brittle material removal processes, characterized by random peaks and valleys which correspond to force build-up followed by sudden fracture occurrence.It is found that the force magnitude for the grooved workpiece is much smaller than for the regular one because the grooves reduce the material strength.It is also reasonable that when the depth of cut is reduced, the cutting force is reduced correspondingly.The average forces are listed in Table 4 for comparison.

Model Validation
In order to validate the simulation model, the orthogonal machining experimental setup is constructed using a vertical Bridgeport milling machine with a cutting tool fixed on the locked vertical spindle column.The soda-lime glass is fixed on the horizontal carriage (Figure 15).The fixture is mounted onto a Kistler three-component dynamometer.A Kistler dual-mode charge amplifier is used to amplify cutting force signals measured by the dynamometer.LabView (National Instruments, Austin, TX, USA) is used to control the computer data acquisition system.The workpiece, fixture, and dynamometer are clamped onto the movable carriage of the Bridgeport milling machine.The horizontal carriage feeds the workpiece on to the stationary tool holder which is mounted in the vertical tool chuck.A 15.875 × 15.875 × 6.35 mm ceramic insert is used as a cutting tool.A special tool holder, tilted to attain a negative 15-degree rake angle is used for the experiment.The cutting conditions are given in Table 5.
The cutting experiments contain four conditions which are exactly the same as in the simulation.For the laser treated samples, the microgrooves are made by laser scanning at the sample surface as shown in Figure 16.The main differences between the samples are the groove depths.For the 0.1 mm depth of cut, the groove depth is 0.1 mm.For the 0.2 mm depth of cut, the groove depth is 0.2 mm.

Model Validation
In order to validate the simulation model, the orthogonal machining experimental setup is constructed using a vertical Bridgeport milling machine with a cutting tool fixed on the locked vertical spindle column.The soda-lime glass is fixed on the horizontal carriage (Figure 15).The fixture is mounted onto a Kistler three-component dynamometer.A Kistler dual-mode charge amplifier is used to amplify cutting force signals measured by the dynamometer.LabView (National Instruments, Austin, TX, USA) is used to control the computer data acquisition system.The workpiece, fixture, and dynamometer are clamped onto the movable carriage of the Bridgeport milling machine.The horizontal carriage feeds the workpiece on to the stationary tool holder which is mounted in the vertical tool chuck.A 15.875 × 15.875 × 6.35 mm ceramic insert is used as a cutting tool.A special tool holder, tilted to attain a negative 15-degree rake angle is used for the experiment.The cutting conditions are given in Table 5.
The cutting experiments contain four conditions which are exactly the same as in the simulation.For the laser treated samples, the microgrooves are made by laser scanning at the sample surface as shown in Figure 16.The main differences between the treated samples are the groove depths.For the 0.1 mm depth of cut, the groove depth is 0.1 mm.For the 0.2 mm depth of cut, the groove depth is 0.2 mm.

Model Validation
In order to validate the simulation model, the orthogonal machining experimental setup is constructed using a vertical Bridgeport milling machine with a cutting tool fixed on the locked vertical spindle column.The soda-lime glass is fixed on the horizontal carriage (Figure 15).The fixture is mounted onto a Kistler three-component dynamometer.A Kistler dual-mode charge amplifier is used to amplify cutting force signals measured by the dynamometer.LabView (National Instruments, Austin, TX, USA) is used to control the computer data acquisition system.The workpiece, fixture, and dynamometer are clamped onto the movable carriage of the Bridgeport milling machine.The horizontal carriage feeds the workpiece on to the stationary tool holder which is mounted in the vertical tool chuck.A 15.875 × 15.875 × 6.35 mm ceramic insert is used as a cutting tool.A special tool holder, tilted to attain a negative 15-degree rake angle is used for the experiment.The cutting conditions are given in Table 5.
The cutting experiments contain four conditions which are exactly the same as in the simulation.For the laser treated samples, the microgrooves are made by laser scanning at the sample surface as shown in Figure 16.The main differences between the treated samples are the groove depths.For the 0.1 mm depth of cut, the groove depth is 0.1 mm.For the 0.2 mm depth of cut, the groove depth is 0.2 mm.The dynamometer is calibrated before conducting the experiments.Force data is collected at a sampling rate of 200 data points per second.Table 6 shows the average cutting force recorded during each experiment.In order to compare the simulation and experimental results for each condition, the average forces for the three replications are taken and plotted in Figures 17 and 18.In general, the experimental forces are 7-26% more than the simulation results.For most of the cases, the difference between the simulation and experimental results are below 20%.Especially for the regular workpiece, the results are more reliable and repeatable because the macro-properties of the EDM model are close to the real condition.However, for the grooved workpiece, the results are not as good as for the regular workpieces.For the 0.2 mm depth of cut of grooved workpieces, the thrust force from the experiment is about 60% larger than the simulation result.A possible reason is that the groove geometry of the model is different from that of the real workpieces.The grooves of the laser-treated sample have an inverted trapezoidal shape, which means the bottom wall of the groove is thicker and stronger.However, for the DEM model, the grooves are made by deleting the particles within certain areas, so that the groove walls are more like rectangular shape, which are weaker at the bottom compared with the trapezoidal-shaped walls for the real workpieces.This effect will become more obvious for the 0.2 mm depth of cut since the wall becomes taller and less force is needed to break the bonds.This argument will be carefully examined in future simulation work.The dynamometer is calibrated before conducting the experiments.Force data is collected at a sampling rate of 200 data points per second.Table 6 shows the average cutting force recorded during each experiment.In order to compare the simulation and experimental results for each condition, the average forces for the three replications are taken and plotted in Figures 17 and 18.In general, the experimental forces are 7-26% more than the simulation results.For most of the cases, the difference between the simulation and experimental results are below 20%.Especially for the regular workpiece, the results are more reliable and repeatable because the macro-properties of the EDM model are close to the real condition.However, for the grooved workpiece, the results are not as good as for the regular workpieces.For the 0.2 mm depth of cut of grooved workpieces, the thrust force from the experiment is about 60% larger than the simulation result.A possible reason is that the groove geometry of the model is different from that of the real workpieces.The grooves of the laser-treated sample have an inverted trapezoidal shape, which means the bottom wall of the groove is thicker and stronger.However, for the DEM model, the grooves are made by deleting the particles within certain areas, so that the groove walls are more like rectangular shape, which are weaker at the bottom compared with the trapezoidal-shaped walls for the real workpieces.This effect will become more obvious for the 0.2 mm depth of cut since the wall becomes taller and less force is needed to break the bonds.This argument will be carefully examined in future simulation work.

Conclusions
In this paper, the orthogonal cutting of soda-lime glass is simulated with PFC2D base on discrete element method.The results of cutting force, subsurface cracks and chips are studied, and the cutting forces are validated through the orthogonal cutting experiment.The following conclusions are obtained from this study.

•
The cause of crack initiation can be identified by the simulation model.There are two types of cracks: one is caused by shear force and the other is caused by tensile force.

•
The simulations show that cutting forces can be reduced and subsurface cracks can be eliminated when cutting grooved workpieces.

•
The geometries of the simulation model highly affect the result.Compared to the grooved specimen, the simulation result for the regular one is more accurate and repeatable.

•
Simulation time is highly dependent on the number of particles for the model and computing capacity.Smaller particle size and larger number of particles can enhance the details of the model.More powerful computing methods are recommended, such as parallel computing and supercomputing to increase the simulation efficiency.

Conclusions
In this paper, the orthogonal cutting of soda-lime glass is simulated with PFC2D base on discrete element method.The results of cutting force, subsurface cracks and chips are studied, and the cutting forces are validated through the orthogonal cutting experiment.The following conclusions are obtained from this study.

•
The cause of crack initiation can be identified by the simulation model.There are two types of cracks: one is caused by shear force and the other is caused by tensile force.

•
The simulations show that cutting forces can be reduced and subsurface cracks can be eliminated when cutting grooved workpieces.

•
The geometries of the simulation model highly affect the result.Compared to the grooved specimen, the simulation result for the regular one is more accurate and repeatable.

•
Simulation time is highly dependent on the number of particles for the model and computing capacity.Smaller particle size and larger number of particles can enhance the details of the model.More powerful computing methods are recommended, such as parallel computing and supercomputing to increase the simulation efficiency.

Figure 2 .
Figure 2.Force and moment on an element.Adapted from[21]

Figure 7 .
Figure 7. Uniaxial tensile test.(a) PFC2D model of uniaxial tensile test; and (b) the stress-strain curve of the uniaxial tensile test.

Figure 7 .
Figure 7. Uniaxial tensile test.(a) PFC2D model of uniaxial tensile test; and (b) the stress-strain curve of the uniaxial tensile test.

Figure 8 .
Figure 8. Unconfined compressive test.(a) PFC2D model of unconfined compressive test; and (b) the stress-strain curve of unconfined compressive test.

Figure 8 .
Figure 8. Unconfined compressive test.(a) PFC2D model of unconfined compressive test; and (b) the stress-strain curve of unconfined compressive test.

J 15 Figure 9 .
Figure 9. DEM model of the cutting simulation.

Figure 10 .
Figure 10.DEM model of the laser-treated sample.

Figure 10 .
Figure 10.DEM model of the laser-treated sample.

Figure 10 .
Figure 10.DEM model of the laser-treated sample.

Figure 11 .
Figure 11.A 0.1 mm depth of cut without microgrooves and recorded forces.

Figure 12 .
Figure 12.A 0.1 mm depth of cut with microgrooves and recorded forces.

Figure 13 .
Figure 13.A 0.2 mm depth of cut without microgrooves and recorded forces.

Figure 11 .
Figure 11.A 0.1 mm depth of cut without microgrooves and recorded forces.

Figure 12 .
Figure 12.A 0.1 mm depth of cut with microgrooves and recorded forces.

Figure 13 .
Figure 13.A 0.2 mm depth of cut without microgrooves and recorded forces.

Figure 12 . 15 Figure 11 .
Figure 12.A 0.1 mm depth of cut with microgrooves and recorded forces.

Figure 12 .
Figure 12.A 0.1 mm depth of cut with microgrooves and recorded forces.

Figure 13 .
Figure 13.A 0.2 mm depth of cut without microgrooves and recorded forces.

Figure 13 .
Figure 13.A 0.2 mm depth of cut without microgrooves and recorded forces.

Figure 14 .
Figure 14.A 0.2 mm depth of cut with microgrooves and recorded forces.

Figure 15 .
Figure 15.Experimental setup of the orthogonal machining test.

Figure 15 .
Figure 15.Experimental setup of the orthogonal machining test.

Figure 15 .
Figure 15.Experimental setup of the orthogonal machining test.

Figure 17 .
Figure 17.Cutting force comparison for the regular workpiece.

Figure 17 .
Figure 17.Cutting force comparison for the regular workpiece.

Figure 18 .
Figure 18.Cutting force comparison for the grooved workpiece.

Figure 18 .
Figure 18.Cutting force comparison for the grooved workpiece.

Table 1 .
Macro-and micro-parameters of the synthetic specimen.

Table 2 .
Comparison of properties between the DEM model and soda-lime glass.

Table 2 .
Comparison of properties between the DEM model and soda-lime glass.

Table 3 .
Cutting conditions for the simulation.

Table 3 .
Cutting conditions for the simulation.

Table 3 .
Cutting conditions for the simulation.

Table 4 .
Average cutting forces in horizontal and vertical directions.

Table 4 .
Average cutting forces in horizontal and vertical directions.

Table 5 .
Cutting conditions for the experiments.

Table 6 .
Average cutting forces for each experiment.

Table 5 .
Cutting conditions for the experiments.

Table 6 .
Average cutting forces for each experiment.