Real-time machine simulation using cutting force calculation based on a voxel material removal model

Machining processes must be adjusted regarding tolerances in dimension and shape to fulfill product requirements. For this purpose, machine simulations are used to allow a preliminary characterization of the given process, thus minimizing the number of physical prototypes and scrap parts. However, these simulations are either extremely specialized for single problems, e.g., dynamic machine behavior, or they are simplified to a kinematic simulation of the machine without considering the machine behavior at all. This article presents a new approach for a real-time machine simulation by combining four types of simulations to close this gap. This proposed approach uses a voxel-based material removal inside a kinematic machine simulation as input parameters for a cutting force calculation. Afterwards, the forces are applied to a multi-body simulation of the static machine behavior. Starting point of the simulation is a hardware-in-the-loop coupling of a real CNC and a real-time visualization of a virtual machine tool. The simulation is experimental verified by comparing the simulated cutting forces and displacements with the measured forces during the process and the resulting shape of the manufactured work piece. The presented conclusions show the general applicability of the proposed method for the simulation of milling processes.


Introduction
Today, there is an increasing demand of virtual prototypes not only for product development but also for production planning. A virtual prototype decreases the start-up time of a new product on a machine tool and allows simultaneous producing while planning the next process. Today, even with state-of-the-art process planning, the final commissioning steps for a new product are performed on the real machine tool to prevent erroneous production and unconsidered collisions. To further reduce the necessary tests on the machine tool, this paper presents a realtime machine simulation combining a multitude of different simulations of the machine tool. The described machine simulation includes a hardware-in-the-loop (HIL) simulation of a real NC unit, a voxel-based material removal, a cutting force calculation, and a multi-body simulation of the machine tool. The focus of this paper lies on the calculation of the cutting forces based on voxel material removal model for a possible real-time application.
Today's state-of-the-art CAD and CAM systems already include a rudimentary kinematic simulation of the machine tool which are still not sufficient for complex machines tools, especially regarding special purpose machines [1]. A more detailed commercial simulation (VERICUT 1 ) includes a number of machine variants, material removal of the work piece, and collision detection. Scientific research projects focus on improving the simulation range by including machine behavior and the calculation of cutting forces. Current research on machine simulation [2] concludes that it is essential to consider the subsystems "machine" and "process" to make a prediction of machine interaction. Looking even deeper into the topic of machine simulations, there are three relevant topics: a cutting force simulation, a material removal simulation, and virtual machine models as kinematic, static, thermal, or dynamic machine simulation models.

Material removal simulations
Material removal simulations have been researched for several decades. The possibility for predicting the work piece quality enables the evaluation of NC programs without manufacturing cost-intensive prototypes. Different approaches have been developed and evaluated, each with their own advantages and disadvantages.
One possibility is to simulate material removal using CAD functionality, e.g., via CAD kernels (ACIS 2 , Open Cascade 3 ). CAD models are based on volume models. Mostly, modern CAD software implements hybrid models based on B-rep [1] and CSG [2]. This kind of modeling relies on analytic mathematics and is very precise, but has the disadvantage of a high computational effort. Real-time applications cannot be based on such a solution.
Another approach is the application of a polygonal method, in particular triangle-based calculations. In computer graphics, it is common to use meshes based on triangles to describe 3D shapes. Although computational times are much shorter compared with the previously described approach, this method is still much too slow for real-time applications. Likewise, most implementations do not consider the nontrivial coplanar case for the intersection calculation, resulting in wrong meshes (e.g., with cracks).
A third kind of methods uses discrete representations to describe virtual objects and apply calculations to it. Two well-known variants are voxel and dexel methods. Voxel means "volume element." Here, a 3D object is described using cubes [3]. Their size determines the resolution of the model. The classic variant uses equally sized cubes. While calculations with cubes are trivial, a great number is necessary to achieve a high resolution. This leads to high memory consumption and high computational cost. Therefore, another approach was developed, called dexel [4], which stands for "depth pixel" or "depth element." It can be seen as a rectangular solid with an edge length for the base area and a length that is the depth of the dexel. Using the length attribute, a 3D solid can be described with much less memory consumption compared with the conventional voxel model. In addition, the computational cost decreases, which makes it attractive for simulation purposes.
Newer methods combine the voxel model with the octree data structure. This is known as sparse voxel octree (SVO) in computer graphics [5]. Voxels are no longer limited to one size. They are represented by the cubes (octants) of the octree. This leads to a lower count of voxels, therefore to lower memory consumption and to a faster computation. Because of these attributes, SVOs are able to work in real time which makes them a proper solution for the simulation approach of material removal in cutting.

Cutting force simulations
Cutting force simulations are a highly researched and still a very complex problem. The first research using an analytic model with empiric data was already carried out in the 1950s by Kronenberg [6] and Kienzle [7], using an exponential function describing the cutting force with increased cutting area. Altintas [8] describes a linear ratio between cutting forces and cutting area. These analytic models of the cutting process use a single cutting edge to determine the cutting conditions. Because of the continuous cut, the turning process can easily be simplified with this model. The evolving cutting area of the milling process must be considered for these simulations. In order to compensate for the change of the cutting area, the analytic models use either the average or the maximal cutting area. Additional adjustments are made for the number of cutting edges and the axial angle. Further research improved the models either for different kinds of end mills or with different material properties [9].
Ernst and Merchand [10] developed a mechanistic model using an approximation of the cutting forces based on the shear stress along the shear planes of the cutting edge. In recent years, FEM simulations have been developed for a more detailed view of the cutting process [11]. As the complexity and the computational effort for such simulations is very demanding, achieving the real-time goal with these methods is unlikely.

Machine simulations
In 2005, Altintas et al. [12] already showed the necessary simulations to obtain better understanding of the machine tool and the process as a whole. Denkena et al. [13] showed the next steps to achieve the merged machine simulations using the subsystems "process" and "machine." Surmann [14] showed a more detailed approach, developing a kinematic machine simulation using CSG material removal for the cutting force calculation and considering the machine statics. Newer developments also include the machine dynamics [15].
Tukora [16,17] developed a cutting force simulation based on a multi-dexel approach. In earlier research [18], the authors have included the real NC unit to indicate the specific NC unit behavior. Other developments include the control loop of the drives [19]. Zhang et al. [20] describe the usage of a milling simulation with augmented reality visualization. A more detailed overview of current research is shown in [21].
Surmann [22] describes the dynamic of cutting processes and simulation of milling and turning processes with a machine simulation. Surmann simulates chattering of different cutting processes based on the used cutting force calculation and a CSG material removal. Wiederkehr et al. [23] showed a point-based tool model to further improve the quality of the simulation by adding effects like runout.

Methodology
2.1 Concept of the modified machine simulation As described above, a plurality of process simulations and kinematic machine simulations has been developed in the recent past. However, most of these simulations are problem specific or insufficient for obtaining a working virtual prototype close to the real process. This paper presents the first steps to improve the simulation described in [24] towards a higher degree of realism regarding the kinematics.
In order to achieve a more realistic kinematic machine simulation, a number of different characteristics have to be considered. According to VDI 3441 [25], the machine characteristics are divided in static, dynamic, and thermal machine behavior. Additionally, the work piece material and the CNC, including the control loops, are influencing the final work piece shape. To implement all these factors, a concept of a real-time machine simulation was developed (Fig. 1).
The concept uses the existing machine simulation with an HIL coupling to a real NC unit. Additional simulation modules are embedded between the existing interfaces. The NC unit communicates the new targeted axis positions to the simulation modules. Each module calculates the resulting changes (e.g., axis position or tool bending) of the machine. These changes are communicated to the kinematic machine model and will result in a new tool center point (TCP). This TCP is then used by the material removal simulation to calculate the new cutting condition. The final resulting axis positions are communicated back to the real NC unit. The cutting parameters are communicated to the cutting force calculation. Figure 2 shows the proof-of-concept model for the later performed verification. As the first verification only considers static

Calculation of cutting forces based on voxel material removal
As shown in Fig. 1, the material removal is part of the kinematic machine simulation. In contrast to the existing machine simulation [18], an overall performance improvement has been made by using a sparse voxel octree. During the cutting process with a cylindrical tool of a prismatic work piece, the initial voxel element is subdivided into eight equally sized sub voxels. For each affected voxel taking part in the cutting process, a new subdivision is reached. This principle is repeated until a predefined recursion depth is achieved. The recursion depth defines the resolution of the voxel model as a function of the work piece size.
The key features of the voxel material simulations are the scalable resolution and a comparably fast and stable cutting algorithm. The disadvantages are comparatively large memory consumption and slow visualization. For the simulation of a continuous milling process, the cutting process has to be discretized, as only time-discrete positions of the work piece can be recorded. There are two possible ways to improve the accuracy of the simulation. Either the time discretion has to be reduced or sweeping has to be performed between the two aligned positions. However, while the time-discrete variant works with simple primitives, the algorithm for sweeping has to include a number of additional forms, which will result in an overall slower performance.
For the cutting force calculation, the voxel material removal provides the cutting volume by summarizing the participating voxels. The cutting force calculation for an end mill is based on Kienzle [7] and is adjusted [26] by the number of simultaneous used cutting edges z ie , the axial angle λ and the product of correction factors ∏K.
Equation (1) uses the basic milling parameter [27]. Analog equations can be made for F CN using the values k f11 and m f from turning. However, in voxel material removal, these values cannot easily be derived. The voxel model will return the minimum radial voxel distance l r from the end mill center, the maximum voxel distance along the end mill l tan , and the volume V. All these values are derived just by counting the number of voxels taking part in the cutting (red voxels), as shown in Fig. 3.   (8).
With N φ ϵℤ and a p ≈l tan ð6Þ

Results
In order to determine the applicability of the presented method, a machine simulation was implemented with real-time voxel material removal, a HIL to a Siemens 840D NC unit, a kinematic and a static machine simulation, including the bending of the end mill. The implanted machine simulation uses Matlab Simulink for the cutting force calculation. Additionally, a multi-body simulation of the machine and the end mill was implemented in Matlab Simulink to calculate the static behavior of the machine. MatlabCom was used as an interface between the kinematic machine simulation and Matlab Simulink.

Test setup
Preliminary tests were executed to determine the applicability of the simulation using a three-axis milling center and three-component dynamometer to record the forces acting on the work piece during the process.
The work piece is an aluminum AlCu4Mg block with the dimensions of 80 × 80 × 60 mm. The milling is performed using an 8-mm solid carbide end mill. The block is milled on one side with an increasing width of cut from 1 to 2.3 mm. The depth of cut, the feed rate, and the spindle speed resulting in the cutting speed are adjusted according to Table 1.
The force calculation uses the specific values of Table 2, which were determined on the same machine with the same material in a previous research project. The used correction values for the rake angle K γ and the tool wear K Ver were set to 1 as the specific cutting values were defined with the same tool. Figure 4 shows the comparison between the implemented force calculation and the measured forces during the process for sample 9. The depicted forces are shown in the X and Y directions of the machine over time. The forces in X direction increased from 150 to 350 Nm on the second edge, as it is expected with the increase of the width of cut. The graph shows conspicuously a fitting of the simulation results in the X direction but a misfit in the Y direction compared with the measured results. This misfit is varying over the samples (Table 1). In general, the misfit is increasing with a reduction of cutting speed. This difference was already stated by Kienzle [7] and is found to be caused by the unconsidered behavior with small areas of cut, e.g., squeezing of material. As the Kienzle model is based on a macroscopic geometric description of the cutting process and empirical data, the general results are better with larger width and depth of cut. With smaller process parameters, microscopic effects induced for example by the tool wear were getting more dominant and have to be considered more carefully. Although the specific  cutting values were measured and defined with the same tool, the tool wear during the process and over all samples is not considered. In general, the process forces will increase with increasing tool wear. Depending on the process parameters, the impact on the forces are different, which can be seen on the force F Y .

Resulting displacement on the work piece
A comparison between a real work piece and a simulated work piece was carried out using the presented machine simulation. The real work piece and the resulting voxel model of the machine simulation were measured. Figure 5 shows the resulting displacements Δx (sum end mill bending u and machine displacement Δx M ) along the milled side of samples 4, 7, and 9, as well as the simulation without any cutting forces.
For samples 7 and 9, the resulting displacement between simulation and the real work piece are nearly identical. Sample 4 shows a bigger difference between the simulated and calculated displacement. A feed per tooth of 0.0078-mm microscopic effect as stated above is more dominant. The shown vibrations in the simulations are not related to chattering of the process and are just a result of the timediscrete record of the milling process and the used material removal. The simulation uses a time-discrete model to calculate the cutting forces. For each step, the cutting force and the displacement are calculated. The calculated displacement is afterwards added to the model and defines the next position of the end mill for the next simulation step. Because of this delayed displacement, the simulation system tents to a time-discrete dependent vibration. With a smaller step size of the recording, the effect would be smaller. Therefore, the simulation has no correlation to the achieved surface of the real work piece.

Conclusion and outlook
This paper describes a new research approach for a machine simulation using voxel-based material removal to calculate the cutting forces during the process. These forces are used in a multi-body simulation of a real machine tool to calculate the resulting displacements of the tool center point. As shown in this paper, it is possible to use this simplified approach to determine the cutting conditions and the resulting displacements with sufficient accuracy. After mending the biggest flaws and performing some additional optimizations, the complete simulation has the potential to achieve real-time capability. Currently, there are several issues to obtain a fully real-time simulation of this kind. The current biggest issue for this solution is the slow communication between the simulation models using MatlabCom, which is not optimized for real-time communication. Although the used voxel model is able to calculate the shape of the work piece in real time, it is not yet possible to achieve real-time visualization, as the voxel model has to be converted into a polygon-based model for fast computational rendering. However, as the visualization is only relevant for the user other solutions, as for example, a delayed visualization is possible. After further optimization of the simulation and its interfaces, a fully fledged benchmark is planned.
Besides resolving the previously stated flaws, the next research focus will lie on the empiric validation of the method on a more complex and realistic work piece and process, and to include working simulation models of dynamic and thermal machine behavior.
Open Access This article is distributed under the terms of the Creative Comm ons Attribution 4.0 International License (http:// creativecommons.org/licenses/by/4.0/), which permits unrestricted use, distribution, and reproduction in any medium, provided you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons license, and indicate if changes were made.