Bastien Toubhans
email: bastien.toubhans@ensam.eu

Guillaume Fromentin
email: guillaume.fromentin@ensam.eu

Fabien Viprey
email: fabien.viprey@ensam.eu

Habib Karaouni
email: habib.karaouni@safrangroup.com

Théo Dorlin
email: theo.dorlin@safrangroup.com

Machinability of inconel 718 during turning: Cutting force model considering tool wear, influence on surface integrity

Keywords:

Machining accuracy can be compromised by elastic workpiece deformation and subsurface residual stress introduction during cutting. In order to anticipate the impact of cutting forces and surface integrity evolutions on finished surface and its geometrical errors, it is necessary to better understand the influence of cutting conditions and tool wear. In this study, machinability of Inconel 718 using a round carbide tool in finish turning conditions is assessed. Cutting forces evolution during tool life are analysed and accompanied by advanced investigations of cutting phenomena. An original mechanistic cutting force model is developed, identified and tested. It includes the effect of tool wear over time in its local formulation. This model allows predicting cutting forces evolution along tool pass for a wide range of finishing cutting conditions. Furthermore, a thorough analysis of residual stress profiles at different tool wear levels is led. It features quantitative results for fresh and worn tools. A study on the influence of cutting parameters and tool wear on residual stress profiles in the machining affected zone is highlighted.

Introduction

Accuracy is one of the major concerns when machining high-added value parts. Among the significant amount of possible sources of geometrical errors encountered in machining of thin parts, this study focuses on two phenomena directly linked to the cutting action. The first one is the workpiece elastic deformation under load during machining. As the workpiece is pushed away from the cutting tool, an undercut defect appears on the finished part. The second one is the alteration of surface integrity during machining. Heterogeneous plastic deformations, arising from multiple mechanisms during cutting, result in residual stresses in a thin layer under the machined surface. [START_REF] Brinksmeier | Residual stressesmeasurement and causes in machining processes[END_REF] find that these residual stresses can generate deformations by disturbing the workpiece equilibrium. Although these two phenomena are often neglected when machining rigid parts, they may have a prominent role when dealing with low stiffness parts. [START_REF] Toubhans | Prediction of form error during face turning on flexible Inconel 718 workpiece[END_REF] quantify their influence. In the case of face turning on thin discs, 80 % of the total geometrical error is attributed to elastic part deformation while 20 % is linked to machining induced residual stresses and stress rebalancing following material removal. These two phenomena depend on the material considered.

The material of this study is Inconel 718, a nickel-based alloy. It is widespread in turbine engines for its remarkable mechanical and corrosion resistance properties at high service temperatures. Although the machinability of this material has been improved over the years, it remains a 'hard to cut' material. Manufacturers are mainly struggling with low attainable cutting speeds and rapid tool wear. The presence of hard carbide particles in the microstructure, the low heat conductivity and the Built Up Edge (BUE) phenomenon are evoked by [START_REF] Polvorosa | Tool wear on nickel alloys with different coolant pressures: Comparison of Alloy 718 and Waspaloy[END_REF]. These three characteristics, inherent to nickel alloys, are identified in the literature as responsible for rapid tool wear. When machining Inconel 718 with carbide tools in finishing conditions, [START_REF] Devillez | Cutting forces and wear in dry machining of Inconel 718 with coated carbide tools[END_REF] establish that the main wear mechanisms are abrasion and adhesion. According to the same authors, on the one hand, wear manifestations due to abrasion are mainly flank wear and notching of the cutting edge. On the other hand, adhesion is responsible for coating flaking and removal of material from the rake face. [START_REF] Arrazola | Correlation between tool flank wear, force signals and surface integrity when turning bars of Inconel 718 in finishing conditions[END_REF] observe that tool wear and surface degradations appear faster above a certain flank wear around 0.15 mm. Tool life is often determined by the amount of flank wear VB which can be measured following the ISO 3685 (1993). A value equal to 0.3 mm is currently considered in literature for carbide finishing tools to deem a tool worn.

In addition, cutting forces tend to rise with tool wear. When using round inserts for finishing Inconel 718 in turning, [START_REF] Arrazola | Correlation between tool flank wear, force signals and surface integrity when turning bars of Inconel 718 in finishing conditions[END_REF] confirm that trend for the three cutting forces in the machine referential. However, the passive cutting force, normal to the machined surface, appears to be much more sensitive to tool wear. Indeed, the passive cutting force level is multiplied by 5 during the tool life. The physical explanation given by [START_REF] Grzesik | Investigation of tool wear in the turning of Inconel 718 superalloy in terms of process performance and productivity enhancement[END_REF] is the ploughing effect favoured by tool wear. Therefore, modelling cutting forces evolution in respect with cutting conditions and tool wear is a challenge. Numerous contributions are available in the field of cutting force modelling. Different approaches emerged with time such as empirical, analytical and numerical as [START_REF] Arrazola | Recent advances in modelling of metal machining processes[END_REF] relate. The present study focuses on the mechanistic approach which is part of the empirical kind. In a mechanistic model, local forces are expressed as a linear function of the cut thickness. Among primary developments, [START_REF] Armarego | Drilling with flat rake face and conventional twist drills-I. Theoretical investigation[END_REF] first introduce the edge discretisation methodology in order to calculate global forces by summing local forces along the cutting edge. [START_REF] Armarego | Computer based modelling of popular machining operations for force and power prediction[END_REF] propose a two components local cutting force model describing what is known as the "cutting effect" and the "edge effect". [START_REF] Kapoor | Development of mechanistic models for the prediction of machining performance: model building methodology[END_REF] state the importance of modelling chip flow direction on accuracy. Indeed, consistent chip flow direction suggests a dependency between the cutting edge segments which is absent in Armarego's formulation. [START_REF] Molinari | A new thermomechanical model of cutting applied to turning operations[END_REF] propose a mechanical cutting model based on the elementary chip equilibrium including a resulting force on its lateral sides, i.e. tangential to the cutting edge. Following this idea, [START_REF] Dorlin | Generalised cutting force model including contact radius effect for turning operations on Ti6Al4V titanium alloy[END_REF] and [START_REF] Chérif | Phenomenological study of chip flow/formation and unified cutting force modelling during Ti6Al4V alloy turning operations[END_REF], with a mechanistic approach, completed the original local force model used by [START_REF] Armarego | Computer based modelling of popular machining operations for force and power prediction[END_REF] by adding a third local force tangent to the cutting edge. It allows better prediction of cutting forces when machining with round insert or with the tool nose.

Several attempts have been made to take tool wear into account in cutting force models. A commonly found approach is to model the VB evolution during tool life and then incorporate it in a cutting force model. This method is used by [START_REF] Zhu | A generic tool wear model and its application to force modeling and wear monitoring in high speed milling[END_REF] with a time dependant wear model in high speed milling. Besides, an analytical model describing the evolution of flank wear is presented by [START_REF] Das | Analytical model of progression of flank wear land width in drilling[END_REF]. [START_REF] Grzesik | Investigation of tool wear in the turning of Inconel 718 superalloy in terms of process performance and productivity enhancement[END_REF] develops a global empirical model considering the removed chip volume and the cutting speed in order to model the flank wear during tool life. Other approaches do not consider tool wear directly but develop a mechanistic model with addition of a time dependent term. The instantaneous force is expressed as A + B(t) where, A is the cutting force level with no wear, and B(t) is a time dependant term representing the effect of wear. Lacalle et al. (2017) present a time dependent additional term considering cutting speed and feed rate in a global cutting force model without considering the nose radius effect. This last approach is considered in this study by means of a local cutting force model. This way, the nose radius effect, which is of prime importance when using round tools, is considered.

Concerning surface integrity, a wide consensus appears to describe the residual stress profiles in the machining affected layer of Inconel 718. [START_REF] Devillez | Dry machining of Inconel 718, workpiece surface integrity[END_REF] depict the typical profile as tensile on the surface, decreasing toward a minimum compressive stress at shallow depth around couple hundreds of micrometres, and increasing back toward the bulk residual stress. Among the multiple parameters influencing the residual stress profile shape, [START_REF] Devillez | Dry machining of Inconel 718, workpiece surface integrity[END_REF] list the most influent as tool geometry, tool wear, cutting conditions and lubrication. [START_REF] Sharman | An analysis of the residual stresses generated in Inconel 718 TM when turning[END_REF] describe that tool wear has the largest influence on surface integrity. Using a worn tool leads to higher tensile peak at the surface and deeper compressive stress. [START_REF] Sharman | An analysis of the residual stresses generated in Inconel 718 TM when turning[END_REF] argue that a worn cutting tool loses its ability to cleanly cut the material and tends to rub against the machined surface which leads to higher temperature and more plastic deformation. These observations are confirmed in a more recent study by [START_REF] Sharman | The effect of tool nose radius on surface integrity and residual stresses when turning Inconel 718 TM[END_REF] using round tools. According to [START_REF] Javadi | Inconel 718 Superalloy Controlled Surface Integrity for Fatigue Applications Produced by Precision Turning[END_REF], there are no linear relationships between residual stress profiles and process parameters. Indeed, [START_REF] Thakur | State-of-the-art in surface integrity in machining of nickel-based super alloys[END_REF] show that it is widely accepted in the community that the final residual stress profile is the result of a competition between thermal and mechanical inputs during the cutting action. Thermal and mechanical loads respectively favour tensile and compressive stresses.

The literature related to the turning of Inconel 718 in finishing condition reveals that cutting forces and surface integrity are significantly influenced by cutting conditions and tool wear. In this study, the focus is set on developing a cutting force model taking tool wear into account during face turning. The second section presents the experimental procedure. In the third section, a cutting force model is improved in order to predict the cutting forces when using a fresh tool. This model is generalised for different cutting speeds. In the fourth section, this model is complemented by an additional term to take tool wear into account. The fifth section focuses on the study of surface integrity depending on tool wear levels and cutting conditions. Finally, residual stress profiles analysis are compared to literature observations.

Study parameters and experimental setup

Material and cutting conditions

The material machined in this study in an Inconel 718 treated to 45 HRC. The cutting tool is a round carbide insert with a 4 mm nose radius. The cutting insert has a PVD micro grain (Ti,Al)N + TiN coating specially designed for the finishing of superalloys. When mounted in the tool holder, the rake, inclination and clearance angle are respectively 0°, 0°and 7°. Flood cooling conditions are applied using a water based emulsion with 5 % oil content.

The study focuses on face turning operation with constant cutting speed. Three cutting speeds are used, V c ∈ {35; 52.5; 70} m/min. Feed f and depth of cut a p are chosen in order to explore a range of maximum cut thickness h max ∈ [0.05; 0.14] mm representative of finishing conditions. A description of the cut section geometry is given in Fig. 1 where the position of the cutting edge is given in the P r plane (cf. ISO 3002-1 (1982)) during one spindle revolution. The engaged portion of the cutting edge is bounded by θ min and θ max , which depends on the depth of cut a p , the feed f and the nose radius r ε . The cut thickness h varies along the cutting edge of a round insert and reaches its maximal value h max close to θ max .

Fig. 1. Cutting edge engagement during face turning with a round cutting tool.

Experimental setup

The cutting tests are performed on a 3-axis CNC lathe. The global cutting forces are measured in-process in the machine coordinate system (MCS) with a piezoelectric dynamometer and a charge amplifier. In-situ measurements of flank wear is performed using a digital microscope to avoid disassembling the tool and the resulting repositioning issues when using round inserts. VB max is measured (cf. ISO 3685 (1993)) and is called VB in this study. The goal is to monitor the cutting forces evolution during the tool life in order to predict it with a cutting force model presented in the next two sections.

3. Generalised cutting force model when using a fresh tool

Local forces model and edge discretisation methodology

This section focuses on the mechanistic approach to model cutting forces when using a fresh tool. This model is generalised for different cutting speeds. The experimental data are gathered through short dressing operations. Nine cutting tests with different engagements, summed up in Table 1, are performed for three different cutting speeds, giving 27 tests. The diameter range during these tests is [115; 84] mm. No significant wear is measured on the cutting inserts after these tests. Indeed, VB does not exceed 60 μm which is consistent with flank wear values encountered during the running-in period. The cutting forces are averaged on three spindle turns during steady state. The measurements during the 27 cutting tests are available in Table 1. Considering the large nose radius and shallow depths of cut used, the edge orientation is such that F p modulus is important compared to F f one. As expected, at similar chip load, cutting forces decrease when the cutting speed increases. Therefore, the cutting speed is taken into account in the model.

The local cutting forces model, displayed in Eq. ( 1), is calculated using the cutting edge discretisation methodology. The engaged cutting edge is discretised into N segments of equal lengths b as depicted in Fig. 2. The local force model used in this study is based on Armarego's early work on cutting force modelling [START_REF] Armarego | Computer based modelling of popular machining operations for force and power prediction[END_REF]. Armarego's model considers two local force components acting on a cutting edge segment f v and f h , respectively axial and radial to the cutting edge as depicted in Fig. 2. Each component is expressed as a linear function of the local cut thickness h. An additional component, tangential to the cutting edge, is added. It ensures that the local forces on each segment contribute to the global chip flow direction. This tangential component f o is proportional to η cf,i , the local gap between the average chip flow direction θ cf and the angular position θ i of the considered segment as shown on Fig. 2. The average chip flow direction θ cf is computed using a geometric formulation shown in Eq. ( 1). This formulation displays convincing correlation with experimental tests in Chérif's study when the inclination angle is equal to zero [START_REF] Chérif | Phenomenological study of chip flow/formation and unified cutting force modelling during Ti6Al4V alloy turning operations[END_REF]. Finally, the effect of the cutting speed is added to the model by multiplying f v and f h by a dimensionless term V c /V o raised to a certain power (n v and n h ) where V o = 52.5 m/min. The cutting speed has a limited influence of the average chip flow direction according to Chérif's observation. As a consequence, the effect of V c is not added to f o .

⏟ ∫ ∫ ⎜ ⎟ ⎜ ⎟ ⎧ ⎨ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎩ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ = + ⎛ ⎝ ⎞ ⎠ = + ⎛ ⎝ ⎞ ⎠ = = - = k k k k k f b h V V f b h V V f b η h η θ θ θ h θ θ dθ h θ dθ ( .
).

( . ).

. . . ( ).

( )

cv
∑ ∑ ⎧ ⎨ ⎪ ⎪ ⎪ ⎩ ⎪ ⎪ ⎪ = = - = + = = = F f F f θ f θ F f θ f θ
.

sin ( ) . cos ( )

. cos ( ) . sin ( )

c i N v i f i N h i i o i i p i N h i i o i i 1 , 1 , , 1 , , MCS MCS MCS | | | (2)
Table 1 Cutting conditions and mean cutting force measurements during short face turning tests. V c is in m/min. Fig. 2. Cutting edge discretisation, local and global forces.

f a p h max F c_mean (N) F f_mean (N) F p_mean (N) mm/rev mm mm V c = 35 V c = 52.5 V c = 70 V c = 35 V c = 52.5 V c = 70 V c = 35 V c = 52.5 V c =

Cutting force modelling: results and discussions

The local forces model contains 7 unknown constants to identify, respectively k cv , k ev , k ch ,k eh , k o , n v and n h . The identified coefficients are given in Table 2. The negative values of the exponents associated to the cutting speed terms are consistent with experimental observations. It can be noted that cutting speed has similar influence on both components f h and f v . Each of the 27 cutting tests gives 3 global forces data samples. As a result, the Residual Degree Of Freedom (RDOF) for f v model is equal to 24 and the one for f h and f o is equal to 50. The RDOF is calculated as the difference between the number of equations and the number of unknown parameters. For example, according to Eq. ( 2), there are 54 equations (27 F p -samples and 27 F f -samples) to compute f h and f o and 4 unknown constants giving a RDOF equal to 50. Cross validations tests are performed by using only a portion of the experimental data-set to identify the model and the rest of the data-set for validation simulations. The identification results are similar when using half of the data-set, two thirds of the data-set and the entire data-set. Consequently, the well distributed data-set combined with the high model RDOF leads to a robust identification.

The errors between simulated and measured forces are displayed in Fig. 3, with black points, for the three components F c , F f and F p . Fig. 3 is a box plot where half of the data is contained in the boxes and the other half between the whiskers (dotted lines) outside the boxes. The whisker length is equal to 1.5 times the interquartile range or box length. Values outside the whiskers are called outliers and are symbolised by red crosses.

Although the maximum relative error on F f seems important compared to the other global forces, the levels of F f are 6 times lower than F c and 5 times lower than F p . Consequently, the absolute error on F f remains low. For the same reason, F f is more sensitive to measurement noise which explains why outliers are found only for F f. This discrepancy is due to the optimisation strategy which minimises the absolute error. Indeed, the global forces with the highest magnitudes are favoured by this optimisation strategy.

The presented three components model is compared to Armarego's two component model. The prediction is unchanged for F c as the f v local force is the same for the two models and only appears in F c calculation. However, when including the third component f o in the model, the average relative error for F p is divided by two and the average relative error on F f is slightly degraded (by 0.3 %). The relative error on F p is shown in Fig. 4 for the 9 cutting conditions used at V c = 52.5 m/min using both models. Both errors are evenly distributed around zero as dictated by the optimization strategy. However, three groups of red dots, corresponding to the two component model, are clearly identified. There is one group for each depth of cut used. Experimentally, changing the depth of cut modifies the length of the engaged cutting edge and hence the chip flow direction. In definitive, Fig. 4 shows that the two component model lacks the description of the influence of chip flow direction on cutting forces. The same observation is made for the two other tested cutting speeds. Considering Eq. ( 1), the contribution of f o component increases when the gap between the angular position of the considered segment and the average chip flow direction increases. This contributes to orientate F f and F p resultant along θ cf in the P r plane, tangent to the rake face.

To conclude, a cutting force model is developed in order to give cutting forces when using a fresh cutting tool. It is valid on a certain range of cutting speeds and tool engagements representative of finishing conditions in Inconel 718 using a carbide cutting tool. The addition of a third component f o in the local force model allows for better prediction of the passive cutting force F p , normal to the machined surface, which is the principal responsible for elastic deformation of the workpiece under load. However, tool wear has the first role in the evolution of cutting forces during a normal tool life. The next section focuses on analysing cutting forces evolution with wear and improving the current model to take tool wear into account.

Cutting force model taking wear into account

In this section, experimental observations are presented concerning tool wear and cutting forces evolution during tool life. A specific paragraph deals with the built-up-edge effect on these two phenomena. Drawing conclusions from the experimental observations, a cutting force model taking tool wear into account is developed and tested

F c f v (N) 24 k cv (N. mm -2 ) 2651 k ev (N. mm -1 ) 5 7 n v (-) -0.136 F f , F p f h (N) 50 k ch (N. mm -2 ) 1726 k eh (N. mm -1 ) 104 n h (-) -0.144 f o (N)
k o (N. mm -2 . rad -1 ) 53807

Fig. 3. Relative and absolute errors between experimental and modelled cutting forces using the fresh tool cutting force model in Eq. ( 1). As the metal is being cut, some material stagnates in front of the cutting edge. It is commonly called a Built-Up Edge (BUE) and can be seen in Fig. 5.a above the red line materialising original cutting edge. Some of it accumulates on the tool rake face forming a Built-Up Layer (BUL). Meanwhile, some of it is evacuated between the cutting tool and the machined surface. This material ends up smeared on the machined surface and the tool flank face as one can observe on the Scanning Electron Microscope (SEM) image on Fig. 6. Built-Up Edge is observed for every tested cutting condition. To sum up, increasing cutting speed and feed tends to favour flank wear unevenness and more dramatic changes of the cutting edge geometry. These changes in cutting edge geometry may cause important variations of cutting forces as analysed in the next paragraph.

Cutting forces evolution during tool life

On Fig. 7, cutting forces are represented as a function of machining time during the test at the following cutting conditions, V c =70 m/min, f =0.1 mm/rev and a p =0.5 mm. The discontinuities on the curves are the transitional regimes between each pass which have been cut out for legibility. This particular test took 5 passes to reach the critical flank wear level VB > 0.3 mm. A typical cutting force evolution is observed with a short running-in period followed by a controlled wear region and finally by the tool end of life. Cutting forces tend to increase with tool wear. However, the passive cutting force F p , normal to the machined surface, is more sensitive to it than F c and F f . It is consistent with observations made by [START_REF] Arrazola | Correlation between tool flank wear, force signals and surface integrity when turning bars of Inconel 718 in finishing conditions[END_REF] in comparable conditions. During a normal tool life, F p is multiplied by 3-7 while F c and F f are multiplied by 2-4 depending on the cutting conditions.

Tool wear is known to be a phenomenon with poor repeatability. To ensure that cutting forces can be modelled, reproducibility tests are performed for certain cutting conditions. On Fig. 8, passive force F p evolution is illustrated during three repetitions of the same wear test. It appears that in the studied case, the running-in and controlled wear regions display repeatability. However, the inflection point towards the tool end of life, materialised by dashed lines, is quite variable. There are several factors able to disturb a stable wear evolution. Indeed, when the tool is wearing, the cutting edge integrity is compromised and prone to catastrophic failures, such as collapsing of the cutting edge, which may cause the inflexion point occurrence to vary. Same remark applies to notching of the cutting edge which creates a starting point for catastrophic failure. Moreover, notching randomly occurs at any time during tool life. Furthermore, entering and exiting the material are critical times for wear evolution as transitional regime may imply Tool wear tends to favour BUE formation. Indeed, as the cutting edge becomes dull, more material seems to stagnate in its vicinity. At low cutting speeds, BUE has a significant effect on cutting forces evolution. This phenomenon is detailed in the next paragraph.

Built-Up Edge effect on cutting forces evolution

When performing wear tests at the lowest cutting speed 35 m/min, perturbations of the cutting forces evolution occur (Fig. 9). Instead of having continuous evolution as observed in Fig. 7, long period oscillations on all cutting forces, discontinuity between passes and peaks when engaging and disengaging the tool appear. The slight discontinuity between the first passes observed on Fig. 9 may be explained by the grain size gradient along the workpiece radius. Indeed, the average grain size is equal to 12 μm close to the outer diameter and 50 μm close to the inner diameter. It results in a moderate hardness gradient along the radius, 45 HRC for fine grain structure and 43 HRC for coarser grain structure. Additionally, as stated by [START_REF] Olovsjö | The effect of grain size and hardness of wrought Alloy 718 on the wear of cemented carbide tools[END_REF], the adhesion related phenomenon are more pronounced in coarser grain structure. It could be another factor explaining the rise of cutting forces observed in the coarser grain structure at the end of the first passes on Fig. 9.

During the long period oscillations, cutting forces drop significantly, up to 45 % for F p , compared to their expected levels. Several possible causes, such as heterogeneous hardness along the radius, poor chip evacuation or temporary loss of lubrication due to chip accumulation in the cutting zone, are discarded due to the magnitude of the oscillations. Process damping might have occurred at lower speed considering the ploughing effect due to BUE formation and the extent of flank wear. However, the low frequency oscillations, under 0.1 Hz, cannot be linked to process damping as it is far from the process frequencies. These phenomena, which manifest above certain wear level, are thought to be adhesion related. Indeed, at lower cutting speeds, the BUE may become stable and act as a substitute cutting edge. As illustrated on Fig. 10, in presence of a BUE, it is supposed that the rake angle becomes more positive and the contact between the flank face and the machined surface is interrupted. The consequence of these two factors is a cutting forces drop. Indeed, the rubbing of the tool flank face on the machined surface is a major contributor to global cutting forces. The oscillations would be resulting from the alternation of evacuation and reformation of the BUE.

This assumption is supported by further analysis. This test is repeated and the machined surface is scanned with a laser profilometre along the workpiece radius after a pass during which perturbations occurred. The passive cutting force during this pass is superimposed with the machined surface topography on Fig. 11. Reliefs are observed on the machined surface. Valleys, 30 μm deep on average, are measured and correlated with forces drops. These observations checks with the BUE hypothesis as BUE acts as a substitute cutting edge made of stagnating material which slightly increases the depth of cut.

From the above observations, it is recommended to avoid cutting conditions where BUE is stable as it will greatly degrade accuracy and surface roughness. The first choice is to increase the cutting speed. During the performed tests, increasing V c from 35 to 52.5 m/min made these phenomena disappear as BUE may have become unstable. Considering the significant increase of all cutting forces during a normal tool life, tool wear influence is added to the cutting force model in the next section.

Cutting force model development

Formulation of the model

In section 3, a mechanistic model designed to give cutting forces when using a fresh tool is presented. Previous observations made it clear that cutting forces increase during tool life. Hence, in this section, the previous model is complemented by an additional term to take wear into consideration. Additionally, only the running-in and controlled wear periods are modelled as the inflexion point toward the tool end of life is unpredictable. Tool flank wear VB is measured during the tests but is not used to model the effect of wear on cutting forces evolution. Indeed, measuring flank wear implies a stoppage of the cutting process. VB measurements are performed after each pass. As a consequence, few measurements are done during the tests to ensure realistic test durations and to limit the number of tool engagement cycles which favour quick wear. In addition, wear manifestations such as smearing on the flank face, notches or wear unevenness (see Fig. 5) make it more complicated to accurately measure VB. Instead of using a parameter directly linked to tool wear, such as VB, the cutting force evolution is modelled in relation to a process parameter which can be continuously measured during the cutting tests.

Parameters such as the cumulated uncut chip length L m and removed volume V m or cumulated machining time t m are considered. L m is the real distance travelled by the generative point of the tool (i.e. an Archimedes' spiral during facing). V m is L m multiplied by the undeformed chip cross sectional area. The idea is to represent the actual amount of work done by the cutting edge during its life. Exploratory developments are made by adding a power function term of L m or V m to the fresh tool model. However, L m and V m are global quantities and there is no bijection between them and the actual cut area represented by f and ap. Indeed, L m is independent of a p and V m is independent of f (cf. Appendix A for supplementary equations). Consequently, these models are discarded.

The most advanced model is presented in Eq. ( 3). An additional term is added to f h and f v formula. So the global forces are expressed as A + B(t) where A is the cutting force level when the tool is fresh and B (t) represents the force increase due to wear over time. The f o component is untouched as chip flow direction is unaffected by tool wear according to experimental observations in similar conditions [START_REF] Chérif | Phenomenological study of chip flow/formation and unified cutting force modelling during Ti6Al4V alloy turning operations[END_REF]). Considering constant lubrication conditions, cutting speed is the most influent process parameter on tool life followed by feed and depth of cut. In this model, the feed and depth of cut effects are encapsulated in the local cut thickness h i which represents the chip load. Finally, the additional term takes the form of a power function of machining time, cutting speed and local cut thickness with 8 new coefficients to identify, respectively k wv , n hv , n tv , n vv , k wh , n hh , n th and n vh . It can be noted that f v and f h local formulations now contain two terms depending on V c /V o with different power exponents. Indeed, cutting speed has significant influence on both cutting forces initial level and evolution over time. (3)

⏟ ⏟ ⏟ ⏟ ⎜ ⎟ ⎜ ⎟ ⎜ ⎟ ⎜ ⎟ ⎧ ⎨ ⎪ ⎪ ⎪ ⎪ ⎪ ⎩ ⎪ ⎪ ⎪ ⎪ ⎪ = ⎛ ⎝ ⎜ ⎜ ⎜ + ⎛ ⎝ ⎞ ⎠ + ⎛ ⎝ ⎞ ⎠ ⎞ ⎠ ⎟ ⎟ ⎟ = ⎛ ⎝ ⎜ ⎜ ⎜ + ⎛ ⎝ ⎞ ⎠ + ⎛ ⎝ ⎞ ⎠ ⎞ ⎠ ⎟ ⎟ ⎟ = k k f b k h k V V h t V V f b k h k V V h t V V f b k η h ( . ). . . . ( . ). . . . . . . wv n n n wh n n n v i cv i ev c o n A fresh tool i c o B t wear effect h i ch i eh c o n A fresh tool i c o B t wear effect o i o cf i i , ( ) ( ) ( ) , ( ) ( ) ( ) ,

Identification of the model and discussions

The seven coefficients which are identified in section 3.2 for the fresh tool model are kept in the new model. The same edge discretisation is used to identify the eight new coefficients. To do so, experimental data are gathered by discretising the cutting forces evolution curves over machining time. The test at V c = 35 m/min and f = 0.1 mm/rev is not used for the identification as the cutting force evolution is highly disturbed by the BUE related phenomenon presented in section 4.1.3.

The identification results are presented in Table 3. The simulated cutting force evolutions are plotted in Fig. 12 against the experimental samples for each test. An average of 20 samples per test and per cutting force (F c , F f and F p ) is selected. As a result, the model is highly overdetermined as evidenced by the residual degree of freedom (cf. Table 3). These samples are represented by discrete markers in Fig. 12. The sampling is denser in the running-in region as the cutting forces evolution is less linear than during the controlled wear region. As it is shown in Fig. 12, the model is able to predict the evolution of the three cutting force components in the machine referential. In the next paragraphs, this model is evaluated on its ability to predict the cutting force evolution, especially the passive cutting force which is the most sensitive to tool wear.

The maximum and median values of absolute and relative errors are given in Fig. 13. While the average error values are quite low, the discrepancy between the maximum and average values and the presence of numerous outliers are questionable.

In order to better visualise the results, a box plot is displayed in Fig. 14 to evaluate the prediction of the passive cutting force F p . The black dots represent the actual distribution of relative error between the experimental and simulated forces for each test. Additional values such as the average error, median, 25th and 75th percentiles are displayed on each box plot.

In Fig. 14, the median value indicates the ability of the model to simulate a correct force evolution over time due to tool wear. If the median is close to zero, the slope is correct (Fig. 12 (b,f)). If the median is over or under zero, the simulated slope respectively overestimate or underestimate the slope as in Fig. 12(e,h). Overall, the model is able to simulate the passive cutting force evolution during the running-in and controlled wear region within 10 % of relative error in the tested range of cutting conditions.

During cutting tests 'a' and 'b', the simulated slopes are correct but the error values are scattered and present outlying values. These isolated errors have two possible origins. The first is linked to discontinuity between passes that appear during the tests at lower cutting speeds, where built-up edge disturbs the cutting forces evolution as mentioned in Section 4.1.3. It means that the samples are on either side of the simulated evolution curves as shown in Fig. 12a after 15 min of machining. It is evidenced by the negative and positive isolated errors on the 'a' distribution in Fig. 14. The second is linked to random effects that manifest as slope modification, named bends, in the experimental cutting force evolution curves. These bends may randomly appear at higher speeds and disturb the nearly linear evolution of cutting forces during the controlled wear period. They are thought to be due to modification of the cutting edge geometry when notches or cutting edge collapse occur. Such a bend appears at the beginning of the second pass during the test at V c = 52.5 m/min and f = 0.35 mm/rev (Fig. 12e) causing the simulated slope to be wrong from the beginning. As it can be seen in Fig. 12(e,g,h), higher cutting speeds and feed are impractical as the controlled wear region duration does not exceed two minutes.

In the studied case, transient regimes do not have a significant effect on cutting force evolutions during the running-in and controlled wear regions. Indeed, no significant cutting force discontinuity is observed between passes (see Figs. 7 and8). For this reason, to a certain extent, the model is valid for different workpiece diameters which would have more or less enter/exit tool periods.

To conclude this section, a local force model able to give the evolution of cutting forces relative to tool wear during most of the tool life is developed. This model is valid on the tested range of finishing conditions during the running-in and controlled wear regions. However, variability in bulk material characteristics and unpredictable catastrophic wear manifestations are significant causes of gaps between experiment and model. In addition to its influence on cutting forces evolution, cutting tool wear has a significant impact on surface integrity which will be the next section focus.

Surface integrity analysis

Experimental approach

In this section, residual stress profiles of the machined surface are presented for various cutting conditions and tool wear levels. The specimens are 35 mm thick cylinders, considered as rigid bodies, to ensure no deformation due to machining induced stress or residual stresses rebalancing following material removal. Face-turning operations are performed until the tool is deemed worn. The diameter range during these tests is [155; 30] mm. Flat steps are machined on the specimen face in order to preserve the machined surface integrity at strategic time as it is depicted on the top view in Table 5. Five cutting conditions are tested based on experimental observations made in Section 4.1. Cutting conditions as well as wear levels and machining times (t m ) associated to each step are listed in Table 4. Compared with the classic cutting force evolution during tool life displayed in Fig. 7, the first step matches the end of the running-in period. The last step represents the tool end of life. The intermediate step corresponds to the end of the controlled wear region before the inflexion point. Stress measurements in the circumferential and radial directions are performed on the machined flat steps by an X-Ray diffraction apparatus (DRX). The axial stresses are not studied. The uncertainty of residual stresses measurements is estimated at ± 40 MPa based on repeatability tests on a standard specimen. Profiles are built by incrementally eroding the machined surface to allow for stress measurements at different depths using electropolishing. Generating a flat hole bottom using this technique is difficult, especially for deep holes. The uncertainty of depth measurements is evaluated at 10 % of the target depth by laser measurements.

Effect of tool wear on residual stress

The residual stress profiles matches the literature observations stated in the introduction [START_REF] Sharman | The effect of tool nose radius on surface integrity and residual stresses when turning Inconel 718 TM[END_REF]. The comparison between the profiles when using a fresh tool and a worn tool demonstrates a strong trend as shown in Fig. 15. This graph represents the circumferential residual stresses for 3 wear levels during the test at V c = 52.5 m/min and f = 0.2 mm/rev. Machining with a worn tool increases the tensile stress on the surface and the compressive stress at the compressive peak. Moreover, the affected depth increases by a factor 5 and the compressive peak is deeper, for every tested cutting conditions and in both measured directions. While the difference between a fresh tool (light green circle curve) and a worn tool (red square curve) is remarkable, the intermediate profile (dark green diamond curve) shows moderate difference relative to the fresh tool profile. This observation is comparable to the cutting forces evolution during the test shown in the top right corner of Fig. 15. The steep increase in cutting forces occurs when the tool reaches a certain dullness which corresponds to a flank wear level VB between 0.15 and 0.2 mm. After this inflexion point, the cutting edge becomes extremely dull as evidenced by the rightmost flank wear picture in Fig. 15.

The worn cutting edge geometry is the major factor influencing stress profile shape. However, it seems that stress profiles are moderately impacted by tool wear until VB reaches the tipping point between 0.15 and 0.2 mm. The effect of cutting conditions combined with tool wear is analysed in the next paragraph.

Effect of cutting conditions on residual stress

The effect of cutting conditions and tool wear on residual stress profiles shape is shown in Table 5. Global trends are described in response to an increase in feed or cutting speed. The latter are given for both tool states: a run-in tool and a worn tool in both measurement directions. Indicators such as the magnitude of the tensile and compressive peaks, depth of the compressive peak and affected depth are considered. The trends are symbolised by blue arrows modifying the typical residual stress profile (solid red line). Due to the discrete nature of the residual stresses measurements and the depth limitations, it may be risky to infer trends. When no trend can be extracted from the data, an orange cross is placed instead of an arrow. For example, the depth of the machining affected zone may be hazardous to analyse for the following reasons. When the stress does not reach equilibrium at the deepest measured point or if the measurement spacing is too large. However, the difference of affected depth when using a run-in tool and a worn tool is noticeable and respectively evaluated at 0.1 mm and 0.5 mm on average.

Residual stress profiles are displayed in Fig. 16 to complete Table 5. The subfigure Fig. 16 a and c depict the effect of cutting speed on residual stresses in the circumferential direction in function of tool wear. The Fig. 16 b andd show the effect of feed on residual stresses in the radial direction. Indeed, residual stresses appear to be more sensitive to cutting speed or feed evolutions in their respective directions, i.e. circumferential for V c and radial for f.

As illustrated by Fig. 16 a, an increase in cutting speed tends to lower the surface tensile stresses when using a run-in tool. This observation is consistent with the literature when using round carbide tools in wet or dry finishing conditions and for moderate cutting speeds under 80 m/min [START_REF] Devillez | Cutting forces and wear in dry machining of Inconel 718 with coated carbide tools[END_REF]; [START_REF] Pawade | An investigation of cutting forces and surface damage in high-speed turning of Inconel 718[END_REF] and [START_REF] Sharman | An analysis of the residual stresses generated in Inconel 718 TM when turning[END_REF]. In addition, when using a run-in tool, the profiles shapes are comparable while V c is greater than 52.5 m/min as shown in Fig. 16 a. However, the profile at V c = 35 m/min is incomparable as it presents a deeper affected zone. The ploughing effect resulting from BUE formation could explain this discrepancy as the material stagnating at the cutting edge tends to create more compressive stress ahead of the tool.

When using a run-in tool at V c = 52.5 m/min, increasing the feed tends to favour tensile stresses on the surface and to reduce the magnitude and depth of the compressive peak in the radial direction (see Fig. 16 b). Same trends are observed by [START_REF] Sharman | An analysis of the residual stresses generated in Inconel 718 TM when turning[END_REF] in comparable cutting conditions except for the compression peak depth which increases with feed rate. [START_REF] Sharman | The effect of tool nose radius on surface integrity and residual stresses when turning Inconel 718 TM[END_REF] discuss that the increase in surface tensile stress is linked to an increase in cutting force and possibly heat generation due to the higher chip volume being removed.

Moreover, compressive surface residual stresses are obtained when using a fresh tool at high cutting speeds and low feed namely V c = {52.5; 70} m/min and f = 0.1 mm/rev. This state of stress is desirable to enhance fatigue strength. However, machining with these conditions when the tool is worn leads to high surface tension as shown in Fig. 16 c and d which is detrimental for fatigue behavior.

Some stress profiles display brutal changes in the first microns under the machined surface as shown in Fig. 16 c (blue curve) and Fig. 16 d (light green curve). These rapid evolutions are consistent with a high heat flux at the surface which can induce metallurgical modifications at the vicinity of the machined surface.

To summarize, tool wear has a moderate impact on residual stresses profiles shape while the flank wear VB is lower than a critical value between 0.15 mm and 0.2 mm. Above this critical value, the profiles are dramatically changed displaying high surface tension and deep affected zones. Using high cutting speeds and low feed tend to favor compressive stress at the surface and shallow affected depth in the controlled wear region.

Conclusions

This study presents a new local formulation of cutting force model considering the effect of tool wear over time. In addition, the effect of tool wear on surface integrity is analysed. The key conclusions are:

• Including a third component, relative to chip flow direction, in the local force model, allows better prediction of cutting forces when using round inserts. In particular, the relative error on the passive cutting force is divided by 2.

• The new local formulation accurately predicts the cutting forces evolution over a wide range of finishing parameters during the running-in and controlled wear periods. Variability of tool wear evolution is still an obstacle to develop fully robust models.

• Tool wear has a critical impact on cutting forces and surface in- tegrity. In this study, the passive cutting force is multiplied by up to 7 over tool life. Concerning surface integrity, the affected depth is multiplied by 5 and reaches 0.5 mm when machining with a worn

  calculated on each segment according to the local force model. The global forces are then computed by projection and summation of the local forces in the MCS according to Eq. (2). Finally, the model coefficients are optimised by least square minimisation of the absolute error of the three global forces (F c , F f and F p ). The Nelder-Mead simplex algorithm is used to optimise the solution. The data-set used to identify the model contains all measured data gathered from the 27 tests described above.

∑

Fig. 4 .

 4 Fig. 4. Influence of the third local force component f o on F p passive force prediction -V c = 52.5 m/min.

Fig. 5 .

 5 Fig. 5. Flank wear manifestations depending on cutting conditions.

Fig. 6 .

 6 Fig. 6. SEM image of a tool with wear manifestations and BUE -V c = 35 m/ min, f = 0.1 mm/rev, a p = 0.5 mm.

Fig. 7 .

 7 Fig. 7. Typical cutting forces evolution during tool life -V c = 70 m/min, f = 0.1 mm/rev, a p = 0.5 mm.

Fig. 8 .

 8 Fig. 8. Evolution of F p passive force during three repetitions -V c = 70 m/min, f = 0.1 mm/rev and a p = 0.5 mm.

Fig. 9 .

 9 Fig. 9. Unstable cutting forces evolution during wear test -V c = 35 m/min, f = 0.1 mm/rev and a p = 0.5 mm.

Fig. 10 .

 10 Fig. 10. Effect of built-up edge formation on cutting geometry and tribology.

Fig. 11 .

 11 Fig. 11. Correlation between cutting forces drops and negative reliefs on machined surface.

Fig. 12 .

 12 Fig. 12. Experimental vs modelled cutting forces during running-in and controlled tool wear region.

Fig. 13 .

 13 Fig. 13. Relative and absolute errors between experimental and modelled cutting forces using the cutting force model in Eq. (3).

Fig. 14 .

 14 Fig. 14. Relative error distribution on F p passive force modelling for each test.

Fig. 15 .Fig. 16 .

 1516 Fig. 15. Circumferential residual stress profile at different wear levels and cutting force evolution -V c = 52.5 m/min, f = 0.2 mm/rev, a p = 0.5 mm.

Table 2

 2 Identified coefficients of the local force model for fresh tool.

	Global forces	Local forces	RDOF	Identified components

Table 3

 3 Identified coefficients of the local force model taking wear into account.

	Global forces	Local forces	RDOF	Identified coefficients
	F c	f v (N)	155	k wv (-)	131
				n hv (-)	0.58
				n tv (-)	0.56
				n vv (-)	1.32
	F f , F p	f h (N)	314	k wh (-)	147
				n hh (-)	0.46
				n th (-)	0.71
				n vh (-)	2.19

Table 4

 4 Machining time and flank wear levels associated to every step.

	V c (m/min)	35	52.5	52.5	52.5	70
	f (mm/rev)	0.1	0.1	0.2	0.35	0.1
	a p (mm)	0.5
	Step 1	t m = 5.2 min	t m = 1.36 min	t m = 0.68 min	t m = 0.24 min	t m = 0.44 min
		VB = 0.13 mm	VB = 0.1 mm	VB = 0.09 mm	VB = 0.07 mm	VB = 0.1 mm
	Step 2			t m = 7.03 min	t m = 3.86 min
				VB = 0.16 mm	VB = 0.12 mm
	Step 3	t m = 47.7 min	t m = 20.2 min	t m = 11.5 min	t m = 6.9 min	t m = 8.83 min
		VB = 0.3 mm	VB = 0.3 mm	VB = 0.3 mm	VB = 0.3 mm	VB = 0.29 mm

Table 5

 5 Residual stress profiles evolution in relation with cutting conditions and tool wear.

Acknowledgements

This research did not receive any specific grant from funding agencies in the public, commercial, or not-for-profit sectors. The authors wish to thank T. Bergey from ArianeGroup for sharing his expertise on X-ray diffraction measurements.

• Wear has limited influence on residual stresses profiles shape while the flank wear VB remains under a critical value between 0.15 mm and 0.2 mm. This study brings physical analysis and model necessary to deal with the challenge of machining workpiece with low stiffness. Moreover, the presented cutting force model offers promising foundation for the development of strategies to compensate workpiece deflexion. Indeed, the update of tool trajectories to accommodate with workpiece elastic deformation is a meaningful perspective in order to improve geometrical quality of thin machined parts. Furthermore, in order to anticipate distortions linked to machining induced residual stresses, a major helpful step would be to establish further model of residual stress based on mechanistic approach as developed by [START_REF] Su | Modeling of residual stresses in milling[END_REF].

CRediT authorship contribution statement

Declaration of Competing Interest

The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.


Laureline L Chaise
email: llchaise@gmail.com

Dominic J Mccafferty

Adélie Krellenstein

Susan L Gallon

William D Paterson

Marc Théry

André Ancel

Caroline Gilbert

Environmental and physiological determinants of huddling behavior of molting female southern elephant seals (Mirounga leonina)

Keywords: Phocids, aggregation behavior, heat loss, habitat type, meteorological variations, stomach temperature Phocids, aggregation behavior, heat loss, habitat type, meteorological variations, stomach temperature

While endotherms can rely on their insulation to reduce heat loss to adapt to cold environments, renewing of fur during molt impairs insulation while they have to perfuse the periphery to support epidermal tissues. The southern elephant seal Mirounga leonina undertakes an annual catastrophic molt while fasting on land in a wet, windy and cold environment. However, southern elephant seals show characteristic aggregation patterns that are predicted to reduce high metabolic costs during the molt. Between 2012 and 2016, 59 female elephant seals were tracked on land during their molt to study their aggregation behavior in relation to molt stage, habitat type and local weather conditions. Infrared thermography and stomach temperature loggers were used to observe variation in body surface and internal temperature in relation to molt stage and aggregation behavior. We found that thermal constraints varied during the molt, with a peak in surface temperature during the mid-stage of the molt. Wallows (mud pools) appear as favorable habitat to aggregate while molting. Indeed, wallows offered a warmer microclimate with greater ground temperature and lower wind speed. Moreover, there was a greater proportion of aggregated seals and larger group size in wallows. These aggregation patterns in wallows were influenced by local weather such that a greater proportion of seals were located in the center of the aggregation, and larger group size occurred during days of unfavorable meteorological conditions. We also observed a higher proportion of seals at mid-stage of molt amongst aggregated seals compared to isolated individuals. This aggregation behavior may reduce the cost of thermogenesis as surface body temperature and stomach temperature were cooler by 1.0°C and 1.5°C, respectively, in aggregated compared to isolated seals. As a consequence, huddling behavior may be thermally advantageous for female southern elephant seals during the molt.

-ΔT h = T h -T a gradient of surface temperature of the aggregation (huddle)

A C C E P T E D M A N U S C R I P T

Introduction

In order to maintain a high and relatively constant body temperature, endotherms living in a cold environment must increase their insulation to compensate for heat loss. This adaptive insulation allows them to have an extended thermoneutral zone, without increasing their metabolic rate [START_REF] Scholander | Adaptation to Cold in Arctic and Tropical Mammals and Birds in Relation to Body Temperature , Insulation , and Basal Metabolic Rate[END_REF][START_REF] Schmidt-Nielsen | Animal Physiology: Adaptation and Environment[END_REF][START_REF] Willmer | Environmental Physiology of Animals, Second[END_REF][START_REF] Clarke | Principles of Thermal Ecology: Temperature, Energy, and Life, First[END_REF]. Most Antarctic and Sub-Antarctic birds and mammals forage at sea and breed or molt while fasting on land, thereby experiencing contrasting periods of energy use and environmental conditions. Living in different physical and thermal environments requires morphological, physiological and behavioral adaptations. For example, pinnipeds are adapted to heat conservation in a cold environment because of low surface area to volume ratio, thick layer of subcutaneous blubber, and vascular peripheral anastomoses regulating cutaneous blood perfusion [START_REF] Schmidt-Nielsen | Animal Physiology: Adaptation and Environment[END_REF][START_REF] Mauck | Thermal windows on the trunk of hauled-out seals: hot spots for thermoregulatory evaporation?[END_REF]. However, behavioral adaptations to minimize heat loss are also widespread in the wild. Indeed, social thermoregulation (huddling), widely used by endotherms, allows metabolic savings by reducing thermoregulatory costs [START_REF] Gilbert | One for all and all for one: The energetic benefits of huddling in endotherms[END_REF]. Reducing heat loss by huddling permits reallocation of energy savings to other physiological processes (e.g. reproduction or growth).

However, during the molt, heat loss increases as blood flow to peripheral tissues must be maintained to support tissue growth and regeneration.

Southern elephant seals Mirounga leonina experience an annual 'catastrophic' molt lasting one month, not only renewing their hair but also their cornified epidermis [START_REF] Ling | The skin and hair of the southern elephant seal, Mirounga leonina (L.): III. Morphology of the adult integument[END_REF]. During this period, molting elephant seals are likely to be more sensitive to heat loss because of their inability to avoid peripheral vascular circulation. The molt is particularly costly in this species: molt metabolic rate in southern elephant seals is estimated to be 2-3 times greater than resting metabolic rate, and appears to be more energetically costly than in northern elephant seals (Mirounga angustirostris; [START_REF] Worthy | Molt energetics of the northern elephant seal (Mirounga angustirostris)[END_REF][START_REF] Boyd | Water flux, body composition, and metabolic rate during molt in female southern elephant seals (Mirounga leonina)[END_REF].

A C C E P T E D M A N U S C R I P T

6 However, decreases in metabolism during the molt have been described in other related species when captive during the same period (harbour seals Phoca vitulina: [START_REF] Ashwell-Erickson | Metabolic and hormonal correlates of molting and regeneration of pelage in Alaskan harbor and spotted seals (Phoca vitulina and Phoca largha)[END_REF][START_REF] Rosen | Correlates of seasonal changes in metabolism in Atlantic harbour seals (Phoca vitulina concolor)[END_REF]. While molting, female body mass loss averages 4-5 kg per day in southern elephant seals, mainly fat from blubber metabolism [START_REF] Boyd | Water flux, body composition, and metabolic rate during molt in female southern elephant seals (Mirounga leonina)[END_REF][START_REF] Hindell | Body mass loss of moulting female southern elephant seals, Mirounga leonine, at Macquarie Island[END_REF][START_REF] Carlini | Mass changes during their annual cycle in females of southern elephant seals at King George Island[END_REF][START_REF] Postma | Spatial variation in female southern elephant seal mass change assessed by an accurate non-invasive photogrammetry method[END_REF]. [START_REF] Paterson | Seals like it hot: Changes in surface temperature of harbour seals (Phoca vitulina) from late pregnancy to molt[END_REF] showed that heat loss of phocids increases during the molt, with an increase in skin temperature due to perfusion by vasodilation through the blubber layer to supply nutrients for epidermis renewal [START_REF] Ashwell-Erickson | Metabolic and hormonal correlates of molting and regeneration of pelage in Alaskan harbor and spotted seals (Phoca vitulina and Phoca largha)[END_REF]. [START_REF] Feltz | Thermal requirements in vitro of epidermal cells from seals[END_REF] demonstrated in vitro that the epidermal cells of phocids require a minimum temperature of 17°C, and an optimal temperature of 37°C, to grow. The high metabolic rate observed in southern elephant seals while molting could be related to the fact that they molt in a cold environment with increased thermoregulatory costs.

Southern elephant seals aggregate in large groups, mostly in mud pools (wallows), while molting on land [START_REF] Laws | The elephant seal (Mirounga leonina Linn.) II. General, social and reproductive behavior. Falkland Islands Dependencies Survey[END_REF][START_REF] Boyd | Water flux, body composition, and metabolic rate during molt in female southern elephant seals (Mirounga leonina)[END_REF][START_REF] Chaise | Local weather and body condition influence habitat use and movements on land of molting female southern elephant seals (Mirounga leonina)[END_REF] and this behavior seems to be influenced by local weather conditions (Cruwys and Davis 1995; [START_REF] Chaise | Local weather and body condition influence habitat use and movements on land of molting female southern elephant seals (Mirounga leonina)[END_REF]. In California sea lions Zalophus californianus, huddling behavior increases in colder weather and allows animals aggregated in the middle of a group to maintain a higher surface temperature than the substrate compared to isolated individuals [START_REF] Liwanag | Thermal benefits of aggregation in a large marine endotherm : huddling in California sea lions[END_REF]. In the same way, huddling individuals of most mammal and bird species maintain a higher and more constant core and surface body temperatures than isolated individuals [START_REF] Gilbert | One for all and all for one: The energetic benefits of huddling in endotherms[END_REF].

However, a higher gradient between body surface temperature and ambient temperature would result in increased heat loss and associated thermoregulatory cost [START_REF] Canals | Energetics and Geometry of Huddling in Small Mammals[END_REF][START_REF] Mccafferty | Estimating metabolic heat loss in birds and mammals by combining infrared thermography with biophysical modelling[END_REF]. In contrast, [START_REF] Gilbert | Body temperature changes induced by huddling in breeding male emperor penguins[END_REF] showed that microclimate created by huddling (up to 37.5°C within tight huddles) allows male emperor penguins Aptenodytes forsteri to save energy during their breeding fast through decreases in core temperature. We therefore predict that aggregated adult southern elephant seals benefit from huddling, either

A C C E P T E D M A N U S C R I P T 7
by the maintenance of a higher body temperature allowing a more rapid molt, or a lowering of thermal gradients between surface and ambient temperature, allowing energy savings during the molt.

The aim of this study was therefore to evaluate thermal consequences of aggregation in molting females, using infrared thermography as a non-invasive technique to measure body surface temperature, concurrently with the recording of internal temperature, of free-ranging pinnipeds [START_REF] Mauck | Thermal windows on the trunk of hauled-out seals: hot spots for thermoregulatory evaporation?[END_REF][START_REF] Mccafferty | Factors influencing the radiative surface temperature of grey seal (Halichoerus grypus) pups during early and late lactation[END_REF][START_REF] Norris | Environment and activity affect skin temperature in breeding adult male elephant seals (Mirounga angustirostris)[END_REF][START_REF] Paterson | Seals like it hot: Changes in surface temperature of harbour seals (Phoca vitulina) from late pregnancy to molt[END_REF][START_REF] Liwanag | Thermal benefits of aggregation in a large marine endotherm : huddling in California sea lions[END_REF][START_REF] Codde | The effect of local weather conditions on the behavior of molting southern elephant seals, Mirounga leonina (L)[END_REF]. We expect aggregation behavior to be more intense in wallows, where elephant seals have been observed aggregating [START_REF] Laws | The elephant seal (Mirounga leonina Linn.) II. General, social and reproductive behavior. Falkland Islands Dependencies Survey[END_REF][START_REF] Boyd | Water flux, body composition, and metabolic rate during molt in female southern elephant seals (Mirounga leonina)[END_REF][START_REF] Chaise | Local weather and body condition influence habitat use and movements on land of molting female southern elephant seals (Mirounga leonina)[END_REF], compared to other habitats (grass and rocky beach), as wallows could be a warmer habitat for molting. We also expect that aggregation behavior increases in the middle of the molt (corresponding to a peak in surface temperature; [START_REF] Paterson | Seals like it hot: Changes in surface temperature of harbour seals (Phoca vitulina) from late pregnancy to molt[END_REF] and when weather conditions are deteriorating [START_REF] Liwanag | Thermal benefits of aggregation in a large marine endotherm : huddling in California sea lions[END_REF]. We therefore predict that aggregation allows elephant seals to modulate their internal and surface temperatures to reduce the cost of thermoregulation and increase their rate of molt in order to minimize mass loss during the molting fast. We defined three molt stages, assessed by the percentage of old hair/skin shed (0 %: no old hair shed to 100 %: all old hair shed; ± 10 %; Fig. 1): this was then divided into initial stage (0-40 % of old hair shed; still largely covered with old hair), mid-stage (50-80 %; most of old hair shed and new hair still not grown) and final stage (90-100 %; new hair growing). We recorded whether the observed seal was aggregated or isolated, where an aggregation was defined when at least two elephant seals were in physical contact. In an aggregation, an elephant seal was considered in a peripheral position (P) when only one of its sides was in physical contact with others, otherwise it was recorded as in a central position (C). We also calculated an aggregation score (number of aggregated seals/total number of seals) for each transect and quadrat, the size of each aggregation (number of aggregated seals) and the C/P

Individual data

Between 2012 and 2016, 59 adult females were captured at the initial stage of molt, then tracked during 4.7 ± 2.9 days (range 1-13 daily observations), and 53 were recaptured at the final stage of molt (Table 1), 8.6 ± 4.3 days after their first capture (range 3-20 days).

Elephant seals were captured on the colony, anaesthetized using tiletamine and zolazepam [START_REF] Mcmahon | Field immobilisation of southern elephant seals with intravenous tiletamine and zolazepam[END_REF][START_REF] Chaise | Implantation of subcutaneous heart rate data loggers in southern elephant seals (Mirounga leonina)[END_REF]) and tagged on one hind flipper with plastic identification tag (Dalton Tags, UK). Females were weighed at capture and recapture (HST Mini-Weigher, 0-1000 kg ± 0.5 kg, HST Scales UK Ltd) to calculate body mass loss (kg.d -1 ).

All captured seals were equipped with VHF transmitters (Series MM300 Marine Mammal Headmount, model MM340B, 7.1 x 3.5 x 2.1 cm; 92 g, Advanced Telemetry Systems, USA)

to track them on land, and nine individuals were equipped with stomach temperature pills and

A C C E P T E D M A N U S C R I P T

time-depth recorders (TDR-STP-207D; MK10-L/SPLASH10-309, 76 x 56 x 32 mm, 125 g, Wildlife Computers, USA) between 2014 and 2016. Stomach temperature pills were placed in the stomach under anesthesia using a lubricated flexible tube, and recorder tags were fixed on the head with epoxy bi-composed glue Araldite®. Stomach pills were set up to record stomach temperature (T stom ) every 10 s. Temperature was recorded during 4.9 ± 2.7 days (range 1-10 days) before the signal was lost, due to natural passage of pill through the gut. For each observation of a tracked elephant seal, digital and thermal images were taken, and molt stage and aggregation status (aggregated or isolated) were determined. In order to compare stomach temperatures between aggregated and isolated individuals, we selected stomach temperature data at the time when individuals were observed and photographed (six individuals out of nine equipped). As we recorded paired data (in both positions: aggregated vs. isolated) for very few individuals (half of our stomach temperature data were unpaired) we assumed independence for statistical analyses. We calculated the individual aggregation rate (number of observations in aggregation/total number of observations). We also recorded meteorological variables (air and ground temperature, relative humidity, wind speed and solar radiation) close to the seal.

2.2.Data analyses

Thermal images

Thermal images were analyzed using the software ThermaCAM® Researcher Pro 2.10 (FLIR Systems, USA). For each image, we specified air temperature (i.e. measured air temperature T a, o C), relative humidity (%), distance (m) and mammal pelage emissivity of 0.98 [START_REF] Humes | Variability of emissivity and surface temperature over a sparsely vegetated surface[END_REF][START_REF] Norris | Environment and activity affect skin temperature in breeding adult male elephant seals (Mirounga angustirostris)[END_REF][START_REF] Mccafferty | Estimating metabolic heat loss in birds and mammals by combining infrared thermography with biophysical modelling[END_REF]. For measurements on caught individuals, mean body surface temperature T b ( o C) was determined by fitting a

A C C E P T E D M A N U S C R I P T

13 polygon around the visible body of the seal, and by measuring mean fore flipper surface temperature T f ( o C) from spot measurement at the axillary, as this is an important thermal window [START_REF] Mauck | Thermal windows on the trunk of hauled-out seals: hot spots for thermoregulatory evaporation?[END_REF][START_REF] Nienaber | Surface temperature patterns in seals and sea lions: a validation of temporal and spatial consistency[END_REF] (Fig. 2A). For transects and quadrat data, mean surface temperature of the aggregation T h ( o C) was determined by fitting a polygon around all visible bodies of aggregated seals and T f from spot measurement from visible axillary of each seal (Fig. 2B). Obvious wet seals and images out of focus were discarded from analysis. We used thermal gradients (ΔT b,f,h = T b,f,h -T a (°C)) for thermal analysis of surface temperatures. was found between grass and beach habitats (N = 73, χ² = 24.92, df = 2, P < 0.0001; Table 2).

A C C E P T E D M A N U S C R I P T

Differences in ground temperature were found between all three habitats, wallows being the warmest (N = 73, χ² = 35.19, df = 2, P < 0.0001) while air temperature was not different between the three habitats (N = 147, χ² = 0.33, df = 2, P = 0.85; Table 2). Wind speed was lower in wallows compared to grass habitat (N = 145, χ² = 6.24, df = 2, P = 0.04; Table 2) but not when compared to beach habitat. Relative humidity and solar radiation were similar between all three habitats (N = 144, χ² = 0.03, df = 2, P = 0.99; N = 139, χ² = 2.4, df = 2, P = 0.31). ). The thermal gradient of body surface temperature was greatest by 2.0 ± 3.6°C during mid-stage compared to initial (estimate ± SD = 0.12 ± 0.06, z = 2.18, P = 0.03) and final stages (estimate ± SD = 0.12 ± 0.05, z = 2.55, P = 0.01) but no significant difference was found between initial and final stage of molt (estimate ± SD = 0.006 ± 0.04, z = 0.17, P = 0.87; Fig. 3).

Habitat

3.2.Influence of molt, habitat and weather on aggregation behavior

Distribution of molt stages varied between aggregated and isolated seals observed on transects and the quadrat (D = 0.18, P < 0.0001). We observed a higher proportion of seals at

A C C E P T E D M A N U S C R I P T 21

Huddling behavior of molting female elephant seals (population based-data) was influenced by habitat type and local weather conditions. The proportion of aggregated seals (i.e. aggregation rate, N = 176; Fig. 5A) was less in grass (mean ± SD = 0.5 ± 0.2; estimate ± SD = -0.58 ± 0.08, z = -7.41, P < 0.0001) and beach habitats (0.6 ± 0.2; estimate ± SD = -0.41 ± 0.07, z = -5.74, P < 0.0001) compared to wallows (0.9 ± 0.2), and less in grass compared to beach habitats (estimate ± SD = -0.16 ± 0.06, z = -2.63, P = 0.009; Fig. 5A). The model also estimated that aggregation rate was marginally higher during days of "bad weather" compared to days of "good weather", for all habitats (estimate ± SD = 0.12 ± 0.07, z = -1.88, P = 0.06) but this result was non-significant. When we compared aggregation rate with weather index for each habitat, we observed that seals aggregated more in grass habitat during "bad weather" days (estimate ± SD = 0.24 ± 0,11, z = 2.15, P = 0.03; Fig. 5A). Indeed, when we analyzed the effect of meteorological parameters, we observed that seals aggregate more in grass habitat with increased wind speed (estimate ± SD = 0.05 ± 0.02, z = 2.30, P = 0.02) and less with increased solar radiation (estimate ± SD = -0.0004 ± 0.0002, z = -2.49, P = 0.01). In the same way, our results showed that elephant seals aggregate less in beach habitat with increased solar radiation (estimate ± SD = -0.0003 ± 0.0001, z = -1.97, P = 0.049) and more in wallows with increased relative humidity (estimate ± SD = 0.005 ± 0.002, z = 2.18, P = 0.03).

The number of elephant seals per aggregation (i.e. aggregation size, N = 754) was less in grass (mean ± SD = 3.1 ± 1.5; estimate ± SD = -1.03 ± 0.15, z = -6.66, P < 0.0001) and beach habitats (4.8 ± 6.9; estimate ± SD = -0.78 ± 0.15, z = -5.36, P < 0.0001) compared to wallows (9.4 ± 8.7), and less in grass compared to beach habitats (estimate ± SD = -0.25 ± 0.10, z = -2.44, P = 0.02; Fig. 5B). The model also estimated that aggregation size was greater during "bad weather" days compared to "good weather" days, regardless of habitat types (estimate ± SD = 0.32 ± 0.16, z = 2.05, P = 0.04). Moreover, the model estimated that the difference in aggregation size between wallows and grass habitat was increased during days

A C C E P T E D M A N U S C R I P T

of "bad weather". Indeed, when we analyzed the effect of the weather for each habitat, aggregation size increased in wallows with increased relative humidity (estimate ± SD = 0.01 ± 0.006, z = 2.33, P < 0.05; Fig. 5B) and the same tendency was observed in beach habitat (estimate ± SD = 0.01 ± 0.006, z = 1.95, P = 0.05). In comparison, aggregation size in grass habitat marginally increased with increased wind speed (estimate ± SD = 0.02 ± 0.01, z = 1.85, P = 0.06), but this result was non-significant.

To study the influence of habitat and weather on C/P ratio, we considered aggregation size as a covariate the model, since the C/P ratio increased with aggregation size (estimate ± SD = 0.07 ± 0.005, z = 13.90, P < 0.0001; N = 689). The model estimated that there were less seals in a central position compared to peripheral positions in wallows compared to grass habitat (estimate ± SD = -1.82 ± 0.23, z = -8.03, P < 0.0001) and beach habitat (estimate ± SD = -0.42 ± 0.17, z = -2.45, P = 0.01) but no difference between grass and beach habitats (estimate ± SD = 0.09 ± 0.13, z = 0.71, P = 0.48). Moreover, the difference in C/P ratio between wallows and grass habitat increased with increased aggregation size (estimate ± SD = 0.29 ± 0.03, z = 8.27, P < 0.0001). There was no general effect of weather index on aggregation C/P ratio (estimate ± SD = 0.03 ± 0.13, z = 0.20, P = 0.84). However, when analyzed per habitat, aggregation C/P ratio increased in grass habitat with increased relative humidity (estimate ± SD = 0.02 ± 0.006, z = 2.38, P = 0.02), and this effect was more important when aggregation size increased (estimate ± SD = 0.006 ± 0.002, z = 2.34, P = 0.02). In comparison, C/P ratio marginally increased in wallows, but non-significantly, with increased wind speed (estimate ± SD = 0.05 ± 0.03, z = 1.77, P = 0.08), and there was no effect of weather in beach habitat.

A C C E P T E D M A N U S C R I P T

3.3.Aggregation behavior and body temperature

A C C E P T E D M A N U S C R I P T

In wallows, gradient of body surface temperature of aggregated individuals (ΔT h = T h -Ta; 15.2 ± 3.4°C) decreased with aggregation size (estimate ± SD = -0.006 ± 0.003, z = -2.01, P = 0.045). In comparison, no significant relationship was found between gradient and aggregation size in grass (estimate ± SD = -0.003 ± 0.02, z = -0.19, P = 0.85) and beach (estimate ± SD = -0.002 ± 0.007, z = -0.31, P = 0.75) habitats.

In grass habitat, gradient of surface flipper temperature of aggregated elephant seals (15.8 ± 4.3°C) observed during transects scans was lower than gradient of isolated seals (17.2 ± 4.5°C; estimate ± SD = -0.08 ± 0.03, z = -2.40, P = 0.02; N = 376). But no difference of temperature was revealed between aggregated and isolated animals in beach habitat (aggregated: 15.5 ± 4.3°C; isolated: 16.2 ± 3.8°C; estimate ± SD = -0.003 ± 0.03, z = -0.10, P = 0.92; N = 419) or wallows (aggregated: 17.6 ± 4.1°C; isolated: 17.7 ± 3.7°C; estimate ± SD = 0.01 ± 0.05, z = 0.24, P = 0.81; N = 189).

No significant difference was found in gradient of flipper surface temperature between central and peripheral seals observed in aggregations (grass: central: 15.6 ± 3.6°C, peripheral:

15.7 ± 4.4°C, estimate ± SD = 0.03 ± 0.08, z = 0.36, P = 0.72, N = 103; beach: central: 16.6 ± 3.4°C, peripheral: 15.3 ± 4.5°C, estimate ± SD = -0.06 ± 0.07, z = -0.90, P = 0.37, N = 134; wallows: central: 17.6 ± 4.3°C, peripheral: 17.7 ± 4.1°C, estimate ± SD = 0.002 ± 0.04, z = 0.05, P = 0.96, N = 156).

Stomach temperatures averaged 36.6 ± 0.5°C (range 30.3-39.9°C). We did not observe major or constant variations in average or extreme values over time before the signal was lost.

Stomach temperatures recorded for individuals observed in aggregation were lower than

A C C E P T E D M A N U S C R I P T 25
stomach temperatures recorded for isolated individuals (N = 12, W = 3, P = 0.02). Both aggregated and isolated individuals were observed in wallows or grass habitat.

3.4.Aggregation, body mass and molt

Females at first capture weighed 316.2 ± 34.4 kg (range 259-410 kg) and their body mass loss during the molt averaged 3.3 ± 0.9 kg.d -1 (range 1.7-6.7 kg.d -1 ). Mean relative individual aggregation rate (number of observations in aggregation / total number of observations) was 0.69 ± 0.21 and molt rate averaged 10.0 ± 3.7 %.d -1 (range 2.0-17.5 %.d -1 ).

We found no significant correlation between individual aggregation rate and initial body mass (N = 35, r = -0.17, S = 8342.7.1, P = 0.33) or body mass loss per day (N = 33, r = 0.11, S = 5356.6, P = 0.56) or with individual molt rate (N = 41, r = 0.29, S = 5039.2, P = 0.09).

Discussion and conclusion

We found an increase in surface body temperature in molting elephant seals during the mid-stage of molt in agreement with previous data on harbor seals [START_REF] Paterson | Seals like it hot: Changes in surface temperature of harbour seals (Phoca vitulina) from late pregnancy to molt[END_REF].

Thus, thermal constraints of elephant seals depend on molt stage and heat loss is greatest during the mid-stage of molt. Moreover, our results showed that aggregation behavior depends on habitat type. Indeed, the proportion of aggregated seals was greatest and seals formed larger aggregations in wallows compared to beach and grass habitats. [START_REF] Laws | The elephant seal (Mirounga leonina Linn.) II. General, social and reproductive behavior. Falkland Islands Dependencies Survey[END_REF] and [START_REF] Boyd | Water flux, body composition, and metabolic rate during molt in female southern elephant seals (Mirounga leonina)[END_REF] already described that elephant seals prefer to aggregate in muddy wallows at South Georgia. This aggregation behavior in mud is not colony-specific, but common for this species, as long as the substrate of the site allows wallow formation. We the weight of the animal and movement destroying the vegetation after a period of time. The grass does not grow back where wallows are formed and the same wallows are used year after year (the same wallows were observed in the quadrat since 2012). Empty wallows dry off but every season molting elephant seals move hundreds of meters from the shore to join a wallow, preferably choosing one already occupied [START_REF] Chaise | Local weather and body condition influence habitat use and movements on land of molting female southern elephant seals (Mirounga leonina)[END_REF]. Mud is likely to contain a mixture of soil, feces, urine and replenished by rainfall and may be influenced by heat generated from seals.

Elephant seals appear to aggregate in wallows either in search of a warmer environment or increase the temperature of the environment through local heating from metabolic heat production. Indeed, this study showed that ground temperature is higher in wallows compared to other habitats, so seals would lose less heat by conduction with the substrate in wallows compared to grass and beach habitats. Moreover, wallows are less exposed to wind than grass and beach habitats, reducing heat loss by forced convection. Our previous work found that wallows were preferably selected at initial and mid-stages of the molt [START_REF] Chaise | Local weather and body condition influence habitat use and movements on land of molting female southern elephant seals (Mirounga leonina)[END_REF] when seals are shedding their old hair and skin but new hair has not yet grown. Thus, molting females experience a warmer habitat (i.e. wallows) when aggregating when heat loss is greatest [START_REF] Paterson | Seals like it hot: Changes in surface temperature of harbour seals (Phoca vitulina) from late pregnancy to molt[END_REF]. The loss of old skin and hair and increased peripheral blood flow for cell growth are two mechanisms that increase heat loss during molt [START_REF] Paterson | Seals like it hot: Changes in surface temperature of harbour seals (Phoca vitulina) from late pregnancy to molt[END_REF]. The role of fur as insulation in seals compared to blubber, is questionable but may have some value in air, mainly for pups and less for adults, but negligible in water [START_REF] Kvadsheim | Thermal function of phocid seal fur[END_REF][START_REF] Paterson | Seals like it hot: Changes in surface temperature of harbour seals (Phoca vitulina) from late pregnancy to molt[END_REF]. The fact that elephant seals spend most of their lifetime at sea has resulted in the evolution of a divingspecialized hair structure or to provide mechanical protection of the skin surface, suggesting that fur is not an important part of their insulation [START_REF] Ling | The skin and hair of the southern elephant seal, Mirounga leonina (L.): III. Morphology of the adult integument[END_REF][START_REF] Ling | Pelage and Molting in Wild Mammals with Special Reference to Aquatic Forms[END_REF]. However, molt stages air until new hair growth), and thus correspond to the visible part of the molting process [START_REF] Boyd | Water flux, body composition, and metabolic rate during molt in female southern elephant seals (Mirounga leonina)[END_REF]. The fact that aggregation behavior is less developed in grass habitat could be related to the use of grass as a transition habitat between beach and wallows, depending on molt stage and weather [START_REF] Chaise | Local weather and body condition influence habitat use and movements on land of molting female southern elephant seals (Mirounga leonina)[END_REF]. The variation of aggregation behavior between habitat types could be linked to the use of habitat depending on molt stage and to the sensitivity to weather conditions (i.e. heat loss) during the molting process [START_REF] Chaise | Local weather and body condition influence habitat use and movements on land of molting female southern elephant seals (Mirounga leonina)[END_REF].

Wallows appear to be a specific and favorable habitat to molt and aggregate and aggregation behavior of molting females in this habitat was influenced by weather. Indeed, in accordance with observations of [START_REF] Liwanag | Thermal benefits of aggregation in a large marine endotherm : huddling in California sea lions[END_REF] in hauled-out California sea lions, our results showed that the proportion of aggregated elephant seals (in all three habitats) and the aggregation size (in wallows) increase when weather is deteriorating. We also observed that elephant seals aggregate more during the mid-stage of molt, when cost of thermoregulation is greatest. Thus, aggregation behavior in female elephant seals appears to correspond to social thermoregulation (i.e. increase of huddling during unfavorable climatic conditions; [START_REF] Gilbert | One for all and all for one: The energetic benefits of huddling in endotherms[END_REF]. Molting in wallows may also facilitate shedding of skin through physical contact/abrasion with other seals and may be an additional benefit of wallow habitats for elephant seals.

Flipper surface temperature (correlated with body surface temperature), and stomach temperature, were greater in isolated individuals compared to aggregated seals. These results suggest that huddling while molting decreases thermoregulatory costs by reducing heat loss through warming of surrounding microclimate and by decreasing body surface area exposed to air, allowing females to reduce their internal thermal set-point and corresponding metabolic heat production [START_REF] Gilbert | One for all and all for one: The energetic benefits of huddling in endotherms[END_REF]). However, we did not find body surface temperature

A C C E P T E D M A N U S C R I P T

differences between peripheral and central individuals (central individuals being less exposed to heat loss). Identification of positions in the aggregation may not have been precise enough to observe an effect of reduced body surface area, or that the thermal benefits of being in the center of an aggregation are related to the aggregation's geometry [START_REF] Gilbert | One for all and all for one: The energetic benefits of huddling in endotherms[END_REF].

Indeed, we found that the number of central positions compared to peripheral ones (C/P ratio) increased with aggregation size, and that surface temperature of aggregated seals decreased with aggregation size in wallows, where mean aggregation size was greater than beach or grass habitats. Moreover, in grass habitat, we showed that C/P ratio of aggregations increases when weather is deteriorating, independently of an increase in aggregation size. Based on these results, it would be interesting to study in greater detail the dynamics of huddling in wallows to determine if elephant seals try to reach the center depending on their molt stage, body condition, dominance or variation in local weather (Cruwys and Davis, 1995;[START_REF] Gilbert | One for all and all for one: The energetic benefits of huddling in endotherms[END_REF]. We observed that there were more seals in a central position than in the periphery in grass and beach habitats compared to wallows and that this difference increased with aggregation size between wallows and grass habitat. This could be related to the fact that on non-wallow habitats (i.e. open area of habitats), seals gather in long linear aggregations resulting in several animals located in the center. In contrast in wallows, seals follow the edges of the muddy hollow, creating a more ellipsoidal huddle.

We expected that aggregation behavior would lead to faster molt in female elephant seals. If elephant seals could renew their hair faster, they could then spend less time fasting on land and save energy (i.e. blubber reserves). However, we found no correlation between aggregation rate and rate of molt. A possible explanation for this may be related to the difficulty of assessing stage of molt and molt completion [START_REF] Ling | The skin and hair of the southern elephant seal, Mirounga leonina (Linn.). IV. Annual cycle of pelage follicle activity and molt[END_REF]. Molt involves not only the shedding of skin and hair but also the regrowth of new hair [START_REF] Ling | The skin and hair of the southern elephant seal, Mirounga leonina (Linn.). IV. Annual cycle of pelage follicle activity and molt[END_REF]. We observed only the first process and therefore may not have been able to fully assess differences in molt We did not find any correlation between individual aggregation rate and body mass loss during the molt (i.e. energy expenditure). The link between aggregation and energysaving could be a more complex process (multifactorial) than the model proposes in this study. Energy-saving from thermal benefits in aggregations could be balanced with other processes such as requirements for cell growth and synthesis of new skin/hair or may be costly in terms of energy expenditure when active in huddles. The individual aggregation rate defined in this study was based on observations usually once per day which may not have been representative of individual aggregation behavior. Aggregation behavior may not be only driven by perceived cooling of the body but also by social factors. Elephant seals show high fidelity to breeding and molting sites [START_REF] Laws | The elephant seal (Mirounga leonina Linn.) II. General, social and reproductive behavior. Falkland Islands Dependencies Survey[END_REF] and we do not yet understand how aggregation behavior may be influenced by relatedness and other social aspects of conspecifics. We also measured body mass loss as an indirect measure of total energy expenditure during the molt as elephant seals were supposed fasting while on land [START_REF] Codde | The effect of local weather conditions on the behavior of molting southern elephant seals, Mirounga leonina (L)[END_REF]. However, some studies have questioned the possibility of fast-breaking during the molt [START_REF] Boyd | Water flux, body composition, and metabolic rate during molt in female southern elephant seals (Mirounga leonina)[END_REF][START_REF] Chaise | Local weather and body condition influence habitat use and movements on land of molting female southern elephant seals (Mirounga leonina)[END_REF]. Other heart rate measurements, accelerometry or doubly labelled water may therefore further complement body mass loss for estimating energy use related to aggregation behavior during the molt.

In conclusion, aggregation behavior in female elephant seals was found to be influenced by environmental factors resulting mainly in differences in surface body temperature which corresponds to previously reported behavioral and physiological aspects of social thermoregulation. The thermal advantages from aggregation mainly occur in wallows that are selected as a specific social habitat where seals can benefit from huddling at the peak of molt, when thermoregulatory costs are greatest. However, it is not clear if female southern elephant seals express huddling as an individual strategy in order to save energy during the molt. The relationship between behavioral and physiological thermoregulatory mechanisms seems more complex as body temperature and heat loss of molting elephant seals change during the molting process in association with their behavior on land. Future studies estimating the metabolic rate of animals may reveal further insights into energy savings associated with aggregation behavior during the molt, which is a relatively understudied phase of the elephant seal lifecycle.

-

  Abbreviations

  aggregation (number of central seals/number of peripheral seals). Thermal (± 0.1°C) and visual images were taken (ThermaCAM® P25, FLIR Systems, accuracy ±2°C;Lumix DMC-FS35 EF-K, Panasonic; TG-4, Olympus) for each aggregation or isolated female.

Figure 1

 1 Figure 1 Molting female southern elephant seals shedding old skin and hair: initial stage (A),

Figure 2

 2 Figure2Analyses of thermal images from observations of tracked individuals (A) and daily

  relative humidity, solar radiation and wind speed) on aggregation behavior. GLMMs were also used to analyze the effect of aggregation behavior (aggregated vs. isolated, and central vs. peripheral individuals) on the gradient of flipper surface temperature (ΔT f ), the effect of A C C E P T E D M A N U S C R I P T aggregation size on surface temperature of aggregated individuals (ΔT h ), and the effect of stages of molt on body surface temperature (ΔT b ). Analyses of aggregation effect on surface temperatures included meteorological variables (relative humidity, solar radiation and wind speed) as covariates, and date as a random effect. Aggregation identity number was added as a random effect for ΔT f analyses, and aggregation size was added as a covariate for comparison of ΔT f between central and peripheral individuals. Models were fitted with a Poisson distribution and final GLMMs were selected based on Akaike Information Criterion (AIC) for removal of non-significant effects and interactions. We used Wilcoxon tests to compare stomach temperature (T stom ) between aggregated or isolated, and central or peripheral elephant seals. Local meteorological variables between habitats were compared using Kruskal-Wallis rank sum tests followed by a multiple comparison test when significant (adjusted pairwise comparisons; kruskalmc in pgirmess package). We used a Kolmogorov-Smirnov test to compare distribution of molt stages between aggregated and isolated elephant seals during transects and quadrat scans. Spearman's rank coefficient tests were used to analyze correlations between body condition variations (initial body mass and daily body mass loss), molt rate, and aggregation behavior (individual relative aggregation rate). Results were expressed as mean ± standard deviation (SD). All statistical analyses were performed with R statistical software (R Development Core Team, version 3.5.1; RStudio Inc., version 1.1.456) and statistical significance was accepted at P < 0.05. 3. Results 3.1.Thermal environment of molting elephant seals: habitats and stage of molt When we compared meteorological variables between transects and quadrat in 2012 and 2014-15, we observed that the difference between ground temperature and air temperature (T g -T a ) was greater in wallows compared to beach and grass habitats, while no difference ACCEPTED MANUSCRIPT A C C E P T E D M A N U S C R I P T 18

Table 2

 2 Local meteorological variables (mean ± SD) recorded at the start of daily transects and quadrat scans(2012 and 2014-15) for each habitat type and significant differences (posthoc test: P < 0.05). We recorded a mean of 5 ± 2 rainy days per year during scans over the field session.We compared body surface temperature between the three stages of molt for individual(2012 and 2014, N = 134) and population surveys data (2012 and 2014-15, N = A C C E P T E D M A N U S C R I P T 19 708

Figure 3

 3 Figure 3 Body surface temperature (ΔT b ) gradient between body surface temperature (T b ) and

  of molt amongst aggregated seals compared to isolated seals, and mainly seals at initial or final stage of molt amongst isolated seals (Fig.4).

Figure 4

 4 Figure 4 Distribution of molt stages amongst aggregated seals (A) and isolated seals (B)

Figure 5

 5 Figure 5 Variation in aggregation rate (number of aggregated seals / total number of seals; A)

  Suzanne that wallows are created by seal aggregations within grass habitat,

  study are based on the percentage of old hair shed (i.e. bare-skin exposed to

  individuals. Measurement of the rate of growth of new hair may therefore provide more precise estimates of molt completion in relation to aggregation behavior of females.

S1

  Figure1

Figure 5

 5 Figure5

Years Number of transects scans (grass/beach) Number of quadrat scans (wallows) Number of tracked/recaptured [equipped] females Number of observations Tracking duration (days)

	2. Material and methods 2015 2.1.2. Population surveys
	2.1.Data collection 22 Jan to These data are based on daily observations of unidentified female elephant seals in
	2.1.1. Study site 27 Feb / defined areas. Two strip transects on beach habitat (49°26'02''S, 70°26'23''E -49°25'59''S, / 12/13 [4] 4.5 ± 1.6 6.1 ± 1.9
	Field work took place during four molting seasons of adult female elephant seals in 70°26'17''E; 150 m long; ± 10 m from the line transect) and grass habitat (49°26'00''S, 2016
	2012 and 2014-16 (Table 1), at the colony of Pointe Suzanne (49°26'S, 70°26'E) in Table 1 Dates of transects and quadrat scans, number of female elephant seals recaptured, 70°26'16''E -49°26'03''S, 70°26'22''E; 150 m long; ± 10 m from the line transect) and one
	M A N U S C R I P T Kerguelen Island (French Southern and Antarctic Lands). Based on the assumption that different areas would offer different thermal environments, the study site was divided into three different habitats based on substrate type and topography: stony beach, grassland and wallows (mud pools without vegetation created by aggregation of molting elephant seals within grassland). M A N U S C R I P T tracked and equipped with stomach temperature pills, mean number of observations and tracking duration (days). quadrat in wallow habitat (49°26'16''S, 70°25'59''E -49°26'20''S, 70°25'45''E -49°26'32''S, 70°25'46''E -49°26'26''S, 70°26'09''E; 0.14 km²) were defined to study the influence of habitat type and local meteorological parameters on elephant seal aggregation behavior during the molt. The transects and quadrat were scanned daily for a total of 220 scans between 2012 and 2014-15 (Table 1). At the start of each scan, meteorological variables, including air temperature (T a , °C), ground temperature (T g , °C; 5 cm depth; from 2014), relative humidity (%), wind speed (m.s -1 ) and solar radiation (W.m -²), were measured M A N U S C R I P T using hand-held devices (Kestrel 3000 Pocket Weather Meter; pyranometer SKS111, Skye
	A C C E P T E D 30/30 14 39/40 33 A C C E P T E D Instruments Ltd, Llandrindod Wells, UK). During each scan, we recorded the molt stage of 18 Jan to 19 Feb 2012 15/12 [0] 5.7 ± 3.3 13.7 ± 3.8 29 Dec 2013 to 1 Mar 25/21 [3] 5.2 ± 2.9 8.6 ± 3.3 A C C E P T E D each seal observed (mainly adult females, with possible presence of few juveniles and males).
	2014
	23 Dec
	2014 to	11/11	12	7/7 [2]	1.3 ± 0.5	4.1 ± 1.0
	15 Jan
						8 9 10

Conflict of Interest: The authors declare that they have no conflict of interest.

Ethical approval

All applicable international, national and institutional guidelines for the care and use of animals were followed. All procedures performed in the study were in accordance with the ethical standards of the institution or practice at which the study was conducted. Highlights :

-Wallows (mud pools) offer a warm and sheltered habitat.

-Heat loss of female elephant seals increases at the peak of the molt.

-A higher proportion of female elephant seals at the peak of their molt was observed amongst aggregated compared to isolated seals.

-Proportion of aggregated female elephant seals, and number of females per aggregation, are greater in wallows compared to open habitat types, and increase when weather is deteriorating.

-Gradient of surface temperature and stomach temperature are cooler in aggregated compared to isolated female elephant seals.