Monica Daniela Iordache

Claudiu Badulescu
email: claudiu.badulescu@ensta-bretagne.fr

Malick Diakhate

Marius Adrian Constantin

Eduard Laurentiu Nitu

Younes Demmouche

Matthieu Dhondt

Denis Negrea

A numerical strategy to identify the FSW process optimal parameters of a butt-welded joint of quasi-pure copper plates

Keywords: Friction stir welding, Optimal process parameters, Cu-DHP copper, Finite element method, Mass scaling, Digital images correlation

ou non, émanant des établissements d'enseignement et de recherche français ou étrangers, des laboratoires publics ou privés.

Introduction

Copper and its alloys are widely used in many engineering applications since they have particular properties such as electrical and thermal conductivity, ductility, mechanical strength, corrosion resistance, etc. However, these materials are difficult to weld using conventional welding processes (fusion melting) because of the high thermal conductivity and the higher oxidation rate at temperature values close to fusion one [START_REF] Hwang | Experimental study on Friction Stir Welding of copper metals[END_REF]. Thus, an alternative and promising solution is friction stir welding (FSW). This joining process, environmentally friendly since it releases no toxic gas or radiations, is solid-state welding based on the generation of heat due to both the friction of the pin (rotating element) at the interface of the two plates to be joined and the plastic deformation of the material being welded. FSW is a metal-working process commonly characterized by a flow progressing through preheating, initial deformation, extrusion, forging, and cool-down processing zones. The rotating FSW pin plunges into the junction of two rigidly clamped plates until the shoulder touches the surface of the materials being welded before traversing along the weld line under the applied normal force [2,[START_REF] Galvão | Critical review on friction stir welding of aluminium to copper[END_REF].

Many research studies have focused on the use of FSW to joint aluminum-based alloys, such as [START_REF] Cole | Weld temperature effects during friction stir welding of dissimilar aluminum alloys 6061-t6 and 7075-t6[END_REF][START_REF] Dourandish | Microstructure, mechanical properties and failure behaviour of protrusion friction stir spot welded 2024 aluminium alloy sheets[END_REF][START_REF] Jacquin | A review of microstructural changes occurring during FSW in aluminium alloys and their modelling[END_REF][START_REF] Zuo | Thermal fields in dissimilar 7055 Al and 2197 Al-Li alloy FSW T-joints: numerical simulation and experimental verification[END_REF], while limited literature on dissimilar FSW of copper is available. The limited research on FSW, and consequently its industrial application, is due to both the high melting point and the good thermal conductivity of the cooper as they lead to increasing the heat input during the welding process to obtain a joint without defects. This also results in higher requirements of materials and design of the welding tools. The effects of pin rotational speed on the microstructure, the mechanical properties, and the fracture localization of copper butt-joints have been investigated [START_REF] Zhou | Effects of rotation speed of assisted shoulder on microstructure and mechanical properties of 6061-T6 aluminum alloy by dual-rotation friction stir welding[END_REF]. The effects of both the shoulder cavity and the welding parameters on the applied torque, the formation of weld defects, and the mechanical properties of deoxidized copper (thin films) butt-joints have been studied [START_REF] Leal | Effect of shoulder cavity and welding parameters on friction stir welding of thin copper sheets[END_REF]. From the industrial point of view, the major challenge today is to able to predict the optimal welding parameters (advancing and rotational speeds of the pin), and many researchers have experimentally investigated this challenge [START_REF] Ramachandran | Performance analysis of dissimilar friction stir welded aluminium alloy AA5052 and HSLA steel butt joints using response surface method[END_REF][START_REF] Shashi Kumar | Identifying the optimal FSW process parameters for maximizing the tensile strength of friction stir welded AISI 316L butt joints[END_REF]. Most recently, some studies [START_REF] Zhang | The optimal temperature for enhanced lowtemperature superplasticity in fine-grained Ti-15V-3Cr-3Sn-3Al alloy fabricated by friction stir processing[END_REF] have shown an optimal welding temperature at which the joint is free from defects, the homogeneity of the grain size within the joint is close to that of the base material, and the mechanical properties of the joint are optimum. This optimal temperature, which is associated with a welding configuration (material and geometry of the pin, the relative positioning of the parts to be welded) is often experimentally evaluated, thus leading to costly and long tests campaign.

Furthermore, this optimal welding temperature seems to be an intrinsic parameter of the material, and the torque, the geometry of the pin, and the kinematics of the process help to achieve this temperature. Another recent study [START_REF] Heidarzadeh | Development of a fuzzy logic based model to elucidate the effect of FSW parameters on the ultimate tensile strength and elongation of pure copper joints[END_REF] used, for the first time, a fuzzy logic model to investigate and optimize the effect of FSW parameters on the tensile properties of the copper joints, such as ultimate tensile strength (UTS) and elongation.

The microstructure of the optimum joint was characterized using different microscopy techniques to elucidate the origin of its behavior.

Several researchers have investigated the relationship between the microstructure of the welded assembly and the resulting mechanical properties. Findings from these studies show that the microstructural evolution of pure copper joint during the FSW process is mainly influenced by welding temperature. At temperature values lower than around 0.5 times Tf (Tf being the material fusion temperature) the microstructure is mainly governed by continuous recrystallization. The grain structure development at welding temperature values greater than around 0.5 times Tf is mainly dominated by discontinuous recrystallization, thus resulting in a relatively coarse grain structure in the stir zone [START_REF] Mironov | Microstructural evolution of pure copper during friction-stir welding[END_REF].

Dynamic and static changes of grain size and texture of copper during the FSW process have been studied by Xu et al. [36]. They have shown that during the plastic deformation stage, the initial coarse grains in the base metal were subdivided with the increase in strain and temperature, and the stir zone exhibits ultrafine grains with a large number of low angle boundaries and a symmetrical simple shear texture. The effect of stacking fault energy on the restoration mechanisms and the mechanical properties of joints has been studied by Heidarzadeh et al. [37]. Authors have shown that in pure copper the continuous dynamic recrystallization was the only restoration mechanism that leads to new grains formation. Moreover, the effect of restoration mechanisms on both the yield strength and the strain hardening of joints has been studied. For a deeper understanding, one can read the review paper recently published [38] which highlights the latest progress made in the microstructure analysis during the FSW process.

Based on the literature review [START_REF] Ramachandran | Performance analysis of dissimilar friction stir welded aluminium alloy AA5052 and HSLA steel butt joints using response surface method[END_REF][START_REF] Shashi Kumar | Identifying the optimal FSW process parameters for maximizing the tensile strength of friction stir welded AISI 316L butt joints[END_REF], one can identify the optimal parameters of the process through numerical simulations, thus leading to a drastic reduction of the experimental campaign. This finding highlights the importance of developing modeling strategies as well as efficient and robust numerical simulations. Today, the modeling of heat generation and material flow during the FSW process are addressed using different modeling techniques such as computational solid mechanics (CSM) and computational fluid dynamics (CFD) methods. The current literature proposes some approaches useful to tackle the mesh distortion. With this in mind, Al-Badour et al. [START_REF] Al-Badour | Coupled Eulerian Lagrangian finite element modeling of friction stir welding processes[END_REF] proposed the coupled Eulerian-Lagrangian (CEL) model. In this method, the plates being welded, governing the equations, are discretized using the Eulerian formulation whereas the welding tool is modeled as a Lagrangian body. These methods are based on explicit integration schemes which imply an extremely reduced step time, a high number of increments, and a high computational time. This is why most of the FSW numerical simulations are restricted to short computational times (penetration phase of the pin).

However, there are very few studies that propose robust numerical strategies and help to estimate the optimal process parameters without having to conduct extremely timeconsuming and costly experimental parametric studies.

The current work aims to develop an efficient strategy to determine the welding process parameters based on numerical simulations with the finite element method and the optimal temperature that is supposedly known for all materials that will be welded. The numerical model developed in this work should be able to simulate within a reasonable computational time (one to two days) an FSW process characterized by a welding time of about a few tens of seconds and a traversing length of the pin of about 100 mm. The mass scaling technique is used to minimize the computational time. Also, the numerical model has enabled us to take into account the thermo-mechanical boundary conditions (heat exchange with the ambient air), the heat dissipation from plastic deformation, and that from the friction between the pin and the materials being welded. The friction coefficient was supposed as a function of the temperature. This strategy has been used to simulate the butt welded joint of 3-mm thick quasi-pure copper and has been validated by comparison with experimental results (identification of defects using x-ray technique, local mechanical behavior of the joint analyzed using digital images correlation, and fracture surface analysis through fractography technique).

Finite element model description

Finite element approach

Numerical simulations have been conducted to assess the effect of different FSW process parameters on the quality of the FSW joints. Modeling the FSW process is a challenging task since it needs to address a multi-physical problem that is characterized by large plastic deformation and high temperature. The more difficult to suitably model is the material flow in the vicinity of the pin. Nowadays, several numerical simulation strategies are proposed to model the FSW process. In the non-flow-based models, excessive deformations appear and thus leads to early termination of the computation. The ALE formulation is often used to ensure a better mesh quality during the simulation [START_REF] Bussetta | Comparison of a fluid and a solid approach for the numerical simulation of friction stir welding with a non-cylindrical pin[END_REF].

However, because of extremely large deformations, this strategy cannot eliminate the mesh distortion and consequently leads to prohibitive computation time [START_REF] Dialami | Material flow visualization in friction stir welding via particle tracing[END_REF]. Recently, coupled Eulerian-Lagrangian (CEL) models have been developed to simulate severe plastic deformations [START_REF] Chauhan | Modeling of defects in friction stir welding using coupled Eulerian and Lagrangian method[END_REF]. CEL models, which were originally used to describe the thermo-mechanical response (forces and temperature distribution), now allows predicting defects (tunnel defects, cavities, and excess flash formation) that may occur [START_REF] Chauhan | Modeling of defects in friction stir welding using coupled Eulerian and Lagrangian method[END_REF].

In this work, a 3D thermomechanical finite element model based on the CEL formulation has been developed using ABAQUS software. The direct integration scheme was adopted, and the resolution was explicit. This integration scheme is well suited for problems that include numerous nonlinearities, dynamic phenomena, and thermal effects.

The finite element model is described below.

Geometry of the model and applied boundary conditions

The two workpieces to be welded are 100 mm in length, 40 mm in width, and 3 mm in thickness. As shown in Fig. 1-a, these workpieces are positioned into the Eulerian domain. This latter has a thickness of 3,5 mm. The geometry of the tool has a conical shape with diameters of 4 mm at its base and 3 mm at its upper part and a length of 2.8 mm (Fig. 1a). The tool is considered to be rigid and the workpieces deformable. All the lateral and bottom faces of the Eulerian domain are embedded. The tool is characterized by three movements. The first one is the rotation at a constant speed (𝑣 𝑟 ) coupled with a penetration speed (𝑣 𝑝 ) of about 2.8 mm/s for all the simulations. Once the penetration step is completed, the tool will begin the actual welding phase over a distance of 90 mm in the x-direction (Fig. 1-a). During this latter, the tool keeps its rotational speed constant while at the same time it is subjected to an advancing speed. Several sets of parameters, which are characterized by different pairs of welding speed (𝑣 𝑟 , 𝑣 𝑎 ) values are simulated.

To save computational time, the tool retraction phase was not modeled. A constant ambient temperature field (22°C) is applied to the whole Eulerian domain and at the beginning of the simulation. The heat exchange during the welding process has been taken into account by using a convection model that is characterized by a convection coefficient (ℎ 𝑓 ) of 30 𝑊 / m 2 . °𝐶 [START_REF] Constantin | An efficient strategy for 3D numerical simulation of friction stir welding process of pure copper plates[END_REF]. The heat exchange through the tool is not taken into account firstly to limit the size of the numerical model and secondly because of the low heat exchange that should be simulated since the thermal conduction coefficient of the copper remains relatively greater than the one of the tool. This step allows defining the boundary conditions at the surface to be welded.

Mesh

The Eulerian domain mesh comprises 6785 thermally coupled 8-node Eulerian elements (EC3D8RT) and 8364 nodes as shown in Fig. 1-b. The mesh was progressively performed by increasing the element width from the joint line to the edge. The tool is a Lagrangian rigid body and was meshed using 785 quadratic tetrahedral elements. The mesh has been adapted based on the following two main reasons. The first one is to ensure that each part of the tool which comes into contact with the plates being welded has the same contact surface by ensuring a relatively homogeneous of both friction and local distribution of the contact pressure, with effects on the temperature field distribution. The second main reason is to ease the modeling of complex geometries such as the tool's one. The number of elements has been chosen in the view of obtaining a good compromise between the computation time and the convergence of the field temperature.

Thermomechanical behavior of the material

In this work, the material being welded is phosphorus deoxidized copper (DHP-Copper) with 99.9% pure copper. This material is commonly used in industry due to its good corrosion resistance, high thermal conductivity, and electrical conductivity [START_REF] Galvão | Influence of tool shoulder geometry on properties of friction stir welds in thin copper sheets[END_REF]. The linear thermo-elastic behavior of the used material is described using the following constitutive equation [START_REF] Lau | Thermo-elastic behavior of a polymeric layer bonded between rigid interfaces[END_REF] :

𝜎 = 2𝜇𝜀 𝑒 + 𝜆[𝑡𝑟(𝜀 𝑒 ) -𝛼(𝑇 -𝑇 0 )]𝑰 (Eq.1)
Where 𝝈 is the stress tensor; 𝜺𝒆 is the linear-elastic strain tensor; 𝝀 and 𝝁 are the Lamé coefficients; T is the real temperature, 𝑻0 is the reference temperature, 𝜶 is the thermal expansion coefficient, and I is the identity matrix.

The thermal field equation is expressed as follows [START_REF] Hammelmüller | Increasing numerical efficiency in coupled eulerianlagrangian metal forming simulations[END_REF]:

-𝑘∇ 2 𝑇 = 𝛼𝜆𝑇 0 𝑡𝑟(𝜺̇𝒆) + 𝜌𝑐 𝑒 𝑇 ̇ (Eq.2)
Where 𝑘 is the thermal conductivity, 𝜌 is the density, and 𝑐 𝑒 is the specific heat.

The first term on the right-hand side of the previous equation (Eq.2) describes the strain rate effect on the temperature field. The mechanical behavior is formulated using Navier's equations for thermoelectricity and is given hereafter [START_REF] Hammelmüller | Increasing numerical efficiency in coupled eulerianlagrangian metal forming simulations[END_REF]:

𝜇𝛻 2 𝑢 + (𝜆 + 𝜇)𝛻𝑡𝑟(𝜺 𝒆 ) -𝛼𝜆∇𝑇 = 𝜌 𝜕 2 𝒖 𝜕𝑡 2 ̇ (Eq.3)
Where u is the displacement vector and t is the time.

The physical properties of the DHP copper at room temperature are given in Table 1.

In this study, the Johnson-Cook's model is used to describe the material flow. The model uncouples the plastic, viscous, and thermal behaviors and describes each of them through three independents terms [START_REF] Johnson | A constitutive model and data for metals subjected to large strains, high strain rates and high temperatures[END_REF].

𝜎 ̅ = [𝐴 + 𝐵 ⋅ (𝜀 𝑝 ̅ ) 𝑛 ] [1 + 𝐶 ⋅ ln ( 𝜀 ̅ ̇𝑝𝑙 𝜀̇0 )] [1 -( 𝑇-𝑇 𝑟𝑒𝑓 𝑇 𝑚𝑒𝑙𝑡 -𝑇 𝑟𝑒𝑓 ) 𝑚 ] (Eq.4)
Where 𝜎 ̅ is the flow stress; 𝜀 𝑝 ̅ the effective plastic strain; 𝜀p 𝑙 the effective plastic strain rate; 0  the normalizing strain rate; A, B, C, n, Tmelt, and m are material constants; Tref the room temperature (22 o C in this study).

In the previous equation (Eq.4), the parameter n takes into consideration the hardening of the material, whereas m depends on its fusion. C is influenced by the strain rate.

The material's constants were determined based on the experimental results from tensile tests that were carried out at different speeds (3 mm•min -1 and 30 mm•min -1 ), which are equivalent to strain rates of 1.67 × 10 -2 s -1 and 1.67 × 10 -3 s -1 and temperatures (22°C, 300°C, and 500 °C). During the FSW process, the applied strain rate can sometimes, depending on the welding configuration to perform, reach 1600 s -1 [START_REF] Arora | Strains and strain rates during friction stir welding[END_REF][START_REF] Andrade | Analysis of contact conditions and its influence on strain rate and temperature in friction stir welding[END_REF].

Nevertheless, such high-performance equipment is not available in our research entities.

The inverse identification method was used. A hydraulic testing machine (INSTRON 1342) equipped with a load cell (+/-100 kN capacity) was used to apply the loading rate and a climatic test chamber (CERHEC 1400), which can generate a controlled temperature for up to 1500°C, helps regulate the testing temperature. For the sake of clarity, only results from tests at 22°C and 500°C are presented in Fig. 2.

The identified values of the constants of Johnson-Cook's model for DHP copper are obtained by fitting the equation (Eq.4) with the experimental results and are shown in Table 2.

Friction model

A critical aspect of the FSW process simulation is the contact condition modeling between the tool and the plates being welded since the Eulerian domain interacts with the Lagrangian one. Many studies have focused on the development of contact models suitable for the FSW process. Most of them have opted for a friction coefficient that is kept constant during the simulation [START_REF] Chiumenti | Numerical modeling of friction stir welding processes[END_REF][START_REF] Chao | Heat transfer in friction stir welding -experimental and numerical studies[END_REF]. However, the friction coefficient depends on speed, temperature, and deformation rate. Recently, Kareem et al. [START_REF] Salloomi | Fully coupled thermomechanical simulation of friction stir welding of aluminum 6061-T6 alloy T-joint[END_REF] have used the Coulomb friction model with a non-linear coefficient that is dependent on both the local temperature of the melted material and the deformation rate. This evolution of the coefficient has been previously proposed by Meyghani et al. [START_REF] Meyghani | Developing a Finite Element Model for Thermal Analysis of Friction Stir Welding by Calculating Temperature Dependent Friction Coefficient[END_REF] through a highly original work that integrates the shear stress of the contact interface (dependent on the temperature), the partial sliding/ sticking condition, and the geometry of the tool.

Very promising results from experimental tests [START_REF] Meyghani | Developing a Finite Element Model for Thermal Analysis of Friction Stir Welding by Calculating Temperature Dependent Friction Coefficient[END_REF], which were carried out using various FSW parameter sets, have validated this evolution of the friction coefficient as a function of temperature. Thus, this methodology was used in this work to describe the relationship between the friction coefficient evolution and the temperature (Fig. 3). At ambient temperature, 𝜇 0 = 0.22 (evaluated by using the inverse identification method).

Mass scaling strategy

In the FEM model, an explicit integration scheme was used for the resolution. One of the major criticisms of this integration scheme is the extremely long computational time that is associated with it, so it is mainly used for dynamic simulations (simulation time relatively short). If the time increment is less than a critical value ∆𝑡 𝑐𝑟𝑖𝑡 , the integration scheme is considered as conditionally stable. The critical time increment ∆𝑡 𝑐𝑟𝑖𝑡 is computed from the mass and stiffness characteristics of the model and is expressed as follows [START_REF] Ducobu | On the introduction of adaptive mass scaling in a finite element model of Ti6Al4V orthogonal cutting[END_REF]:

∆𝑡 𝑐𝑟𝑖𝑡 = 𝑚𝑖𝑛 ( 𝐿 𝑐 𝑖 𝐶 𝑑 ) (Eq.5)
Where

𝐶 𝑑 = √ 𝐸 𝜌
is the wave propagation velocity within the material and 𝐿 𝑐 𝑖 is the characteristic length of each element 'i' of the mesh.

Mass scaling is a way of reducing the computational time, and it increases artificially the masses of the elements and can be applied even though there is rate dependency. The mass (Eq.3) is scaled by replacing the density term 𝜌 with the fictitious density 𝜌 * = 𝜅 𝑚 • 𝜌, with 𝜅 𝑚 > 0 [START_REF] Hammelmüller | Increasing numerical efficiency in coupled eulerianlagrangian metal forming simulations[END_REF]. The mass scaling factor 𝜅 𝑚 has to be chosen in such a way that the inertial forces, the right-hand side of the equation (Eq.3), remain small. The substitution of the density 𝜌 for a fictitious density 𝜌 * leads to a change in the thermal time constant (Eq.2). This effect can be compensated by introducing the fictitious specific heat 𝑐 𝑒 * = 𝑐 𝑒 𝜅 𝑚 -1 . Thus, we obtained the two following scaled thermo-elastic equations (Eq.6 and Eq.7). Mass inertia effects can be seen explicitly on the right-hand side of the equation (Eq.6) [START_REF] Hammelmüller | Increasing numerical efficiency in coupled eulerianlagrangian metal forming simulations[END_REF].

-𝑘∇ 2 𝑇 = 𝛼𝜆𝑇 0 𝑡𝑟(𝜺 𝒆 ̇) + 𝜌 * 𝑐 𝑒 * 𝑇 ̇ (Eq.6)

𝜇∇ 2 𝒖 + (𝜆 + 𝜇)∇𝑡𝑟(𝜺 𝒆 ) -𝛼𝜆∇𝑇 = 𝜌 * 𝜕 2 𝒖 𝜕𝑡 2 (Eq.7)
To achieve a reasonable accuracy of simulation results, the ratio of kinetic energy to internal one must be less than 2% of the simulated model. According to [START_REF] Constantin | An efficient strategy for 3D numerical simulation of friction stir welding process of pure copper plates[END_REF], a value of 𝜅m= 1000 was chosen. Thus, the temperature error is less than 10% and the computational time is reduced by 25 times. As a consequence, the use of the mass scaling method combined with a reasonable computational time leads to a significant reduction in the increments number as well as the inherent numerical errors.

Validation of the finite element model

To allow validation of the FEM, two DHP copper plates with dimensions of 100 mm (length) x 100 mm (width) x 3 mm (thickness) have been welded with two sets of welding parameters. The assembly process is performed by using a welding machine FSW-4-10 that is characterized by a rotating speed in the range of 300 to 1450 rpm and an advancing speed between 10 and 480 mm/min. Throughout the welding process, this machine also allows both controlling the displacement and recording the force in the 𝑧 ⃗ direction (Fig. 1b). The temperature has been measured using an infrared camera (FLIR A40M) with an accuracy of +/-2°C, and at the interface between the tool and the plates being welded (precisely at 1 mm behind the tool and pointing to the weld bead). The infrared camera moves with the tool.

The first welded assembly, which is labeled W-90-800, was obtained with an advancing speed of 90 mm/min and a rotating speed of 800 rpm. The second welded assembly, which is labeled W-90-1000, was obtained with an advancing speed of 90 mm/min and a rotating speed of 1000 rpm. The rotational speed value of the tool was modified while the advancing one was kept constant. This choice was motivated by the fact that the welding temperature increase is mainly influenced by the rotational speed than the advancing one and this temperature value could be less affected by measurement errors.

The evolution of maximal welding temperature is plotted against the position of the tool during the welding of W-90-1000 (Fig. 4). This first result helps to evaluate the initial value of the friction coefficient 𝜇 0 (Fig. 3) by minimizing the difference between values that are predicted by the model and those obtained from the experimental measurements.

Once the value of 𝜇 0 identified, the validity of the numerical model can be evaluated by comparing the temperature distribution measured within the second welded assembly during the tool advance with the one predicted by the finite element model (Fig. 5).

Moreover, the numerical axial force is compared with the experimental one (Fig. 5-b).

These findings indicate that the experimental results are in good agreement with those obtained from the numerical simulation, the force axial error is less than 6 %.

Parametric study

It has been proven [START_REF] Hwang | Experimental study on Friction Stir Welding of copper metals[END_REF] that it is crucial to reach the optimal welding temperature for obtaining a welded joint of high quality characterized by a mechanical strength close to that of the base material. This temperature can be experimentally evaluated and is about 0.4 to 0.5 times Tmelt for quasi-pure copper materials [START_REF] Hwang | Experimental study on Friction Stir Welding of copper metals[END_REF]. The optimal welding temperature of DHP copper materials [START_REF] Hwang | Experimental study on Friction Stir Welding of copper metals[END_REF][START_REF] Constantin | Experimental investigations of tungsten inert gas assisted friction stir welding of pure copper plates[END_REF] is about 550°C. At this temperature, the base material is in a pasty state, so this allows a homogeneous melting while avoiding defect formation. If this temperature is exceeded and approaching that of melting, the material becomes too fluid, which will result in both voids formation within the joint and inhomogeneous melting. This latter also leads to excessive burr of the weld bead and, under tensile loading, fracture at the welded joint.

This parametric study aims at carrying out numerous simulations which subsequently will help to identify the optimal welding parameters. Twelve simulations are carried out using the welding parameters specified in Table 3. Each simulation is individually labeled, indicating the speeds of both advancing and rotating.

For instance, Fig. 6-a shows the computation field temperature obtained during the simulation of the welding configuration S-90-1000 and precisely when the tool position is at 85 mm from point A. As shown in Fig. 5-a, both experimentally and numerically, the welding temperature takes time to reach its optimal value. In these welded areas where the temperature is not stabilized, defects likely have appeared. Consequently, for each numerical simulation, the temperature is recorded at a tool position greater than 85 mm from point A. These temperature values are shown in Fig. 6-b by the red points. For a clear presentation of these temperature results, a polynomial surface interpolation was performed. This highlights the effects of the FSW process parameters on the stabilized welding temperature. This finding is in good agreement with experimental observations [START_REF] Padhy | Friction stir based welding and processing technologies -processes, parameters, microstructures and applications: A review[END_REF]. At a given advancing speed 𝑣 𝑎 , the stabilized temperature increases with the rotating speed 𝑣 𝑟 . However, when this latter is fixed, a decrease in the advancing speed leads to an increase in the stabilized process temperature.

Simulation Results

From the obtained results one can compute a thermal efficiency surface indicator (𝑻 𝑬 ) that is defined in the following formula (Eq.8):

𝑻 𝑬 (𝑣 𝑎 , 𝑣 𝑟 ) = [1 -

𝑎𝑏𝑠[𝑻 𝑺 (𝑣 𝑎 ,𝑣 𝑟 )-𝑻 𝑶 ]

𝑻 𝑶 ] × 100 (Eq.8)

Where 𝑻 𝑺 (𝑣 𝑎 , 𝑣 𝑟 ) represents the stabilized temperature obtained from the simulation with the speeds 𝑣 𝑎 and 𝑣 𝑟 , and 𝑇 𝑂 represents the optimal welding temperature identified from the literature [START_REF] Constantin | Experimental investigations of tungsten inert gas assisted friction stir welding of pure copper plates[END_REF][START_REF] Hwang | Experimental study on Friction Stir Welding of copper metals[END_REF]. The maximal value of the indicator 𝑇 𝐸 might give a valuable indication of the optimal area of the welding process.

The thermal efficiency surface (𝑇 𝐸 ) mapping is shown in Fig. 7.

This mapping clearly shows an area, the one colored in red, that is characterized by values of 𝑇 𝐸 close to 100%, in other words, by stabilized temperature values close to the temperature 𝑇 𝑂 . This result means that, inside that particular area, one can identify FSW process parameters that help to obtain welded joints with mechanical strengths close to that of the base material. Ideally, it is recommended to perform the welding at relatively high speeds (𝑣 𝑎 = 120 𝑚𝑚 𝑚𝑖𝑛 ⁄ and 𝑣 𝑟 = 1150 𝑟𝑝𝑚) to avoid heat loss in the weld bead, but also to increase productivity. It should be noticed that these values could be limited due to both the kinematics and rigidity of the welding machine. Conversely, selecting low optimal values leads to the use of both moderate kinematics and less effective welding machines, on the condition that the heat losses did not become too high.

In this parametric study, only the heat exchange with the ambient air is integrated into the model. However, evaluating the effects of more realistic thermal boundary conditions remains possible.

From Fig. 7, the optimal dependency between the rotating speed and the advancing one might be expressed as follows (Eq.9).

𝑣 𝑟 = 𝑘 𝑇 × 𝑣 𝑎 + 𝑘 𝑟 (Eq.9)

With 𝑘 𝑇 = 5.52 𝑟𝑜𝑡 𝑚𝑚 ⁄ and 𝑘 𝑟 = 483 𝑟𝑜𝑡 𝑚𝑖𝑛 ⁄ . These two constants depend on both the tool geometry and the material being welded. Their values are identified from the equation of the black-colored line plotted in Fig. 7.

The equation (Eq.8) allows evaluating the optimal value of the tool rotating speed based on an imposed advancing speed value. Some experimental tests are then needed to both validate this strategy and evaluate the difference between the stabilized temperature and the optimal one (𝑇 𝑂 ), as well as the joint strength.

Experimental validation

Experimental setup

In addition to two welding configurations (labeled W-90-800 and W-90-1000) that were used for both identification and validation of the finite element model, five other FSW assemblies have been manufactured using the same material (Cu-DHP) and geometry (length, width, and thickness). The labels (process parameters) of the five other welding configurations are specified in Table 4. The non-destructive testing using penetrating radiation (X-ray radiography) was performed on ANDREX X-ray Equipment model CMA357 using high contrast, very fine grain KODAK INDUSTREX T200 Film. The focal distance was 50 cm, current intensity 2 mA, working voltage 120 kW, and the exposure time 1 min and was performed on four specimens (labeled W-90-1200, W-150-1200, W-90-800, and W-150-800). This investigation aims at qualitatively evaluating the presence or absence of defects having a dimension greater than 300µm within the joint.

For each FSW specimen, the investigated volume size is graphically represented by the blue-colored parallelepiped shown in Fig. 8 The experimental setup used to measure the displacement field on the above-mentioned faces is shown in Fig. 9 The testing machine used to investigate the mechanical behavior of the welded joint is the same as that for the thermo-mechanical characterization of the base material, namely INSTRON 1342. Monotonic tensile tests were performed at a loading rate of 3 mm/min.

The displacement field measurement on the two perpendicular faces (frontal and lateral) of the sample, was performed by simultaneously using two CCD cameras that are fully synchronized with the INSTRON machine: an Aramis-GOM system equipped with a 2448 x 2050 pixels CCD sensor and a Retiga 6000 equipped with a 2758 x 2208 pixels CCD sensor. The resulting strain field is computed by deriving the displacement field.

Before installing the specimen in the testing machine, its two perpendicular faces were speckled and on an area of interest which is 45 mm long and centered on the welded joint (Fig. 8-b). This length has been selected so that a wide area around the welded joint can be monitored and recorded up to failure. A polarized lighting device was also set up for recording images with high contrast at a frequency of 1 Hz.

The mechanical loading is applied along the 𝑦 ⃗-direction (Fig. 8). For image analysis, the spatial resolution is set to 19 pixels, and the standard deviation of the displacement field to +/-0.5 µm.

To complete this experimental investigation, a digital microscope (KEYENCE VHX) was used to analyze the fracture surfaces of three specimens (labeled W-90-1000, W-120-1000, and W-150-1000) at 100x magnification to highlight the failure scenario through different zones of FSW joints.

Results from X-ray analysis

The results from the X-ray analysis of the four specimens are shown in Fig. 10.

As shown in Fig. 10, false-color images of the above-mentioned faces (frontal and lateral) allow detection of the defects along both the axis (y-direction) and width (z-direction) of the joint.

From image analysis, it can be noticed that there is no detectable defect within the specimen labeled W-90-1000 (Fig. 10-a). Consequently, this FSW specimen should exhibit higher mechanical properties (maximum values of stress and strain at failure). The defects detected within the specimen labeled W-150-1000 (Fig. 10-b) seem to be a tunnel type defect, are characterized by a variation in its width along the welding direction and are located in the root of the weld (retreating side). The defects within the two other specimens labeled W-80-800 and W-150-800 (Fig. 10-c and Fig. 10-d, respectively) have a size bigger than that of those detected within the specimen labeled W-90-1000, extend over the whole width of the specimen, are located in the root of the weld (retreating side) and characterized by a 1 mm in depth. These specimens are likely to exhibit low mechanical strength.

Results from monotonic tensile tests

The macroscopic mechanical responses of the welding configurations that were tested under tensile load (Table 4) are shown in Fig. 11. In this figure, the true stress (𝜎 𝑣 ̅̅̅) is plotted as a function of the logarithmic strain (𝜀 𝑣 ̅ ).

Results from tensile tests (Fig. 11) highlight the effects of the welding parameters on the macroscopic mechanical behavior of the welded joint. For some process parameters, the values of both the logarithmic strain and true stress at the specimen failure are close to the ultimate strength of the base material. For a suitable analysis of the results, the mechanical efficiency 𝐸 𝑀 of the joint is computed by using the following formula:

𝐸 𝑀 (𝑣 𝑎 , 𝑣 𝑟 ) = 𝜎 ̅ 𝑣 𝑚𝑎𝑥 𝐹𝑆𝑊 (𝑣 𝑎 ,𝑣 𝑟 ) 𝜎 ̅ 𝑣 𝑚𝑎𝑥 𝐵𝑀 × 100 (Eq.10)
Where 𝜎 ̅ 𝑣 𝑚𝑎𝑥 𝐹𝑆𝑊 represents the maximum value of the true stress resulting from the tensile test on an assembly welded at the speeds (𝑣 𝑎 , 𝑣 𝑟 ) and 𝜎 ̅ 𝑣 𝑚𝑎𝑥 𝐵𝑀 represents the maximum value of the true stress of the base material. Based on the results shown in Fig. 11 and Fig. 2, 𝐸 𝑀 values are plotted on Fig. 12-a as a function of the two welding speeds (𝑣 𝑎 , 𝑣 𝑟 ). This result (Fig. 12-a) highlights an area, around the white-colored line, within which the selected welding speeds lead to a joint characterized by a mechanical efficiency value close to 90%. It is also clear from Fig. 12-a that one can notice a very good correlation between the mapping of the 𝐸 𝑀 values and that of the thermal efficiency surface (𝑇 𝐸 ) values. Moreover, Fig. 12-b shows a good agreement between the welding speeds (𝑣 𝑎 , 𝑣 𝑟 ) selected from the mapping of the 𝑇 𝐸 values and the resulting strain value at specimen failure (𝜀̅ 𝑣 𝑚𝑎𝑥 𝐹𝑆𝑊). These findings validate the strategy for identifying the optimal process parameters.

The maximum strain (𝜀̅ 𝑣 𝑚𝑎𝑥 𝐹𝑆𝑊) value is around 0.3 when the welding is performed at the optimal process speeds, and that strain value is close to the one obtained when the base material fails. The key parameter for manufacturing a joint of high quality is mainly governed by the welding temperature obtained in the vicinity of the tool. Optimizing this welding temperature prevents the generation of both tunnel and kissing bond defects.

Strain maps results

Differences between the macroscopic behavior of FSW joints (Fig. 11) could be better explained by analyzing the local behavior of the different zones of the joint. For this purpose, the strain fields (𝜀 𝑦𝑦 ) obtained on both faces (frontal and lateral) of the specimens W-90-1000 and W-150-1000 are shown in Fig. 13. Local strain maps (in the loading direction) resulting from seven different values of macroscopic stress (Fig. 13-a) are compared. As stated previously, the analysis of X-ray results allows concluding that the defects detected within the specimen W-150-1000 are of tunnel type. This explains why the yield strength of that specimen is lower than that of the specimen labeled W-90-1000. From Fig. 13-b, one can see on the strain maps of the lateral face (W-150-1000) that the defect starts appearing at a stress value of around 130 MPa (from point C, Fig. 13-a).

This localization reflects the tunnel defect propagation from the root of the joint to the welded area.

This propagation scenario is highlighted on strain maps labeled d and e. The sudden decrease in the macroscopic stress value of the specimen W-150-1000 is because of the sudden propagation of the defect to the opposite face.

When analyzing the strain field on each frontal face of the two specimens, one can see that the defect starts being detected from the stress level labeled d. These localizations are undoubtedly associated with defect propagation. From these strain maps, one can also notice that defects are located on the retreating side and their appearance leads to high strain values in their vicinity. Under the same loading rate level, the highest strain value is about 0.2 for the specimen W-150-1000 and 0.02 for the specimen W-90-1000.

To better understand the mechanisms that govern the fracture of joints, fracture surface analysis was performed on three specimens (labeled W-90-1000, W-120-1000, and W-150-1000) using a high-resolution digital microscope. The results are shown in Fig. 14.

From Fig. 14, one can notice a cross-section reduction of about 67.5%, 30.6%, and 15.9%

for the specimen W-90-1000, W-120-1000, and W-150-1000 respectively. The specimen W-90-1000 exhibits a ductile fracture. Periodic striations are observed on the fracture surface of specimens W-120-1000 and W-150-1000. Similar striations have already been observed [START_REF] Zhou | Influence of the kissing bond on the mechanical properties and fracture behaviour of AA5083-H112 friction stir welds[END_REF] on the fracture surface of an aluminum AA5083-H112 alloy FSW joint.

Zettler et al. [START_REF] Zettler | Material Deformation and Joint Formation in Friction Stir Welding[END_REF] have observed machining grooves that can be associated with the tool advancing along both the cross-section and welding direction of an aluminum AA6063-T6 alloy FSW joint.

In the present study, the striations at the fracture surface of W-120-1000 are localized in its lower side of the FSW joint and at a periodic interval. The localization of these striations in the lower half of the cross-section highlights the fact that the insufficient temperature at the root of the weld does not ensure melt quality, thus leading to the appearance of kissing bond defects. By comparing the specimens W-120-1000 and W-150-1000 (having the same rotating speed but different advancing speed), one can observe striations over the entire cross-section of the specimen W-150-1000 and are distributed at low density compared to the specimen W-120-1000. This pattern was more pronounced at the root of the weld.

Microstructure, microhardness, and fracture surfaces

Microstructure and microhardness of joints have been analyzed and compared to those from specimens obtained under extreme welding conditions, i.e. based on thermal and mechanical efficiency surface indicators 𝑇 𝐸 and 𝐸 𝑀 . Thus, the configuration W-90-1000 was obtained with the best indicators and the configuration W-150-1000 with the worst one (Fig. 16). The microstructures of both the base material and the weld nugget were analyzed for the two above-mentioned welding configurations. For each of these, the specimen dedicated to the microstructure analysis was the one located in the vicinity of the specimen mechanically tested. The lateral face of the specimen was first polished to a depth of 3 mm using 600 to 1200 grit papers before it received a diamond polishing (3µm-granularity) to achieve a smooth and even surface. Then, a chemical attack, which was formulated with 50% of distilled water and 50% of 𝐻𝑁𝑂 3 , was performed for 30 seconds and followed by rinsing with ethanol. One can see that the microstructure of the base material (Fig. 16-a) is quite uniform and polyhedral, and is composed of grains having an average size of about 35µm. The microstructure of the weld nugget of the specimen W-90-1000 (Fig. 16-b) is less uniform and is composed of grains having a smaller average size (21µm) and strongly deformed with diffuse edges compared to those of the base material. In Fig. 16-d, one can see that the weld nugget of the specimen W-150-1000 is composed of grains having an even smaller average size (about 15 µm) than those of the weld nugget of W-90-1000. The welding temperature determines the dynamic recrystallization of the material in both configurations. However, the lower welding temperature of the configuration W-150-1000 led to the appearance of grains having an average size smaller than the one of grains from W-90-1000. This has also influenced the microhardness of both configurations (Fig. 16-c). The microhardness of W-150-1000 is higher than the one of W-90-1000. However, the insufficient heat induced by the welding of W-150-1000 did not allow the material to be in its optimal plastic state, thus generating weld defects of type "kissing bond" (Fig. 16-f). This type of defect has not been detected in the configuration W-90-1000 (Fig. 16-e).

One last result is about the failure localization referring to the joint axis. Fig. 15 shows the observations from each welding configuration tested. It can be seen that among the seven tested specimens only the welding configuration W-90-1000 exhibits a fracture path localized outside the welded joint. This indicates that for this specimen, the process parameters (welding speeds) are optimal.

Moreover, one can notice that for each of the other FSW configurations, the resulting fracture path is localized within the weld bead and always at the retreating side (RS).

These findings are in good agreement with the prediction of the optimum parameters.

Conclusion

This work has enabled us to study the friction stir welding process applied on a quasipure copper material. A robust and efficient numerical strategy is proposed and aims at predicting the optimal welding parameters for which a butt joint welding of 3-mm thick quasi-pure copper exhibits the maximum values of both mechanical strength and strain to fracture.

The simulation tool, which is based on the CEL (Coupled Eulerian-Lagrangian) model, takes into account some key aspects that make it robust and efficient: firstly, a friction coefficient that depends on the temperature has been used and its initial value was computed based on the ambient temperature and using the inverse identification method;

secondly, the computational time has been reasonably optimized (two days for the welding configurations simulated in this work) thanks to the use of the mass scaling technique, which has a low impact of the accuracy of the results since the temperature field is lowly overestimated (less than 10%); thirdly, the values of the constants of Johnson Cook's law were determined experimentally.

The key parameter of this simulation strategy lies in evaluating the optimum value of the welding temperature, which is intrinsic to both the welding configuration and material.

This optimum value of welding temperature can be often found in the literature. However, for welding configuration that has never been studied, this optimum welding temperature could be identified through extremely long and costly experimental investigation.

  -b. The x-rays pass through the specimen, and two of its faces (Fig.8-b) are projected for the defects analysis: the frontal one bounded by points A, B, C, and D, and the lateral one bounded by points B, F, E, and C.For each welding configuration to be investigated through monotonic tensile tests, three specimens have been extracted from the welded stabilized area and at the following positions (from point A, Fig.5-a): 55 mm, 70 mm, and 85 mm respectively. The geometry of the extracted specimens is shown in Fig.8-a. Additional details about the specimen geometry are 𝑙 0 = 100 𝑚𝑚 and 𝐿 = 200 𝑚𝑚.During tensile tests, the Digital Images Correlation (DIC) method was used to monitor the local strain field 𝜀 𝑦𝑦 simultaneously on both the frontal face (area bounded by points G, H, I, and J in Fig.8-b) and lateral face (area bounded by points H, K, L, and I in Fig.8b).

Figure 1 -Figure 2 -Figure 4 -

 124 Figure 1 -Geometrical model, boundary conditions, and mesh of the FEM

Figure 8 -

 8 Figure 8 -Geometry and monitored faces of the sample

Figure 9 -Figure 10 -Figure 11 -Figure 16 -Figure 1 -Figure 8 -

 910111618 Figure 9 -Tensile test set-up to investigate the mechanical behavior of the FSW joint

Table 1 .

 1 Thus, the simulation model developed in this work allows identifying the best combinations of welding speeds (𝑣 𝑎 , 𝑣 𝑟 ). An extensive testing campaign involving tensile test, X-ray technique, analysis of local strain fields, analysis of fracture surfaces highlights the robustness of the simulation strategy proposed in this work. The experimental results were in good agreement with the FE simulations and then have enabled us to determine the suitable set of FSW parameters for the studied material. Xu, N., Ueji, R., Fujii,H. Dynamic and static change of grain size and texture of copper during friction stir welding, J Mater Process Technol, 232, 90-99(2016). http://dx.doi.org/10.1016/j.jmatprotec.2016.01.021[37] Heidarzadeh, A., Saeid, T., Klemm, V et al. Effect of stacking fault energy on the restoration mechanisms and mechanical properties of friction stir welded copper alloys, Mater Des, 162,185-197, (2019). https://doi.org/10.1016/j.matdes.2018.11.050.[38] Heidarzadeh, A., Mironov, S., Kaibyshev, R. et al. Friction stir welding/processing of metals and alloys: A comprehensive review on microstructural evolution, Prog Mater Sci, 117, (2021). https://doi.org/10.1016/j.pmatsci.2020.100752 Physical properties of the DHP Copper[START_REF]Physical and Mechanical Properties of Pure Copper[END_REF]

	8. Acknowledgments
	This work was supported by a grant of the Romanian Ministry of Research and
	Innovation, CCCDI-UEFISCDI, project number PN-III-P3-3.1-PM-RO-FR-2019-
	0048/01.07.2019 and Campus FRANCE, France
	9. Declarations
	Ethical Approval
	Not applicable
	Consent to Participate
	Not applicable
	Consent to Publish
	Not applicable
	Authors Contributions
	Monica Daniela IORDACHE -conceived of the presented idea and
	supervised the project
	Claudiu BADULESCU -developed the theory and performed the
	computations
	Malick DIAKHATE -wrote the manuscript in consultation with Claudiu
	BADULESCU

Table 2 .

 2 Constants of Johnson-Cook's model for DHP copper

	Material Tmelt( o C) Tref ( o C) A (MPa) B (MPa)	C	n	m
	DHP-Cu	1083	22	250	250.4	0.0137	0.81	0.73

Funding

The research reported was founded partially by "Romanian Ministry of Research and Innovation, CCCDI-UEFISCDI, project number PN-III-P3-3.1-PM-RO-FR-2019-0048/01.07.2019" and "Campus FRANCE", France

Availability of data and material

Not applicable

Code availability

Not applicable

Marius Adrian CONSTANTIN -designed and performed the experiments

Eduard Laurentiu NITU -involved in planning and supervised the work, the analysis of the results, and the writing of the manuscript Younes DEMMOUCHE -aided in interpreting the results and worked on the manuscript Matthieu DHONDT -aided in interpreting the results and worked on the manuscript Denis NEGREA -performed the x-ray radiography measurements All authors discussed the results and commented on the manuscript

Conflicts of interest/Competing interests (include appropriate disclosures)

We know of no conflicts of interest or personal relationships that could have appeared to influence the work reported in this paper


Gérard Lassibille

Earnings distribution among Spanish engineers: research vs. non-research occupations

Based on individual survey data, we estimate separate earnings equations for engineers who hold a research occupation and engineers who work in other fields using a switching regression model. From these results, we compare the wage-generating process for each subgroup and evaluate the overall wage gaps across the working life of each category of engineers. On average, researchers receive lower wages compared to their peers. However, the estimation of the earnings models indicate that, holding everything else constant, engineers holding research positions have higher salaries at the beginning of their working life. Due to the recent growth of R& D activities in Spain, we explain these earnings gap by a transitory cohort effect. However, in this new context, we do not exclude that this differential might reflect also a deep change in the working conditions that private sector firms offer to scientists.

Introduction

Economic literature has devoted numerous studies to career choice decision of engineers and scientists, and to their earnings and turnover over the lifecycle (see, e.g., [START_REF] Park | The effect of schooling, postschool investments, and unemployment on the earnings of engineers: a human capital approach[END_REF][START_REF] Azevedo | Scientists, engineers and the reservation wage[END_REF][START_REF] Cooke | Turnover and earnings: the scientist and engineer case[END_REF][START_REF] Wilson | The Rate of Return to Becoming a Qualified Scientist or Engineer in Great-Britain, 1966-1976[END_REF][START_REF] Biddle | Private sector scientists and engineers and the transition to management[END_REF][START_REF] Ferrall | Levels of responsibility in jobs and the distribution of earnings among US engineers, 1961-1986[END_REF][START_REF] Lavoie | The occupational dynamics of recent Canadian engineering graduates inside and outside the bounds of technology[END_REF]Trem-blay et al., 1998). Due to the lack of appropriate data, there is no such study for Spain. This article attempts to fill this gap by analyzing the wages of engineers who hold research and non-research positions in private sector firms. Based on survey data conducted by a Spanish professional association of engineers, we estimate separate earnings functions for both categories of engineers using a switching regression model. From these estimations, we compare the underlying wage-generating process for each category of engineers and evaluate their experience-earnings profiles over their working life. Results show that researchers receive, on average, lower wages compared to their counterparts. However, the regression results indicate that, all else remaining the same, engineers in research path have higher salaries at the beginning of their working life than their peers. We explain this wage gap by the scarcity of new entrants in the labor market for scientists, and thus by a transitory cohort effect. However, given the recent growth of R& D activities in Spanish private sector firms, we do not exclude that this differential might reflect a deep change in the working conditions that business firms offer to scientists. The paper proceeds in the following manner. Section 1 presents briefly the data and the characteristics of engineers who hold research and non-research positions. Section 2 describes the methodological approach we use to estimate the wage gaps between the two categories of engineers. Section 3 presents the empirical results.

The data

Empirical results are based on data from a postal survey conducted by the Professional Association of Industrial Engineers in Catalonia (see Enginyers Industrials de Catalunya, Associació/ Col-legi, 1991Col-legi, , 1992)). The survey was administered to all male engineers working in the region; among those, about 30% answered the questionnaire.1 These cross-sectional data contain valuable information for the analysis attempted here. They detail the main activity and the job of each individual at the time the survey was conducted. They describe also their professional experience, their wages, and their level of responsibility within the firm. Apart from these personal characteristics, the survey gives also information on the firm where the individuals work (sector, size, national/multinational capital controlled firms). Note that the survey does not report any information on the educational background of respondents. This is because all Spanish engineers receive similar train-ing.

After completing secondary education, they enroll in Escuelas Tecnicas Superiores2 (5 years of higher education) which award all the same kinds of degrees. Differences in the quality of schools and students in such a non-differentiated system might theoretically exist. This said, in Spain, these differences are probably small, because the number of schools is very limited3 and candidates for a place in an engineering school are highly selected. In order to obtain more observations, we used data from both the 1991 and 1992 surveys. For individuals who responded to both surveys, only the 1992 questionnaire was used. The selected sample is restricted to individuals who are working full time. Self-employed engineers as well as engineers working in the public sector are excluded from the analysis. After deleting missing values, these exclusions left us with about 2700 observations. Table 1 presents selected statistics for the sample of engineers who hold research and non-research occupations (management, production, project). In this table, wages refer to annual labor earnings and are not defined by the hourly wage rates, because hours of work are not identified in the survey. The results show that engineers in non-research occupations receive, on average, higher wages compared to their counterparts; the earnings gap between both categories is about 20%. Researchers have about 4 years less professional experience than others. On average, engineers working in research occupations have less responsibility than their peers: 26% compared to 50% is assigned to the highest job level (level 1); 16% holds occupations which entail no supervision of other engineers (level 4) compared to only 10% for their counterparts. They are working more frequently in very large business firms (61% vs. 48%). Finally, engineers in research occupations are shown to work mainly in the electrical and manufacturing sectors, i.e. in high technology industries. These observations are in line with more general findings on the distribution of R&D activities in Spain (see Instituto Nacional de Estadística, 1993).

The conceptual framework

To evaluate earnings gaps between both categories of engineers, we follow the tradition by postulating a switching regression model (see, e.g., [START_REF] Lee | Unionism and wage rates: a simultaneous equation model with qualitative and limited dependent variables[END_REF][START_REF] Maddala | Limited Dependent and Qualitative Variables in Econometrics[END_REF]. Let us denote the research and non-research occupations as occupations 1 and 2, respectively, and express the corresponding wage equations as:

ln w 1 = Xβ 1 + u 1 , (1)
ln

w 2 = Xβ 2 + u 2 , (2)
where w stands for annual labor earnings, X is a vector of characteristics, β 1 is a parameter vector, and u 1 is the random error term. We assume that the expected benefits to the engineer are equal, or proportional to the differences in the wage rates between the two occupations. Thus, an individual will work in a research occupation if :

(ln w 1 -ln w 2 ) > Z 1 γ 1 + ε 1 , ( 3
)
where Z is a vector of explanatory variables associated with the probability of holding a research job, γ is a parameter vector, and ε is an error term.

Substituting Eqs. 1 and 2 into Eq. 3 and combining terms yields:

I* = Zγ + ε (4) I = 1 (research occupation), if I* ≥ 0 I = 0 (non-research occupation), if I* < 0
In this specification, I* is a partially observed dichotomous variable that describes the selection process. We observe the outcome (i.e. wages when holding a research or a non-research occupation) depending on whether I* is positive or negative. This equation summarizes a two-step process. First, the expected wage differential must be large enough for the individual to make it worthwhile to try to obtain a research job. The choice to apply is explained by the expected wage gap (net of the cost of applying for the job) and variables related to the preference for research occupation. Second, the employer determines whether the person is chosen for the job; explanatory variables are the person's attributes that matter to the employer. For this reason, the factors affecting the selection process (i.e. the vector Z) include both supply and demand variables. The switching regression model above is estimated using FIML, under the hypothesis that the three error terms (u 1 , u 2 , ε) follow a trivariate normal distribution with non-zero covariances. While OLS estimations of Eqs. 1 and 2 represent the unconditional wages in the two occupations (i.e. the expected earnings of engineers prior to the job choice process), the switching regression model gives an estimation of the wages conditional upon being selected for a particular occupation. The conditional wages of each category of engineers are given by:

E [ln w 1 | I* ≥ 0] = Xβ 1 + ρ 1ε σ 1 λ 1 , (5)
E [ln w 2 | I* ≥ 0] = Xβ2 + ρ 2ε σ 2 λ 2 , ( 6
)
Where γ 1 = [Ø(Zγ)/Φ(Zγ)] and λ 2 = [-Ø(Zγ)/{1-Φ(Zγ)}], Ø and Φ being, respectively, the standard normal density and cumulative distribution functions, ρ 1ε is the covariance between u i and ε, and σ i is the variance of u i .

The empirical results

Table 2 displays the estimation results of the switching regression model described above. It shows also the OLS estimations of the wage equations for each category of engineers. The variables included in the occupation choice model and the earnings equations are: professional experience, size of firm, sector of activity and whether the engineer is working for a multinational company. 4We discuss the coefficients of these variables only briefly, because for our purpose, the empirical results are best summarised in the earnings profiles presented further on.

Results displayed in Table 2 show first that the correlation between the error terms of the earnings equations and that of the selection equation is statistically significant; otherwise, stated OLS estimations of the earnings equations for each category of engineers are biased. 5In what concerns the choice of occupation, the main noteworthy feature in the results is that, all else remaining the same, the probability of holding a research job declines quite strongly as professional experience increases. Simulations not reported to save space, indicate, e.g., that the probability of holding a research job is around 6% for an individual who has 5 years of professional experience and 3% for 15 years of professional experience. The negative impact of professional experience may arise for several reasons: (a) as shown, e.g., by [START_REF] Biddle | Private sector scientists and engineers and the transition to management[END_REF] or [START_REF] Dasgupta | Information disclosure and the economics of science and technology[END_REF], mobility to production or management positions is the natural evolution for many researchers as these positions offer, on average, higher compensations than R& D occupations;6 (b) in Spain, R&D activities have grown substantially during the last two decades; in this new context, the negative impact of professional experience on the probability of holding a research occupation may only reflect a cohort effect (see the discussion below).

Turning to the earnings equations, results in Table 2 show that, holding everything else constant, the mean earnings of researchers increase more slowly over the life cycle compared to their peers.

Depending on the initial level of experience, the marginal returns to professional experience is between 0.3% and 1.1% lower for engineers working in research occupations (Table 3). According to the human capital theory (see [START_REF] Becker | Human Capital[END_REF], a possible reason for such differentiated returns to professional experience is a faster obsolescence of specific human capital in R& D occupations. We also found that, all else remaining the same, the type of firm a researcher works in makes no significant difference to his earnings. Conversely, engineers in non-research occupations working for a multinational company earn significantly more than their counterparts.

Using the regression results in Table 2, Fig. 1 shows the experience-wage profiles of both categories of engineers. These profiles refer to conditional wages; they are computed according to formulas 5 and 6 for representative individuals defined at the sample means.

The first noteworthy feature in the results is that starting salaries of researchers exceed slightly those of engineers who hold non-research occupations. The second feature to underline is that after a relatively short period of time (about 10 years), non-research engineers have higher wages compared to their peers, and the earnings gap between the two categories of engineers increases widely over the rest of the life cycle.

A potentially important reason for this pattern has to do with the opportunities of promotion in research and non-research jobs. The survey used in this paper offers some interesting insights on this question as it describes the responsibility level of engineers within their firms. Based on these data, we adjust the probability of holding a certain level of responsibility for each sample of engineers by postulating anordered probit model (see, e.g., [START_REF] Maddala | Limited Dependent and Qualitative Variables in Econometrics[END_REF]. Using these regression results, 7 Table 4 reports the predicted probability of being assigned to various levels of responsibility according to the occupation of the engineers and their professional experience; these predictions are computed in a standard way (see, e.g., [START_REF] Greene | Econometric Analysis[END_REF]. The results show that, other things being equal, engineers holding research occupations are more likely to be assigned to low levels of responsibility. For example, their probability of getting an executive position (level 1) after 10 years of professional experience is only 19% compared to 39% for their peers; at the end of their working life, their probability of holding occupations which entail no supervision of other engineers (level 4) is about three times larger than their non-research counterparts. Clearly, these results suggest that over their working life, researchers have less opportunity to move into better-paid jobs than engineers who hold non-research occupations. This conclusion is in line with other empirical findings (see [START_REF] Beltramo | Les rémunérations et les carrières des chercheurs en entreprise[END_REF]. The lower opportunity of promotion in research occupations may be probably explained by the internal organization of R&D activities, and in particular by a less pyramidal structure of occupations.

Lower career expectations and outcomes in research occupations, which lead many researchers to switch to management or production positions offering better opportunities, may explain the wage gaps we observe at mid-career between both categories ofengineers. However, they do not explain why researchers receive an earnings premium at the beginning of their working life (see Fig. 1). In Spain, it is more likely the changes in the labor market conditions for scientists in recent years (and thus, cohort effects) that seem to explain the wage differentials between the two categories of engineers early on in the working life.

7 See regression results in Appendix A.

As indicated in Table 5 below, R&D activities in Spanish firms grew rapidly during the 1980s: the total number of research related jobs in enterprises increased by about 130% and the number of scientists multiplied by almost four. These increases are much larger than those registered in other R&D sectors. As a matter of fact, in public administration and higher education, the number of workers employed in R&D activities increased, on average, by only 50%, and the number of scientists doubled at the most. One consequence of the rapid growth of R& D activities in the private sector has been a substantial shortage of science-related manpower. 8 During the decade under consideration, the number of graduates from engineering schools was still small, and as a result, the supply for scientists increased less rapidly than the demand. The gap in earnings early on in the working life between research and non-research en-gineers (see Fig. 1) probably reflects these time-related market conditions and consequently, the relative scarcity of scientists at the beginning of the 1990s. In other words, the earnings premium researchers receive at the early stage of their career is explained partly by cohort effects. The fact that senior engineers holding research positions earn less than their non-research counterparts is in accordance with such an explanation: in the past, R&D activities in the private sector were relatively scarce and demand for scientists was small compared with the supply. Obviously, the scarcity of R&D activities in the past explains partly the earnings gap between senior researchers and senior engineers who hold non-research occupations. Part of it is due also to selection issues, as in mid-career, many talented researchers switch to management or production po-sitions offering higher salaries. 9

Conclusion

Using a unique data set of engineers working in private sector firms, we focused on wage gaps between engineers holding research and non-research positions. We found that, on average, (i) researchers earn less than their peers; (ii) the rate of returns to professional experience is lower in 8 For an analysis of shortages in the labor market for engineers and scientists, see, e.g., [START_REF] Bosworth | The demand for qualified scientists and engineers[END_REF], [START_REF] Forbes | Growth, technology and the demand for scientists and engineers[END_REF], [START_REF] Hansen | The economics of scientific and engineering manpower[END_REF]. 9 Only panel data would make it possible to sort out the respective influences of each element in explaining earnings gap between each category of engineer at this stage of their working life.

R&D positions; (iii) opportunity of promotion within research occupations is smaller. These results are apparently in line with other empirical findings. However, in the Spanish context, substantial differences do exist. In fact, results from wage models that control for selectivity have shown that the two categories of engineers have different experience-earnings profiles. In particular, we found that at the beginning of their working life, researchers receive higher salaries compared to their counterparts. However, after a short period of time, the wages of engineers working in non-research paths exceed those of researchers. Part of the earnings gap observed at this stage is consistent with the idea that many engineers, and among them probably the most talented, start their careers as researchers and switch to management or production positions offering better-paid jobs.

The fact, that engineers employed in research occupations receive an earnings premium at the early stage of their working life, is apparently in line with the recent development of R&D activities in Spanish private sector firms; it reflects partly cohort effects and the scarcity of new entrants in the labor market for scientists. It might suggest also a deep change in the work conditions offered to scientists by business firms. Additional analysis is required to explore this issue. However, due to the recent growth of the labor market for scientists, more insight is needed to really sort out the influences on the wage-generating process of researchers of, on one hand, these possible new working conditions and, on the other, of the aforementioned generation effect.

Appendix A

As is usual in postal surveys, this response rate is low. Although we are not able to evaluate the bias induced by such a response rate, if we compare the characteristics of the sample with keys aspects of the distribution of R&D activities in Spain (see below), it seems that this bias is probably small.

ETS/Engineering Schools.

In some fields, e.g., like aeronautics, there is only one engineering school in the whole country.

We exclude this variable from the occupation choice model because engineers are equally distributed in each type of firms (see Table1). In our model, all variables that affect the selection process are wage determinants. As shown by[START_REF] Olsen | A least squares correction for selectivity bias[END_REF], this may increase the sensitivity of the results. However, we find no theoretical reasons to exclude from the wage equation variables included in the occupation choice model.

The bias is larger for engineers holding a non-research job since σ 1 ρ 1 is much larger than σ 2 ρ 2 .Earnings distribution among Spanish engineers

Due to data limitations, we are not able to analyze the switchfrom a scientific track to a managerial one.Earnings distribution among Spanish engineers


Mikhail Durnev

Alexey Kavokin

come