On Modelling of Cutting Force and Temperature in Bone Milling

: Cutting force and temperature are the key factors to be controlled during the orthopaedic surgery which could result in mechanical damage and necrosis of the bone tissue. Mechanistic modelling of the bone cutting process is expected to be an efficient method to understand and control these process challenges. However, due to the special structure and properties of the bone tissue (consist of osteon fibres and interstitial lamellae matrix), the conventional metal cutting models are not applicable in bone cutting process. This paper presents a novel cutting force and temperature mechanistic models for milling of bone. A cutting stress model of bone material was developed which takes into account its anisotropic characteristics based on the orthogonal cutting data. The cutting force coefficients are predicted incorporating the osteon orientation, tool geometry and edge effect with unified mechanics of cutting approach. Furthermore, a model of the induced cutting temperature based on heat flux developed during the process was proposed to predict the temperature distribution on bone cut surface. The experimental results showed a better consistency with the proposed model compared with the conventional Johnson-Cook model under different cutting conditions. A necrosis (potential cell injury from thermal effect) penetration depth was also proposed to evaluate the extent of thermal damage of bone tissue by the developed models. The proposed model can be used to assist the robotic surgery, to optimize the cutting parameters as well as to guide the orthopaedic tool design.


Introduction
Bone cutting is an important procedure in orthopaedic surgery that can span from amputations to finely controlled tissue removal to match the implant surfaces. Nevertheless, this is a challenging process due to bone unique structure, which is consisted of osteons (fibres) and interstitial lamellae (matrix) where the osteons showing a higher strength than interstitial lamellae and inherently results in an anisotropic property.
As an anisotropic material with semi-brittle behaviour during machining, the bone workpiece and its surrounding tissue are easily to be damaged under the high cutting force, which could lead to slow recovery rate or even the secondary injury of the patient as reported by Liao et al. (2017). Noordin et al. (2015) reported that improper cutting of bone could result in severe surface damage which is adverse to the post-surgical recovery, while Denis et al. (2001) also reported that improper cutting can lead to inappropriate contact between bone and prosthesis situation that yields non-uniform stresses on the bone surface. Moreover, the high cutting force would lead to tool failure and breakage during the bone surgery. Liao and Axinte (2016a) reported that cutting with the worn and broken tool would result in secondary trauma of the bone for the avoidance of which they developed a force and acoustic emission-based process monitoring to detect the tool malfunction during the bone cutting. On the other hand, high cutting force could lead to high cutting temperature which would result in thermal necrosis of bone tissue. Moritz et al. (1947) showed that in general when the temperature is above 44°C for 1min the bone could be damaged by necrosis. Further, Hillery and Shuaib (1999) reported that when the cutting temperature exceeds 55°C and for more than 30 seconds the bone tissue is damaged irreparably with necrosis. These two thresholds have been accepted by researchers as the critical temperatures during the orthopaedic surgery. Thus, to minimize the mechanical and thermal damage of the bone tissue during the cutting process, the cutting force and temperature have to be critically controlled.
Predicting the force and temperature in bone cutting is an efficient way to assist robotic surgery, to optimize the cutting parameters, and to guide the orthopaedic tool design. While most of the researchers/surgeons analysed the cutting process of bone based on empirical models or using finite elements (FE) methods as reported by Alam et al. (2011) and Sezek et al. (2012), only a few studied the mechanistic modelling of bone cutting process. Davidson and James (2003) were the first to couple cutting theory with a heat transfer FE simulation to predict the temperature in bone drilling. Lee et al. (2012) and Sui et al. (2014) then developed a mechanistic model of bone-drilling force by considering the influence of drill geometry while the bone drill temperature was also predicted by Lee et al. (2011) and Sui et al. (2015) with utilising finite-difference method. However, these models are only applicable in bone drilling and none of them take into consideration the ploughing phenomena that occur in bone cutting. Moreover, while most of the reports focus on the drilling process, very limited research reports on the bone milling process although it is one of the main processes in bone surgery. To the authors' knowledge, up to now, there has not been a significant attempt to propose mechanistic models for milling forces on bone structures. Most researches dealing with force and temperature in bone milling are focused on experimental investigations to optimise the cutting parameters. Mitsuishi et al. (2004) developed an empirical model for predicting bone milling force considering the density verified for both pig and human bones. Denis et al. (2001) investigated the milling process on defrosted human tibiae with different cutting parameters which showed that while milling force increases with the feed per tooth, the cutting temperature decreased. Sugita et al. (2014) reported a similar conclusion in milling of bovine femurs. While these studies are very useful to practitioners (e.g. surgeons), they provide limited insight on the particular phenomena that govern the bone cutting process.
On the other hand, the current models are either considering the bone shear stress as a constant or by employing Johnson-Cook constitutive model. However, this might not be appropriate since the bone presents a specific structure, i.e. combining the osteon fibres and internal lamella matrix, which leads to an anisotropy of the mechanical properties. Schwiedrzik et al. (2014) showed different deformation modes and mechanical properties (yield strength, elastic modulus etc.) in axial and transverse directions by micro-compression test. These differences inherently lead to crack initiation and propagation in different directions of the bone structure as reported by Koester et al. (2008). Liao and Axinte (2016b) reported that in different cutting directions the values of bone cutting stress change significantly hence their influence of the cutting phenomenon could not be ignored. Sugita et al. (2009) also found that different cutting directions influence the chip morphology as well as the cutting forces. Yeager et al. (2008) found the cutting force in transverse direction is nearly double of the value in parallel cutting direction which leads to largest surface roughness and most severe sub-surface damage of bone structure. Thus, in modelling the bone cutting process, its structural anisotropy has to be taken into consideration.
In this paper, a novel mechanistic model for both force and temperature of bone milling is proposed.
A model of bone cutting stress and friction coefficient was developed considering its anisotropic properties. The orthogonal cutting experiments were employed to evaluate the cutting stress in varying osteon cutting angle. The osteon orientation, tool geometry and edge effect were considered in the model along with the bone material specificity. An arc shape moving heat source was employed for the temperature field prediction of bone workpiece. A series of experiments were conducted to validate the proposed model and comparisons with the results from a conventional Johnson-Cook model were carried out. Based on the proposed model, the necrosis penetration depth could be also evaluated.

Modelling on bone cutting stress under different cutting angle
Since the bone structure mainly consists of osteon fibres and interstitial lamellae matrix which mechanical properties depends on the osteon orientation, during the bone milling process the relative angle between cutting and osteon directions is one of the main factor affecting the mechanics of machining. As shown in Fig. 1, the inclination angle of cutting direction to osteon orientation (fibre cutting angle) changes continuously, the cutting stress in the bone structure is also changing with the osteon orientation cyclically, thus the constant strength assumption and Johnson-Cook model is not applicable in these circumstance. Karpat et al. (2012) calculated the cutting force coefficient of carbon fibre reinforced polymer (CFRP) milling process as a simplified function of fibre cutting angle from experimental milling force data. However, since the cutting force coefficient is also mechanistically related to cutting stress, tool geometry and cutting mechanics data (e.g. friction coefficients, shear angle), Karpat's model is not essentially applicable when applied in bone material with different cutting tools. Thus, to allow modelling the bone milling process considering different cutting mechanisms, i.e. by taking into account the influence of cutting conditions, edge effect, friction as well as the inherent variations of material properties, in this paper a bone cutting stress model is proposed, which considers both the fibre cutting angle (ߠ ௨௧ ) and bone inherent shear strength (߬ ௦ ) as: where ‫ܥ‬ ଵ ǡ‫ܥ‬ ଶ ǡ‫ܥ‬ ଷ and ߠ are constants to be calibrated. The interaction angle between the cutting direction and osteon orientation, can be defined as osteon cutting angle in relation to the instantaneous tool position angle ߠ and the initial osteon angle ߠ according to Jamal and Rahul (2005), expressed as: steons According to Davidson and James (2003) the viscoelasticity of the bone is reflected in the relationship between shear strength and strain rate as in Eq. (3).
where the shear plane angle ∅ can be obtained from Ernst-Merchant cutting theory (Shaw 2005). Thus, considering the osteon cutting direction in bone cutting its cutting stress can be express as: According to our previous study (Liao and Axinte 2016b) the friction coefficient between bone and cutting tool is not changing with the variation of osteon cutting angle but varies with the uncut chip thickness as different chip formation mechanisms occur. Thus, based on the preliminary experimental results, the friction coefficient in bone cutting could be expressed as the function of uncut chip thickness as: where ߤ and C4 are constant friction coefficient and constant exponent respectively that can be calibrated.
In this respect, the true shear stress and friction coefficient in bone machining can be obtained from the orthogonal cutting according to cutting theory (Shaw and Vyas 1998) Thus, with the true shear stress (Eq. (6)) and friction coefficient (Eq. (7)) from orthogonal cutting data the proposed cutting stress (Eq. (4)) and friction coefficient (Eq. (5)) models can be calibrated 8 accordingly. This will then allow the consequent bone milling force and temperature prediction to take into consideration of the material and tool influences (e.g. osteon orientation, tool geometry and edge effect etc.) from cutting mechanics point of view.

Modelling on bone milling force
During the milling process, the cutting forces in each direction can be divided into small incremental forces acting on a differential element (dz) of the milling cutter (Altintas 2000) and can be expressed as: Because of the inherent roundness of the cutting edge ( Fig. 2) coupled with viscoelastic nature of the bone structure, in bone milling process significant ploughing effect has been found ( Fig.2), phenomenon which cannot be ignored when modelling bone cutting process. In conventional metal cutting when modelling ploughing force it is always assumed that there is a dead metal zone (DMZ) formed in front of cutting edge thus to allow the slip-line field analysis as reported by Karpat and Özel (2008). However, as the bone exhibits more of a viscoelastic behaviour the DMZ assumption is applicable when considering ploughing phenomenon of bone cutting. Thus, as shown in Fig. 3, the ploughing force ‫ܨ(‬ ) in bone cutting then can be modelled directly by integrating the force along round tool edge ݈ , and can be expressed as where the friction factor η=0.5cos -1 (µ) (Jin and Altintas 2011) and ploughing angle ߮ expressed as Since the shearing layer (hc) is formed only when the uncut chip thickness is bigger than minimum cutting thickness ℎ , the depth of ploughing layer (hp), as shown in Fig. 3, which is formed below the shearing layer, can be expressed within the ℎ as where the minimum cutting thickness ℎ can be calculated from (Yuan et al.1996): Thus the ploughing force coefficients can be calculated as where the effective rake angle ߛ in tool edge, as shown in Fig. 3, can be calculated as Considering the relative orientation between tool rotating direction and osteon orientation, the shear force coefficient can be calculated from Armarego and Brown model (1969): The uncut chip thickness can be determined by simplified cosine function of cutting angle and feed rate or estimated from cutting trajectory based models. To get a more accurate estimation, the uncut chip thickness is modelled based on the cutting point trajectories according to Spiewak (1995).
Considering the tool run out effect the uncut chip thickness then can be expressed as Taking into consideration of the ploughing depth of Eq. (11), the shear depth then can be calculated Thus, by combining the cutting stress model developed in previous section the bone milling force can be estimated with various cutting conditions (e.g. different osteon orientations and tool geometries).
This will on one hand allow the surgeons to simulate the cutting force of bone milling process before the surgery thus to select the optimized cutting parameters and cutting tools, on the other hand could be also applied in robotic surgeries as well as to help the design of novel cutting tools.

Modelling the bone milling temperature
Since the bone structure can be damaged by necrosis once the temperature exceeds 44°C it is of importance to understand and model this process output thus to allow the optimization of cutting condition or controlling the cutting process. To model the cutting temperature distribution of bone milling process, a moving heat source method can be applied to the tool-workpiece contact area, as shown in Fig. 4, where the heat source is considered as an arc distribution in as the tool-workpiece contact area in milling process is represented as a circular surface (Jin et al. 2002, Richardson et al. 2006). Fig. 4 Moving heat source in milling (Jin et al., 2002) According to Richardson et al., (2006), the relationship between global heat flux (ܳ ௧௧ ) and local heat flux ‫ݍ(‬ (ߠ)) in the arc heat source can be expressed as： where ߠ ௗ is the tool full immersion angle. Since the cutting heat transferred into the workpiece is mainly from the heat energy generated in shearing (ܳ̇௦(‫))ݐ‬ and ploughing (ܳ̇(‫))ݐ‬ area, the global heat source can be integrated from these two heat energies as： Combine Eq. (18) and (19), the local heat flux in the bone milling process then can be expressed as: (20) According to the cutting theory, the heat generated in shearing area is mainly formed from the shear energy and part of it flows into the workpiece, which can be calculated from the of shear force ‫ܨ(‬ ௦ ) and shear speed (ܸ ௦ ) as: where the proportion of heat flows into workpiece, ‫ܥ‬ ௦ , can be expressed according to Shaw (2005) Similarly, the heat generated in ploughing area can be calculated from the where the proportion of heat flows into workpiece, ‫ܥ‬ , can be expressed according to Berliner and Krainov (1991): Thus, the heat flux flow into the workpiece during bone milling process can be expressed as Assuming the bone workpiece is an infinite body, the generation solution of a moving heat source (q) is given as Eq. (26) according to Hahn and Zisik (2012): where the heat source ‫ݍ‬ can be expressed as (27) When calculating the temperature in bone workpiece at point M (x, z), as shown in Fig. 4, the Xaxis projection distance ‫ܦ‬ from M to any point (R, ߠ ) of arc heat source is ‫ݔ‬ − ܸ ‫−ݐ‬ ܴ ‫݅ݏ‬ ݊ߠ , and the spatial distance D is ට ൫‫ܸ-ݔ‬ ‫-ݐ‬R sin ߠ ൯ 2 +൫z-R(1-cos ߠ )൯ 2 , thus the temperature rise at pint M can be expressed as Combining Eq. (25) and (28), the instantaneous temperature rise at point M(x, z) on time t can be get from Thus, with the proposed model the surgeons could predict temperature field in bone material during the milling process thus to optimize the cutting parameters and select the optimized cutting tools to allow restraining the cutting temperature to be below the criterion without damaging the bone tissue.

Experimental setup
In order to obtain real bone cutting stress for calibrating the proposed cutting stress model (Eq. 4) with different osteon cutting angles, orthogonal cutting trials have been performed using an in-house developed 4-axis miniature machine tool to collect the cutting force data with varied osteon cutting angle from 0 to 180 degree under different uncut chip thickness (changing from 0.2µm to 140µm) while the cutting speed was fixed as 100mm/min to reduce the thermal effect during the modelling of bone strength. A custom-made solid carbide cutting tool with rake angle γ=8°, clearance angle α=8° and edge radius re=2µm was employed in orthogonal cutting test. The cutting forces were measured with a 3-component miniature dynamometer (Kistler 9317B) under the sampling rate of 10 kHz. The bone sample was collected from mid-diaphysis of bovine femurs which was kept in physiological saline environment to maintain its mechanical properties. More information of bone orthogonal cutting could be found in previously research (Liao and Axinte 2016b). To validate the cutting force and temperature model, milling tests for measuring cutting forces and temperature close to the cutting zone have been performed on a Makino A55 machine tool (Fig. 5). A 2 flutes solid carbide milling tool with diameter of 4mm, rake γ =10°, clearance angle α =10°, edge radius re=8µm and helical angle λ=30° was used for the milling tests; the tool specifications could be considered relevant for small milled surfaces that are required for medical implants in real surgery application. A series of milling tests were carried out with the axial and radial cutting depth of 2mm under different cutting speeds (6.28-25 m/min), feed rates (30-180µm/tooth) and initial osteon angles (0°and 90°). Since the long axes of femur and its osteons are in the same direction, here 0°and 90°of initial osteon angles could be referred to the parallel and transverses direction of the bone long axes.
The details of cutting parameters can be seen in Table 1 where tests 1-9 and 10-21 were employed for validating milling force and temperature models respectively. K-type thermocouples with diameter of 100µm were embedded into the bone workpiece (1.9mm from the free surface to junction) to allow the real time temperature measurement using a 16-channel data logger (GW instruments Model 100) operating at 10 kHz sampling rate. . Results and discussion

Milling force validation
To calibrate the proposed bone cutting stress and friction model (Eq. (4) and (5)) the orthogonal cutting tests were employed, where the real cutting stress and friction coefficient can be obtained from Eq. (6) and (7). By changing the uncut chip thickness from 0.2µm to 140µm and fitting the cutting force with a regression algorithm, it can be seen from Fig. 6 that when the uncut chip thickness is of small values (<20µm) both the cutting force (Fc) and thrust force (Fth) show a more linear distribution while in zero uncut chip thickness non-zero forces present. However, at this stage there is no chip formed leading to zero shearing force; thus, the intercept of the cutting force distribution vs. UCT could be associated to the ploughing phenomenon.   Osteon cutting angle (°) are bigger than the proposed model while in tangential direction are smaller than the proposed model.
Moreover, the proposed model is more sensitive to the uncut chip thickness since the viscoelastic properties of the bone structure has been considered. The bone milling forces predicted by the proposed model are well in accordance with the experimental results in both amplitude and trends. Fig. 10 and 11 show the predicted and experimental cutting force with 0°and 90°initial osteon cutting angle under test 6 (Vc=6.28m/min, fz=180μm), where the maximum predicted errors of amplitude are 9% and 17% respectively indicating a better agreement in 0°initial osteon cutting angle. This is because the fracture cutting occurs at the tool entrance into the bone structure due to high cutting force while the crack under 90°osteon angle propagates in a severe condition against the cutting direction due to the higher cutting stress and fracture toughness according to Liao and Axinte (2016b). More evidently, it can be seen that under 90°initial osteon cutting angle the experimental cutting force also shows higher values of both amplitude and dynamic component compared with 0°initial osteon cutting angle.

Milling temperature validation
Different from its cutting stress, the thermal conductivity of bone can be considered as isotropic since the osteon orientation does not show an obvious influence on its thermal property (Sui et al. 2015). The thermal conductivity coefficients of bone and milling tool are listed in Table 3. It is also worth to note that in real surgery application the tool is usually made by medical steel (e.g. 316L) which the thermal diffusivity is slightly smaller than the applied tool in this experiment thus the heat transferred to the bone workpiece from ploughing area will be smaller. To validate the modelling results of bone milling temperature, the thermocouple was embedded into the cutting layer of bone sample, as shown in Fig. 13, where the thermocouple would be cut off by the cutting edge thus to collect the real temperature from the cutting area. The average temperature within the cutting time was then calculated to be considered as the cutting temperature of bone milling. It is also needed to be noted from Fig. 13 that the ambient temperature of bone workpiece is 20°C while during the real surgery process the tissue temperature will be warmer and real surgery temperature will be higher. For better comparison between the outputs of the model and the experiments, the surface temperature from modelling is defined as the average value within the circular area of highest temperature at 100µm diameter (equal to thermocouple diameter). On comparing prediction errors, as shown in Table 4, the modelled results are in good agreements with the experiments with maximum and minimum errors of 27.25% and 6.37% respectively. As the thermocouple response time is 120ms some errors could be associated with the experimental setup itself.
Examples of predicted workpiece temperature fields of bone milling (test 10-15) are shown in Fig. 14, from which an explicit image can be drawn on how the bone material suffers from the heating problem during cutting process based on the isotherms. Interestingly, it can be seen clearly the high temperature area (i.e. above 44°C when necrosis occurs) penetrated, depending on the cutting parameters, into the bone structure to 300-1000µm under the free surface, explaining well why the necrosis phenomenon happens in the orthopaedic surgery. These results are consistent with the histochemical observation of in vivo rabbit surgery from Eriksson et al. (1984), where the bone necrotic zone around the cut surface was around 500µm. Since the bone milling temperature decreases from the cut surface to the inner area while there would be a high risk of thermal damage on bone material when the cutting temperature is above the necrosis threshold (44°C, Moritz et al. 1947), to evaluate the risk of necrosis during bone milling process, a necrosis penetration depth could be defined as the distance from the machined surface to the 44°C point of the isotherms, as shown in Fig. 15.  Fig. 16 and 17 show the surface temperature (measured in the bone cutting surface as shown in Fig.   13), global heat flux and necrosis penetration depth as a function of cutting speed and feed rate. It is interesting to see that the heat flux increases with the cutting speed (Vc) and feed rate (fz) while the surface temperature from both experiments and modelling showing a converse trend. The tendency of increasing heat flux is the same with metal cutting as the cutting force increases with the cutting speed and feed rate. However, the decrease of temperature with increasing cutting speed and feed rate is fundamentally different from the conventional metal cutting phenomenon. This is due to the low thermal conductivity of bone compared with the tool material, with which when increasing the moving speed of heat source the amount of cutting heat transmitted into to bone workpiece reduces significantly as the dwell time for the heat transmitting is shorten while most heat is brought away by tool and chips. Moreover, it also can be seen that when increasing the cutting speed to 25m/min or the feed rate to 180µm the necrosis penetration depth reduces to 200µm which means only one osteon (around 0.2mm in diameter) is damaged in the depth direction while in low cutting speed or feed rate several osteons are damaged due to the high necrosis penetration depth. Thus, with the understanding of temperature distribution and necrosis penetration depth under different cutting conditions, the proposed model can be used to assist the robotic surgery, to optimize the cutting parameters as well as to guide the orthopaedic tool design. However, it is also needed to be noted that the modelled process is under dry cutting condition while when the coolant is applied the heat will be reduced and the necrosis penetration depth will be reduced.

Conclusions
In this paper, a novel mechanistic model on bone milling process were developed for predicting the cutting forces and temperature. Firstly, the bone cutting stresses were modelled by considering its anisotropy in mechanical properties. The orthogonal cutting experiments were employed to evaluate the cutting stress in varying osteon-tool relative engagement angle and uncut chip thickness. Secondly, the force model in milling was developed while the osteon orientation, tool geometry and edge effect were taken into consideration in the model with its cutting mechanics. Thirdly, an arc moving heat source with cutting energy based heat flux was employed in the temperature field prediction of bone workpiece in milling process. The model predicted results agrees well with the experimental results, indicating a potential to be used to assist the robotic surgery, to optimize the cutting parameters as well as to guide the orthopaedic tool design.
The main finding of the paper can be summarised as follows: (a) The bone cutting stress can be modelled as a circular function of the osteon cutting angle in which the weakest and strongest cutting strengths of bone are at 30°and 90°osteon cutting angle respectively. The friction coefficient between bone and cutting tool is only related to the uncut chip thickness while the influence from osteon strengthening effect is not significant. The cutting stress and friction coefficient model can be calibrated by orthogonal cutting experiments under different osteon cutting angles. The ploughing effect is significant and needs to be considered in the modelling process due to the tool edge radius and viscoelasticity of bone material.
(b) The milling force from experimental results show a higher amplitude and dynamic component at 90°initial osteon cutting angle compared with 0°initial osteon cutting angle. This effect was captured accurately with the proposed model while the conventional Johnson-Cook model showed its limitation to capture this phenomenon since the osteon strengthening effect is not considered.
(c) A maximum prediction error of 27.25% and minimum error of 6.37% on cutting temperature were reached with proposed model. The predicted temperature distribution is also consistent with the histochemical observation of in vivo rabbit surgery on previous study. With proposed model of bone milling process, the surgeons could simulate the cutting force and temperature under different cutting conditions (e.g cutting and feed speeds, tool geometries) and directions (e.g. parallel, inclined or perpendicular to the bone long axis) before the surgery thus to select the optimized cutting parameters and cutting tools to allow reducing the cutting force and temperature. Moreover, this model could be also applied in robotic surgeries as well as to help the design of novel cutting tools.