Adel Tayeb

Jean-Benoît Le Cam
email: jean-benoit.lecam@univ-rennes1.fr

Michel Grédiac

Evelyne Toussaint

Eric Robin

| Xavier Balandraud

| Frédéric Canevet

Identifying hyperelastic constitutive parameters with sensitivity-based virtual fields

Keywords: Sensitivity-based virtual fields, identification, hyperelastic materials, full-field measurements

This work deals with the identification of hyperelastic constitutive parameters using the virtual fields method. The choice of the virtual displacement fields is a cru-cial aspect of the method, typically for reducing the sensitivity to the measurement noise. A first and simple option is to generate the virtual displacement fields ran-domly. Nevertheless, in case of hyperelastic models for which the stress is not a linear function of the constitutive parameters, improving the choice of the virtual displace-ment fields is not trivial and an alternative strategy has to be found. In the present study, the sensitivity-based virtual fields approach is applied and compared to the randomly-generated virtual displacement fields approach. Two material models were considered: the Mooney model, which describes quite well the behavior of hyperelas-tic materials for small and moderate strains, and the Ogden model, which accounts for the stress hardening phenomenon observed at higher strains. The full kinematic fields are measured by using the digital image correlation technique during an equibiaxial tensile test performed on a cruciform specimen. Identification results are discussed through their capability to predict the external force measured during the test. The sensitivity-based virtual fields approach is found to improve significantly the predic-tion compared to the randomly-generated virtual displacement fields approach.

INTRODUCTION

The characterization of the mechanical behavior of rubber materials is classically carried out within the framework of hyperelasticity [START_REF] Holzapfel | Nonlinear solid mechanics: a continuum approach for engineering science[END_REF]. In this representation, the behavior describes only rate-independent effects without hysteresis, with the assumption of isotropy and incompressibility. In this case, the constitutive parameters depend upon the loading conditions applied in the identification procedure [START_REF] Guo | Application of a new constitutive model for the description of rubber-like materials under monotonic loading[END_REF]. This is the reason why they are generally identified from three homogeneous tests, namely the uniaxial tension (UT), the pure shear (PS) and the equibiaxial tension (EQT) 1 . An alternative approach consists in performing only one heterogeneous test in which the three loading states (UT, PS and EQT) are present along with various intermediate states.

Moreover, the values of the identified parameters are unique for each heterogeneous test, see [START_REF] Guélon | A new characterisation method for rubber[END_REF] and [START_REF] Promma | Application of the virtual fields method to mechanical characterization of elastomeric materials[END_REF]. Several approaches have been developed for the identification of constitutive parameters from experimental heterogeneous tests such as the finite element model updating technique, the equilibrium gap method or the virtual fields method, see [START_REF] Avril | Overview of identification methods of mechanical parameters based on full-field measurements[END_REF] and [START_REF] Pierron | General framework for the identification of constitutive parameters from full-field measurements in linear elasticity[END_REF] for further details. The latter method was applied before in the case of hyperelastic models and gave promising results, see [START_REF] Promma | Application of the virtual fields method to mechanical characterization of elastomeric materials[END_REF]. This is the reason why it will be used in this work for the identification of hyperelastic constitutive parameters. The Virtual Fields Method (VFM) was first introduced in [START_REF] Grédiac | Principe des travaux virtuels et identification[END_REF] for the identification of constitutive parameters of linear elastic materials. It was applied over the years to several types of material behaviors. Within the framework of the small strain assumption, the method was used to characterize orthotropic materials in [START_REF] Grédiac | Special virtual fields for the direct determination of material parameters with the virtual fields method. 2--application to in-plane properties[END_REF], [START_REF] Chalal | Experimental identification of a nonlinear model for composites using the grid technique coupled to the virtual fields method[END_REF] and [START_REF] Grédiac | A t-shaped specimen for the direct characterization of orthotropic materials[END_REF], and elasto-plastic materials in [START_REF] Grédiac | Applying the virtual fields method to the identification of elasto-plastic constitutive parameters[END_REF] and [START_REF] Avril | Stress reconstruction and constitutive parameter identification in plane-stress elasto-plastic problems using surface measurements of deformation fields[END_REF], among others. Dynamic properties were studied in [START_REF] Giraudeau | Simultaneous identification of stiffness and damping properties of isotropic materials from forced vibrating plates[END_REF] for instance. In the recent past, the VFM was extended to the finite strain framework. From a historical point of view, the first application of the VFM for hyperelastic materials was presented in [START_REF] Promma | Application of the virtual fields method to mechanical characterization of elastomeric materials[END_REF]. Mooney [START_REF] Mooney | A theory of large elastic deformation[END_REF] and Yeoh [START_REF] Yeoh | Some forms of the strain energy function for rubber[END_REF] models were considered using a biaxial tensile test of a 4-branch cruciform specimen in a moderate strain range. Other hyperelastic models such as the Ogden model [START_REF] Ogden | Large deformation isotropic elasticity-on the correlation of theory and experiment for incompressible rubberlike solids[END_REF] was employed in [START_REF] Palmieri | Virtual fields method on planar tension tests for hyperelastic materials characterisation[END_REF] using a planar tension tests on a flat rubber specimen subjected to the same moderate strain range as in the previous reference. Furthermore, in addition to hyperelastic parameters, the VFM was used for the determination of viscoelastic properties of rubber, see [START_REF] Yoon | High strain-rate tensile characterization of epdm rubber using nonequilibrium loading and the virtual fields method[END_REF] and [START_REF] Yoon | Application of the virtual fields method to a relaxation behaviour of rubbers[END_REF], where dynamic tests were exploited for the identification. This extensive use of the VFM for different mechanical behaviors and applications rose the challenge of the best choice of the virtual displacement fields, which is a crucial aspect of the method. Several techniques were used for addressing this problem when dealing with hyperelastic materials. In general, the virtual displacement fields are either chosen to cancel the actual loading term in the equation of the virtual work equation, see [START_REF] Yoon | High strain-rate tensile characterization of epdm rubber using nonequilibrium loading and the virtual fields method[END_REF] and [START_REF] Yoon | Application of the virtual fields method to a relaxation behaviour of rubbers[END_REF], or randomly, with a criterion on the conditioning of the linear system arising from the VFM [START_REF] Promma | Application of the virtual fields method to mechanical characterization of elastomeric materials[END_REF]. In a recent work, a new procedure in the choice of independent virtual displacement fields has been introduced in the case of anisotropic plasticity with small deformations [START_REF] Marek | Sensitivity-based virtual fields for the non-linear virtual fields method[END_REF]. The method is based on the sensitivity of the stress field to changes of the constitutive parameters. The virtual displacement fields are then computed proportional to the stress sensitivity fields. This method was extended to finite strains in [START_REF] Marek | Extension of the sensitivity-based virtual fields to large deformation anisotropic plasticity[END_REF]. In the present work, the VFM is applied to the identification of hyperelastic constitutive parameters from an equibiaxial experiment performed with a 4-branch cruciform shaped specimen. Two strategies in the generation of the virtual displacement fields are applied and compared. First, the random virtual displacement fields strategy that was applied to the Mooney and the Yeoh models in [START_REF] Promma | Application of the virtual fields method to mechanical characterization of elastomeric materials[END_REF], is extended to the Ogden model. Contrary to the former ones, the Ogden model leads to a non-linear stress relationship with respect to the constitutive parameters. This property makes it challenging to choose the best virtual displacement fields generated randomly. To overcome this issue, the sensitivity-based virtual displacement fields strategy has been developed within the framework of hyperelasticity by using the Mooney and Ogden models. This is the main contribution of the paper. The paper is organized as follows. In the first section, the theoretical background of the VFM is recalled and the identification problem is formulated for the two hyperelastic models under consideration. The second section is devoted to the experimental setup, including the mechanical test and the digital image correlation technique. In the third section, the results of the identification procedure are discussed in terms of virtual displacement fields to be used in the identification procedure and the ability of the obtained parameters to accurately predict the external force applied during the test. In the final section, the sensitivity-based virtual fields method is applied to four different hyperelastic materials. Concluding remarks close the paper.

THEORETICAL BACKGROUND

In this section, a literature survey is proposed about the use of the VFM for the identification of hyperelastic constitutive parameters. First, the basics of the method are recalled. Then the randomly-generated and sensitivity-based virtual displacement fields strategies are presented. Finally, the minimization problem arising from the application of the method is detailed. The Mooney [START_REF] Mooney | A theory of large elastic deformation[END_REF] and Ogden [START_REF] Ogden | Large deformation isotropic elasticity-on the correlation of theory and experiment for incompressible rubberlike solids[END_REF] models are used for this purpose.

Virtual fields method

Considering a hyperelastic solid subjected to a quasi-static loading and assuming that the volume forces can be neglected, the principle of virtual work, can be written as follows

-∫ ( , ) ∶ * ( , ) + ∫ ⋅ * ( , ) = 0, (1)
where is the volume of the solid, is the Cauchy stress tensor, * and * are the virtual displacement and strain, respectively. is the external force per unit deformed area. In general, the stress, virtual strain and displacement fields are all heterogeneous fields at any time . The stress field can be deduced from the measured displacement or strain fields by using the corresponding constitutive equations and a given set of material parameters. The virtual displacement * is a continuous function, which is kinematically admissible in the volume . A sufficient condition regarding the admissibility of the virtual displacement field is to consider a continuous function, which satisfies the boundary conditions, (see [START_REF] Pierron | The virtual fields method: extracting constitutive mechanical parameters from full-field deformation measurements[END_REF]). The virtual strain field is deduced by using the strain-displacement relationship, namely * = 1∕2 ∇ * + ∇ * . It should be noted that these two functions are also called test functions or weighting functions. They are completely independent from the real strain and displacement fields. In Eq.

(1), the first integral is the contribution of the internal virtual work due to deformations. The second integral accounts for the contribution of the external loading applied to the solid. Given that most full-field measurement techniques provide information on the specimen's surface, a two dimensional problem should be considered. Typically, a plane stress problem is considered in the case of thin specimen. Under this assumption, Eq. ( 1) becomes

-∫ ( , , ) ∶ * ( , , ) + ∫ ⋅ * ( , , ) = 0, (2)
where is the thickness of the solid, is the surface of the solid in the normal direction to the thin dimension and is its boundary. Full-field measurements, especially when performed with digital image correlation are based on the decomposition of the surface of the specimen into a set of subsets. The displacement at the centre of these zones is then deduced by minimizing the optical residual (see [START_REF] Bornert | Assessment of digital image correlation measurement errors: methodology and results[END_REF] and [START_REF] Vacher | Bidimensional strain measurement using digital images[END_REF]). The displacement, strain and stress fields are discrete functions, which are evaluated only at these discrete points. The internal virtual work part of Eq. ( 1) becomes

∑ =1 ⋅ * , ( 3
)
where is the number of data points in the region of interest ROI. , * and are the stress, virtual strain and area surrounding each data point, respectively. The external virtual work contribution is computed using the second integral of Eq. [START_REF] Avril | Overview of identification methods of mechanical parameters based on full-field measurements[END_REF]. Given that only the resulting force and not its distribution is measured along the edge of the specimen , two methods can be used for calculating this part of the virtual work. The first one consists in considering a constant virtual displacement along the edge. The contribution is then obtained by multiplying the resulting force by the constant virtual displacement. The second method consists in assuming a uniform force distribution along the edge and using a non-constant virtual displacement. The product between these two quantities is then integrated over . The stress in Eq. ( 3) is a function of the set of constitutive parameters to be identified. This set is denoted here by , and the strain by . The VFM is based on the application of the principle of virtual work and on the minimization of the difference between internal and external virtual works. For materials with linear dependency of the stress upon the constitutive parameters, the principle of virtual work is applied at least as many times as the number of constitutive parameters contained in . In this case, such as for linear elastic materials, this minimization leads to a linear system which gives the constitutive parameters after inversion. However, for other models such as those used with elasto-plastic materials, for which the stress is nonlinear with respect to the constitutive parameters, a nonlinear least-square minimization procedure must be employed. One or several specimen strain levels are used in the minimization. Two forms of the cost-function are proposed accordingly:

( , ) = ∑ =1 ∑ =1 ( , ) ⋅ * ( ) - * 2 , ( 4
)
for one specimen strain level and

( , ) = ∑ =1 ⎡ ⎢ ⎢ ⎣ ∑ =1 ∑ =1 ( , , ) ⋅ * ( ) ( ) - * ( ) ( ) 2 ⎤ ⎥ ⎥ ⎦ , ( 5
)
for several specimen strain levels. In Eqs. ( 4) and ( 5), and represent the number of independent virtual fields and time steps, respectively. To apply this method to hyperelastic materials, Eqs. (1) to (5) should be expressed in the finite strain framework. As in [START_REF] Promma | Application of the virtual fields method to mechanical characterization of elastomeric materials[END_REF] and [START_REF] Pierron | The virtual fields method: extracting constitutive mechanical parameters from full-field deformation measurements[END_REF] among others, the principle of virtual work is expressed as a function of the first Piola-Kirchhoff stress tensor , which represents the current force per unit of undeformed area. In this case, the principle of virtual work can be expressed as follows

-∫ 0 ( , ) ∶ * ( , ) 0 + ∫ 0 ( ⋅ ) ⋅ * ( , ) 0 = 0, (6)
where subscript • 0 , , and designate quantities given in the reference configuration2 , the coordinates, the thickness and the vector normal to the edge in the same configuration, respectively. Eqs. ( 4) and ( 5) become

( ) = ∑ =1 ∑ =1 ( ) ⋅ * ( ) ( ) - * ( ) 2 , ( 7
)
when considering one specimen strain level and

( ) = ∑ =1 ⎡ ⎢ ⎢ ⎣ ∑ =1 ∑ =1 ( , ) ⋅ * ( ) ( , ) - * ( ) ( ) 2 ⎤ ⎥ ⎥ ⎦ , ( 8
)
when considering several specimen strain levels. The identification of the constitutive parameters is carried out by using the fminsearch function of Matlab for nonlinear models, and by inverting the linear system for linear models. It should be noted that this function can lead to a local minimum. This issue was tracked by reconstructing the force-displacement curve with set of identified parameters. Indeed, curves plotted with parameters extracted from local minima are not in good agreement with the experimental one. In this case, these parameters are discarded and a new identification is performed with different initial values. The procedure can be repeated until a good agreement between the curves is obtained. Note that in all the cases reported here, launching the procedure once was sufficient. Further details on the procedure are presented in the remainder of this section.

Hyperelasticity

The mechanical behavior of elastomers is generally described within the framework of hyperelasticity. The stress is obtained by deriving the strain energy density with respect to the corresponding strain measure. This strain energy density relates the stress to the strain through the principal stretches or the principal invariants of the Cauchy-Green strain tensors. In the following, the material is assumed to be incompressible. Under this assumption, the first Piola-Kirchhoff stress tensor is given by

= --+ , ( 9
)
where is a Lagrange multiplier due to incompressibility (see [START_REF] Truesdell | The non-linear field theories of mechanics[END_REF]), is the deformation gradient tensor and • designates the transpose of a second-order tensor. For the Mooney model [START_REF] Mooney | A theory of large elastic deformation[END_REF], the strain energy density reads as follows

= 1 1 -3 + 2 2 -3 , ( 10
)
where 1 and 2 are the first and second invariants of the right Cauchy-Green strain tensor, respectively. 1 and 2 are the constitutive parameters to be identified. Combining Eqs. ( 9) and ( 10) and replacing by its expression in Eq. ( 6) lead to the following expression of the principle of virtual work for this material

1 ∫ 0 ∶ * 0 + 2 ∫ 0 ∶ * 0 = ∫ 0 ( ⋅ ) ⋅ * 0 , ( 11
)
where and are two functions of the principle stretches defined in Appendix A. Using Eq. ( 11) with two independent virtual displacement fields leads to the following system

= with = ∫ 0 ∶ * (1) 0 ∫ 0 ∶ * (1) 0 ∫ 0 ∶ * (2) 0 ∫ 0 ∶ * (2) 0 = 1 2 and = ∫ 0 ( ⋅ ) ⋅ * (1) 0 ∫ 0 ( ⋅ ) ⋅ * (2) 0 . ( 12
)
This linear system gives the two constitutive parameters 1 and 2 after inversion. The second model considered in this work is the Ogden model [START_REF] Ogden | Large deformation isotropic elasticity-on the correlation of theory and experiment for incompressible rubberlike solids[END_REF] for which the strain energy density is expressed by the following quantity

= ∑ =1 2 2 1 + 2 + 3 -3 , ( 13
)
where 1 , 2 and 3 are the principal stretches and , ; = 1.. are the constitutive parameters to be identified. The stretch is defined as the ratio between the current and the initial length in a given direction. Note that expressing the strain energy density with respect to the principal stretches implies that all tensors are expressed in the principal basis of the strain tensor. Hence, the first step is to determine the principal basis for each data point, and then to express the displacement and strain fields in these bases. From Eqs. ( 9) and ( 13) the eigenvalues of the Piola-Kirchhoff stress tensor are given by

Π = --1 . ( 14
)
The indeterminate coefficient is identified assuming a plane stress state (Π 3 = 0). Replacing Π in the principle of virtual work leads to the following equation

-∫ 0 Π 1 ⋅ * , + Π 2 ⋅ * , 0 + ∫ 0 ( ⋅ ) ⋅ * ( , ) 0 = 0, (15)
where ( , ) is the principal basis for the strains. In this basis, the cost function defined in Eq. ( 8) becomes

( ) = ∑ =1 ⎡ ⎢ ⎢ ⎣ ∑ =1 ∑ =1 Π 1 ( ) ⋅ * ( ) , + Π 2 ( ) ⋅ * ( ) , - * ( ) 2 ⎤ ⎥ ⎥ ⎦ . ( 16
)
In this work, two orders were considered for the Ogden model, namely are the first and the second orders. The identification of the constitutive parameters was performed by minimizing the cost function . At least 2 or 4 independent virtual displacement fields were used in the identification procedure for the first and second order Ogden models, respectively. It should be noted, however, that only the second-order Ogden model satisfactorily describes the stress-hardening phenomenon observed in the stress-strain curves.

Choice of the virtual displacement fields

There is an infinite number of virtual fields * that satisfy the principle of virtual work in Eq. ( 6). The choice of a set of independent virtual fields remains a critical issue. A first attempt was to consider arbitrary expressions for the virtual fields defined over the whole specimen. Examples are available for isotropic and anisotropic elastic materials, see [START_REF] Grédiac | Special virtual fields for the direct determination of material parameters with the virtual fields method. 1--principle and definition[END_REF], [START_REF] Pierron | The virtual fields method: extracting constitutive mechanical parameters from full-field deformation measurements[END_REF] and [START_REF] Grédiac | A t-shaped specimen for the direct characterization of orthotropic materials[END_REF].

Piecewise displacement virtual fields were also used in [START_REF] Toussaint | The virtual fields method with piecewise virtual fields[END_REF]. In the case of hyperelasticity, virtual displacement fields can be chosen randomly. This option is detailed in section 2.3.1. In the present work, the sensitivity-based virtual displacement fields strategy introduced in [START_REF] Marek | Sensitivity-based virtual fields for the non-linear virtual fields method[END_REF] is extended to hyperelasticity. It is presented in section 2.3.2. It should be noted that the hyperelastic model has a significant effect on these strategies and their robustness. This point is addressed in what follows.

Randomly-generated virtual displacement fields

To the best of the authors' knowledge, the VFM was firstly applied for hyperelastic materials in [START_REF] Promma | Application of the virtual fields method to mechanical characterization of elastomeric materials[END_REF]. Motivated by a noisesensitivity reduction, a set of randomly-generated virtual displacement fields was used. The procedure used in [START_REF] Promma | Application of the virtual fields method to mechanical characterization of elastomeric materials[END_REF] relies on the division of the (ROI) into several quadrangular sub-domains, and on the generation of random values for the virtual displacements at the nodes of these sub-domains. Then, the displacement is interpolated in the whole domain using four-noded quadrangular finite element formulations [START_REF] Zienkiewicz | The finite element method[END_REF]. These displacement fields take the following expressions as functions of the nodal

virtual displacements * = 4 ∑ =1 ( , ) * ( ) * = 4 ∑ =1 ( , ) * ( ) , ( 17
)
where are the shape functions, and are the coordinates in the reference element. These displacements are given in the reference element coordinate system ( , ). They should therefore be given in the global coordinate system ( , ). This is achieved by defining the transformation between these two coordinate systems. Note that in the case of a regular element grid, this transformation may be omitted and the global coordinate system ( , ) can be used in Eq. 17. In [START_REF] Promma | Application of the virtual fields method to mechanical characterization of elastomeric materials[END_REF], the randomly-generated virtual displacement fields were used for the Mooney [START_REF] Mooney | A theory of large elastic deformation[END_REF] and the Yeoh [START_REF] Yeoh | Some forms of the strain energy function for rubber[END_REF] hyperelastic models, for which the application of the VFM leads to a linear system. In this case, to ensure the independence of the virtual fields, a good conditioning of the linear system (12) must be obtained. However, for models for which the VFM does not lead to a linear system such as the Ogden model, no criterion was found to generate independent virtual fields with this approach.

Sensitivity-based virtual displacement fields

In a recent work [START_REF] Marek | Sensitivity-based virtual fields for the non-linear virtual fields method[END_REF], a new procedure for generating independent virtual displacement fields was employed for the identification of the constitutive parameters of an anisotropic plastic model in the small strain domain. The method is based on the sensitivity of the stress to changes of the constitutive parameters, typically between 10 and 20 %. The virtual displacement fields are then generated proportionally to the stress sensitivity fields by using an approach inspired from the finite element formulation. These stress sensitivity fields being computed from time-dependent strain fields, the former are time dependent as well as their corresponding virtual displacement fields.The method was then extended to the finite strain domain in [START_REF] Marek | Extension of the sensitivity-based virtual fields to large deformation anisotropic plasticity[END_REF] for anisotropic plasticity. In the finite strain domain, the stress sensitivity field is defined by

( ) ( , ) = + , -( , ) , ( 18
)
where 0.1 ≤ ≤ 0.2 is the range of variation of the ℎ parameter from a mean value ̄ obtained from a literature survey. Note that the stress sensitivity of Eq. ( 18) gives the influence of each constitutive parameter to the global response of the material at each point since the stress field is heterogeneous for this experiment. The virtual displacement fields were related to the stress sensitivity fields by the relation inspired from the finite element formulation, which gives

( ) ( , ) = * ( ) ( ), (19)
where is the global strain-displacement matrix obtained from a virtual mesh generated a priori. This matrix is obtained from the assembly of the elementary strain-displacement matrix obtained directly from the differentiation of the shape functions in Eq. ( 17) with respect to the coordinates in the ( , ) plane. * ( ) in Eq. ( 19) designates the virtual displacement field corresponding to the ℎ constitutive parameter. It is worth remembering that this displacement field is a test function, which is not related to the actual displacement field. In practice, matrix should be modified to account for the boundary conditions of the ROI. Typically, for edges where the external loading is unknown, a null displacement should be imposed to account for this constraint. A new matrix ̄ is therefore obtained from the original matrix . The virtual displacement field is then expressed as follows

* ( ) ( ) = ̄ ( ) ( , ) , (20)
where designates the pseudo inverse operator. Once the virtual displacement field is obtained, its gradient involved in the virtual work principle is deduced by using the following classic relationship inspired from the finite element method * ( ) ( ) = * ( ) ( ).

The contribution of each constitutive parameter to the response of the material is specific and may be very different in magnitude.

A scaling in the cost function should therefore be added, see [START_REF] Marek | Sensitivity-based virtual fields for the non-linear virtual fields method[END_REF] and [START_REF] Marek | Extension of the sensitivity-based virtual fields to large deformation anisotropic plasticity[END_REF]. Eqs. ( 7) and ( 8) become

( ) = ∑ =1 1 ( ) 2 ∑ =1 ( ) ⋅ * ( ) ( , ) - * ( ) 2 , ( 22
)
when one specimen strain level is considered and

( ) = ∑ =1 ⎡ ⎢ ⎢ ⎣ 1 ( ) 2 ∑ =1 ∑ =1 ( , ) ⋅ * ( ) ( , ) - * ( ) ( ) 2 ⎤ ⎥ ⎥ ⎦ , ( 23
)
when several specimen strain levels are considered. In Eqs.( 22) and ( 23), ( ) denotes the mean internal virtual work IVW obtained from the virtual fields corresponding to the ℎ parameter, which is expressed as follows

( ) ( ) = ∑ =1 ( , ) ⋅ * ( ) ( , ) . ( 24
)
In Eq. ( 24), the first Piola-Kirchhoff stress tensor is evaluated by using the set of reference values for the constitutive parameters. In the identification procedure, the number of images considered (i.e.

) strongly depends on the dependence of the stress upon the constitutive parameters. In general, a large distribution of the global stretch is needed in order to successfully perform the identification. Therefore, a heterogeneous experiment is needed.

EXPERIMENTS

The aim of this section is to present the experimental setup and to briefly recall the main characteristics of the DIC measurement system used in this study.

Material and specimen geometry

The material used in this study is a carbon black filled natural rubber. The specimen is shown in Figure 1. It is a 105 long and 2 thick cruciform specimen with 6 cylindrical ends to avoid slippage in the testing machine grips. A similar specimen geometry was firstly used in [START_REF] Promma | Application of the virtual fields method to mechanical characterization of elastomeric materials[END_REF] for the characterization of natural rubber. This specimen shape gives various states of strain when an equibiaxial load is applied. The classic states presented previously, namely the UT, the PS and the EQT, are induced as well as various intermediate states of strain. Therefore, the single heterogeneous test used in this study provides an interesting alternative to these multiple homogeneous tests classically used for the identification of hyperelastic constitutive parameters. This procedure gives a complete response of the specimen for the various strain states obtained with a wide distribution of the biaxiality coefficient.

Loading conditions

The experimental setup is presented in Figure 2. It is composed by a home-made biaxial testing machine and a digital camera. The machine is formed by four independent RCP4-RA6C-I-56P-4-300-P3-M (IAI) electrical actuators controlled by PCON-CA-56P-I-PLP-2-0 controller and four PCON-CA (IAI) position controllers. All these components are controlled by an in-house LabVIEW program. The biaxial testing machine is equipped with two cell loads with a capacity of 1094 , which store the force variation in the two perpendicular directions. In this work, an equibiaxial load was applied to the cruciform specimen. The four independent actuators were linked to have the same movement. This movement is such that the specimen center was motionless during the test. Hence, a reference point was obtained at the center of the specimen with respect to the correlation procedure.

A displacement of 70 was applied to each branch at a loading rate of 150 ∕ for five cycles. This choice was made in order to identify the stabilized behavior of the specimen since the material accommodates during the first cycles [START_REF] Mullins | Softening of rubber by deformation[END_REF]. The maximum applied displacement corresponds to a global stretch3 of 2.33.

Figure 2 Experimental setup

Full-field kinematic measurement

During the mechanical test, images of the specimen surface were stored at a frequency of 5 using an IDS camera equipped with a 55 telecentric objective. The charge-coupled device (CCD) sensor of the camera has 1920 × 1200 joined pixels. The displacement field at the surface of the specimen was determined by using digital image correlation technique. It consists in correlating the grey levels between two different images of a given zone at two different levels of strain [START_REF] Sutton | Image correlation for shape, motion and deformation measurements: basic concepts, theory and applications[END_REF]. Before the test, the specimen was sprayed by a white paint in order to improve the image contrast. Hence, a black and white random grey field was obtained. During the test, the specimen was subjected to a uniform cold lighting thanks to a home-made LED lighting system. This ensures the uniform distribution of the lighting even at the highest strains reached. The correlation process was performed by the SeptD software [START_REF] Vacher | Bidimensional strain measurement using digital images[END_REF], which is dedicated to both finite and small strain domains. Following the recommendation of the DIC guide [START_REF] Jones | A good practices guide for digital image correlation[END_REF], both DIC hardware and analysis parameters are given in tables 1 and 2 respectively. Note that the size of the subset in could vary from one test to another depending on the positioning of the camera and its parameters. Thanks to the symmetry of the test, a rectangular region on one branch of the specimen is sufficient to apply the identification procedure described in Section 2. This choice ensures a trade-off between:

• a large distribution of the loading cases within the ROI,

• a large strain level for each loading case,

• the convergence of the DIC calculation for the chosen ROI during the whole experiment, including the case of large strains.

The rectangular ROI is represented in Figure 3. It corresponds to a zone from the specimen center to the cylinder at the end of its branch. The gauge block shown in Fig. 3 is used for converting pixels to millimeters. The result of the correlation is Figure 3 Region of interest with a 4 pixels step size the displacement field over the ROI at each data point. Strain and displacement gradient tensors fields are then derived from these displacement fields. The displacement fields obtained from the correlation process were stored by the SeptD software as two matrices whose components correspond to each data point of the ROI. Due to the large displacement applied during the experiment, the correlation was not achieved in some subsets (less than 4% of the total number of subsets). The displacement in these subsets was approximated through a polynomial interpolation from all the subsets nearby. The displacement fields are then smoothed using a mean filter in order to reduce the experimental noise, especially where significant gradient occurred. This filter was applied before and after differentiation of the displacement fields in the computation of the displacement gradient tensor. Further details are given in Section 4.

Camera

RESULTS

In this section, the results of the identification carried out with two procedures, namely the sensitivity-based and the randomlygenerated virtual displacement fields strategies, are highlighted. A comparison between the force obtained experimentally and from finite element simulations using the identified constitutive parameters is used to assess the efficiency of each procedure.

Experimental displacement fields

The displacement field obtained for a global stretch of 2.33 with the SeptD software is given in Figure 4. Note that the dark blue zones in the two maps correspond to the subsets where correlation could not be successfully performed. Results of the displacement fields smoothed by the procedure described above are reported in Figure 5. The displacement gradient fields were obtained from the displacement fields by classic differentiation with respect to the spatial coordinates. First, the displacement fields were interpolated over the whole subset using the displacement of the four nodes defining the corners by the following expression ( ,

) = ⋅ + ⋅ + ⋅ ⋅ + ( , ) = ⋅ + ⋅ + ⋅ ⋅ + ℎ . ( 25
)
The displacement gradient tensor was then deduced from Eq. ( 25). Thus

= + ⋅ + ⋅ + ⋅ + ⋅ . ( 26
)
Note that the eight constants ⋯ ℎ in Eqs. ( 25) and ( 26) are different from one subset to another. They are obtained by inverting the linear system of eight equations arising from Eq. ( 25). The displacement gradient fields are presented in Figure 6. These data were smoothed using the same filter as that used for the displacement fields. The distribution of the loading cases and the maximum principal stretches are reported in Figures 7(a) and 7(b), respectively. Fig. 7(a) highlights that the ROI chosen in this work contains a wide range of biaxiality coefficient defined by log min max . Indeed, it covers all loading cases lying between EQT and UT. Furthermore, a wide distribution of loading cases is observed in the ROI chosen (i.e. the points in the ( 1 , 2 ) plane cover a significant range). Figure 8 illustrates this distribution, in which both invariants 1 and 2 were plotted for the classic UT, PS and EQT loading cases, as well as for the present experimental results. Note that this distribution is given for the maximum applied displacement. For lowest applied displacements, the points in Fig. 8 are located from the origin of the ( 1 , 2 ) plane ( 1 = 2 = 3) to the distribution corresponding to the maximum applied displacement. This means that the maximum principal stretch distribution is large enough whatever the loading case considered. The maximum value of the principal stretch reaches 3 in the UT zone for the maximum value of the displacement equal to 70 . This value is large enough for activating the stress hardening phenomenon. Among the models considered in this study, this phenomenon is only described by the second-order Ogden model. The other two models can be therefore suitable for describing the behavior of the material up to a strain level corresponding to the occurrence of the stress hardening. The data obtained experimentally were used in the identification of the constitutive parameters for randomly-generated and sensitivity-based virtual displacement fields.

Identification from randomly-generated virtual fields

The identification procedure described in section 2.3.1 was applied herein for the determination of the hyperelastic constitutive parameters. As described above, the virtual displacement was randomly generated at the nodes of the correlation grid (data points). Then, the virtual displacement field at any point was obtained by interpolating the nodal displacements. For the Mooney model, the virtual displacement fields were chosen in such a way that the conditioning of the matrix in Eq. ( 12) was above The two virtual displacement fields used in the identification of the Mooney parameters are presented. For the Ogden model, no criterion was found for the choice of the virtual displacement fields in Figure 9. Hence, a set of 100 virtual displacement fields was randomly generated and used in the identification procedure, (i.e. the 100 randomly-generated virtual displacement fields were used in the cost-function of Eq. 16). Note that for both models, the identification procedure was performed several times in order to obtain an objective criterion which could assess the efficiency of the method. Furthermore, the initial values for the fminsearch function (denoted 0 ) were set to the reference ( ), minimum ( ) and maximum ( ) values reported 4. The parameters identified using this approach are reported in Table 3. The parameters for the three models appear to be in good agreement with the values found in the literature for natural rubber, (see [START_REF] Marckmann | Comparison of hyperelastic models for rubber-like materials[END_REF]). However, no significant change on the ability of the identified parameters to predict the mechanical response has been found for this model. This is due to the non-linearity of the stress upon the constitutive parameters for this model which led to the lack of an objective criterion in the choice of the virtual fields. The final values of the objective function for the Ogden models were 0.42 and 0.20 for the first and second orders, respectively.

Identification from sensitivity-based virtual fields

The identification procedure presented in section 2.3.2 was then used to obtain the hyperelastic constitutive parameters from the same test data. First, a virtual mesh had to be generated (it can be different from the correlation grid). Then, the virtual fields were generated proportionally to the stress sensitivity fields. The reference values for the parameters used in this work are given in Table 4. These values are found in the literature for similar materials, see [START_REF] Marckmann | Comparison of hyperelastic models for rubber-like materials[END_REF]. As an example, the two virtual fields used with the Mooney model are shown in Figure 10, where * (1)

1 , .., * (2)
2 designate the components of the two sensitivity-based virtual 26) displacement fields 4 . The values of the identified parameters are given in Table 5. Contrary to results obtained with randomlygenerated virtual fields, the identified parameters for all the models considered are in good agreement with the reference values of a natural rubber. The identification took in average about 75 seconds for each model. Note that the parameters reported in Table 5 are obtained for several simulations with different sensitivity parameters, i.e. with different virtual fields. Furthermore, the mean parameter values used in the generation of the sensitivity-based virtual displacement fields do not affect the final result of the identification. In fact, the mean values for each parameter could change within the reference range without affecting the final result of the identification. For the Ogden model, the mean least squared error (the value of the objective function of Eq. 23 at the end of the identification procedure) is about 1.5 10 -5 and 2.5 10 -6 for the first and the second orders, respectively. This values are negligible compared to the ones obtained with randomly-generated virtual fields reported in the previous section.

Comparison of the results

To evaluate the accuracy of the identified parameters, the biaxial experiment used in this work was simulated using Abaqus software for a plane stress problem using the parameters reported in Tables 3 and5. The finite element model is presented in Figure 11. The element type used was the CPS4 element, which is a four-noded bilinear plane stress quadrilateral. The four edges of the cruciform specimen were subjected to a tension displacement of 70 each. The displacement was blocked in the transverse direction along the edges of the specimen. For each set of parameters, the resulting force in every branch was compared to the experimental force obtained during the experiment. It is worth noting that the force was measured in the two perpendicular directions. Since the loading was symmetric, only one load cell was used in each direction. The results of the comparison are shown in Figure 12, in which SBVF and RVF refer to sensitivity-based virtual fields and randomly-generated virtual fields, respectively. For the RVF method, the Mooney model appears to fairly predict the response for a global stretch up to 2, which is just before the inflection point of the force-displacement curve. This is the usual validity range for the Mooney model. On the other hand, both Ogden models either underestimated or overestimated the force in the branch for the whole test and the different initial values for the fminsearch function. This is due to the procedure for the choice of the virtual displacement fields, which were randomly generated in this case. With this approach, no criterion was found in their selection for the Ogden model. For the SBVF method, the Mooney model also correctly predicted the experimental force for a global stretch up to 2.

Better results were obtained with the first-order Ogden model for a wider strain range, corresponding to a global stretch up to 2.15. The results given by these two models were very satisfactory since they do not take into account the stress hardening phenomenon. The second-order Ogden model predicts very well the experimental force for the whole test. Especially, the stresshardening phenomenon, which is observed by the inflection point in the experimental force curve for a global stretch around 2.15. Such a result could not be obtained with the RVF method. This shows the good ability of the SBVF method to generate optimized virtual displacement fields for hyperelastic behavior. For the Ogden model, only four loading levels were considered in the identification procedure. They corresponded to an applied

GENERALIZATION CAPABILITY OF THE IDENTIFICATION PROCEDURE

The section aims at evaluating the generalization capability of the proposed methodology by characterizing other hyperelastic materials. In order to illustrate the potential of the methodology, four hyperelastic materials that differ significantly in stiffness and stress hardening have been considered. They are denoted by M1, M2, M3 and M4. As in the previous section, the accuracy of the identified parameters is highlighted the fact that these parameters enable us to correctly predict the force measured during the test. Figure 13 shows these results. It can be concluded that, regardless of the load level obtained with the material and the amount of the stress hardening, the SBVF is a well suitable method to generate virtual displacement fields for hyperelastic models

CONCLUSION

In this work, the VFM was applied in order to identify the constitutive parameters of two hyperelastic models. Two methods were employed for the generation of the virtual displacement fields. The first method consisted in randomly generating the virtual displacement at all nodes of the correlation grid. Then, the virtual displacement was interpolated around every point using a finite element formulation. The second method is based on the sensitivity of the stress fields to changes of the constitutive parameters. The virtual displacement fields were then computed from the stress sensitivity fields following a finite element formulation. A sufficient number of virtual fields was generated to be used in the identification procedure. Two hyperelastic models were investigated in this work: the Mooney model and two orders of the Ogden model. For the Mooney model, good results were found using both the RVF and the SBVF methods in the corresponding valid deformation range. Indeed, for both methods, the force in the branch of the cruciform specimen obtained from numerical simulations using identified parameters was very close to the experimental force for a global stretch up to 2. For the Ogden model, the RVF method did not give good results for both the first and the second orders. The parameters did not accurately predict the experimental response of the specimen. On the other hand, the parameters obtained from the SBVF method were in good agreement with the experimental results in terms of the force in the cruciform branches. Furthermore, the stress hardening phenomenon was well described with the parameters identified for the second-order Ogden model. As in the experiments, the inflection point in the force-displacement curve was found around a global stretch of 2.15. Hence, for a global stretch greater than 2.15, the force-displacement curve slope increases significantly. This feature is very important in the design of rubber component undergoing large deformations in industrial applications. Consequently, for the models considered in this work, the SBVF method constitutes a rational procedure for the selection of the virtual fields when applying the VFM in the identification of hyperelastic constitutive parameters. This method could be applied for the identification of other parameters characterizing other physical phenomena observed for rubber materials such as viscoelasticity, Mullins effect or Payne effect. To this end, an adaptation of the method has to be done with respect to the loading conditions as well as the numerical implementation of the method. This will be the subject of future work.

Figure 1

 1 Figure 1 Specimen geometry (dimensions in )

Figure 4

 4 Figure 4 Experimental displacement fields in the ROI for a global displacement of 70

Figure 5

 5 Figure 5 Smoothed experimental displacement fields in the ROI for a global displacement of 70

Figure 6

 6 Figure 6Displacement gradient fields obtained from Eq.[START_REF] Sutton | Image correlation for shape, motion and deformation measurements: basic concepts, theory and applications[END_REF]

Figure 7

 7 Figure 7 Loading cases and maximum principal stretches for a global displacement of 70

Figure 8

 8 Figure 8 Loading states distributions over the ROI for a global stretch = 2.33, each point in the graph corresponds to a data point in the ROI

Figure 9

 9 Figure 9 Randomly-generated virtual fields for the identification of the Mooney model

Figure 10

 10 Figure 10 Sensitivity-based virtual fields for the identification of the Mooney model

Table 5 Figure 11 Figure 12 Figure 13 Figure A1 Figure A2

 5111213A1A2 Figure 11 Finite element model for the simulation of the biaxial experiment

Table 2

 2 DIC Analysis parameters

	IDS UI-3160CP Rev. 2

Table 3

 3 Identified hyperelastic constitutive parameters using randomly-generated virtual displacement fields

			-3	-	-
	Ogden 1	1	0.38	0.29	0.36
		1	2.39	2.85	2.31
	Ogden 2	1	0.073	0.378	0.145
		1	2	1.34	2.76
		2	0.375	0.24	0.227
		2	1.99	1.34	2.77
	Model	Parameter Reference	Max	Min
	Mooney	1	0.4	1	0.1
		2	0.04	0.1	0.01
	Ogden 1	1	0.768	0.84	0.68
		1	1.26	1.63	0.8
	Ogden 2	1	0.114	0.81	1.4 10 -4
		1	7.297	11.94	1.58
		2	0.772	1.6	-9.6 10 -3
		2	-0.295	1.31	-7.81

Table 4

 4 Reference values for the constitutive parameters displacement of 17.5, 35, 52.5 and 70. Moreover, the identification was carried out with the fourth loading, for which the Mullins effect was eliminated.

It should be noted that the EQT can be compared to the uniaxial compression (UC)[START_REF] Treloar | The physics of rubber elasticity[END_REF]

chosen here to be the undeformed configuration

the ratio between the current and the initial lengths of the specimen in a given direction

These displacement fields are reported in Appendix A for the other two models.

ACKNOWLEDGEMENTS

The authors thank M. Marek and Prof. Pierron from Southampton University for providing algorithms dealing with the sensitivity-based virtual fields approach. The authors thank the National Center for Scientific Research (MRCT-CNRS and MI-CNRS), Rennes Metropole and Region Bretagne for financially supporting this work, as well as the Cooper Standard company for providing the specimens.


Minh-Quang Pham
email: pham@limsi.fr

Josep Crego
email: crego@systrangroup.com

Franc ¸ois Yvon

Multi-Domain Adaptation in Neural Machine Translation with Dynamic Sampling Strategies

Building effective Neural Machine Translation models often implies accommodating diverse sets of heterogeneous data so as to optimize performance for the domain(s) of interest. Such multi-source / multi-domain adaptation problems are typically approached through instance selection or reweighting strategies, based on a static assessment of the relevance of training instances with respect to the task at hand. In this paper, we study dynamic data selection strategies that are able to automatically re-evaluate the usefulness of data samples in the course of training. Based on the results of multiple experiments, we show that our method offer a generic framework to automatically handle several real-world situations, from multi-source or unsupervised domain adaptation to multidomain learning.

Introduction

A typical setting in machine translation (MT) is to collect the largest possible collection of parallel data for the chosen language pair, with the intent to achieve optimal performance for the task of interest. In such situations, the training data distribution is opportunistic, while the test data distribution is chosen and fixed; a key aspect of training is then to mitigate the detrimental effects of a mismatch between these distributions. Single-source and multisource1 domain adaptation (DA) is a well-studied instance of this setting (see [START_REF] Chu | An empirical comparison of domain adaptation methods for neural machine translation[END_REF][START_REF] Saunders | Domain adaptation and multidomain adaptation for neural machine translation: A survey[END_REF]) for a review), and so is multi-domain (MD) learning [START_REF] Chu | Multilingual and multi-domain adaptation for neural machine translation[END_REF][START_REF] Zeng | Multidomain neural machine translation with word-level domain context discrimination[END_REF][START_REF] Jiang | Multi-domain neural machine translation with word-level adaptive layer-wise domain mixing[END_REF][START_REF] Pham | Revisiting multi-domain machine translation[END_REF]. A related situation is multilingual MT [START_REF] Firat | Multi-way, multilingual neural machine translation with a shared attention mechanism[END_REF][START_REF] Ha | Toward multilingual neural machine translation with universal encoder and decoder[END_REF][START_REF] Johnson | Google's multilingual neural machine translation system: Enabling zero-shot translation[END_REF][START_REF] Aharoni | Massively multilingual neural machine translation[END_REF], where the diversity of training data not only corresponds to variations in the topic, genre, or register but also in language.

This problem is often approached by static instance selection or re-weighting strategies, where the available training data is used in proportion to its relevance for the testing conditions [START_REF] Moore | Intelligent selection of language model training data[END_REF][START_REF] Axelrod | Domain adaptation via pseudo in-domain data selection[END_REF]. Finding the optimal balance of training data is however, a challenging task due, for instance, to the similarity between domains/languages, or to the regularization effects of out-of-domain data (Miceli [START_REF] Barone | Regularization techniques for fine-tuning in neural machine translation[END_REF]. A static policy may also be suboptimal when some target domains or languages are easier to train than others. Finally, improving the performance of the MT system in one domain will often hurt that of another (van der [START_REF] Van Der Wees | Dynamic data selection for neural machine translation[END_REF][START_REF] Britz | Effective domain mixing for neural machine translation[END_REF] and improving model generalization across all domains [START_REF] Koehn | Findings of the WMT 2018 shared task on parallel corpus filtering[END_REF] may not achieve optimally for any particular domain.

Several recent proposals explore ways to instead consider dynamic data selection and sampling strategies: van der Wees et al. (2017) and [START_REF] Zhang | Curriculum learning for domain adaptation in neural machine translation[END_REF] construct a static curriculum, while Wang et al. (2020a) and Wang et al. (2020b) build curricula that automatically adapt to the training data. In this paper, we contribute to this line of research in several ways.

• First, we propose a novel framework (Multi-Domain Automated Curriculum, MDAC for short), a variant of Differentiable Data Selec-tion (DDS) of Wang et al. (2020b), initially applied to multilingual NMT, that simultaneously accounts for the domain adaptation and the multidomain adaptation problems.

• We show that MDAC achieves performance that compare to fine-tuning strategies for DA ( § 5.1) and outperform some static data sampling strategies for multidomain settings (5.3).

• We show that our variant MDAC mitigates some failures of DDS in multidomain training.

• We illustrate the generality of differentiable data selection frameworks (both MDAC and DDS) on less common situations such as DA using unsupervised clustering ( § 5.5); DA using out-of-domain training data and small in-domain validation data ( § 5.4); and twodomain adaptation where the test distribution only mixes two of the training domain ( § 5.2).

Learning with multiple data sources

We conventionally define a domain d as a distribution D d (x) over some feature space X that is shared across domains [START_REF] Pan | A survey on transfer learning[END_REF]: in machine translation, X is the representation space for input sentences; each domain corresponds to a specific source of data, and may differ from other data sources in terms of textual genre, thematic content [START_REF] Chen | Guided alignment training for topic-aware neural machine translation[END_REF][START_REF] Zhang | Topic-informed neural machine translation[END_REF], register, style [START_REF] Niu | Multi-task neural models for translating between styles within and across languages[END_REF], etc. Translation in domain d is formalized by a translation function h d (y|x) pairing sentences in a source language with sentences in a target language y ∈ Y. h d is usually assumed to be deterministic (hence y = h d (x)) but may differ across domains.

It is usual in MT to opportunistically collect corpora from several domains, which means that training instances are distributed according to a mixture

D s such that D s (x) = n d d=1 λ s (d)D d (x), with {λ s (d), d = 1 . . . n d } the mixture weights satis- fying d λ s (d) = 1.
In the sequel, boldface λ denotes a vector with λ(d) the d th component of λ.

The main challenge in this situation is to make the best of heterogeneous data, with the aim to achieve optimal performance for the target test conditions. These might correspond to data from just one of the training domains, as in standard supervised domain adaptation; a more difficult case is when the test data is from one domain unseen in training (unseen domain adaptation); in multido-main adaptation finally, the test distribution is itself a mixture of domains, some of which may also be observed in training. We thus assume that the test distribution takes the form D t (x) = d λ t (d)D d (x) -with only one non-null component in the case of domain adaptation (see Figure 1). These situations have been amply documented from a theoretical perspective (Mansour et al., 2009b;Mansour et al., 2009a;[START_REF] Hoffman | Algorithms and theory for multiple-source adaptation[END_REF]. A general recommendation in the DA setting is to adjust the sampling distribution used to optimize the system so as to compensate for the mismatch between D s (x) and D t (x). This can be approximated by reweighting instances, or more conveniently domains, which are selected during training with a probability λ l (d), with λ l (d) = λ s (d).

A widely-used approach to supervised DA is fine-tuning [START_REF] Luong | Stanford neural machine translation systems for spoken language domain[END_REF], where λ l varies during learning. With our notations, this approach first learns an initial parameter value with all the data (∀d, λ l (d) = λ s (d)), then continues training with only batches from the test domain d t (λ l (d) = I(d = d t )), with I(A) the indicator function for predicate A. This strategy is potentially suboptimal as some out-of-domain samples may contribute to the final performance due to e.g. domain similarity. Optimizing the learning distribution in multidomain settings is even more challenging as the learner needs to take advantage of possible domains overlaps and also of the fact that some domains might be easier to learn than others.

3 Multi-Domain Automated Curriculum MDAC builds an adaptative training distribution λ l that optimizes the data selection policy along with the training of the model. We parameterize λ l by a differentiable function λ l (ψ), which is described in § 4.4. We divide the training into many short sessions; in each session t, the model is trained with a static data distribution λ l (ψ t ). After one learning session, we update the data distribution using the REINFORCE algorithm of [START_REF] Williams | Simple statistical gradientfollowing algorithms for connectionist reinforcement learning[END_REF]. The evolution of ψ is thus defined by:

ψ t+1 = ψ t + lr 1 * n d d=1 R(d) * ∂λ l (d; ψ t ) ∂ψ ,
where the reward R(d) is computed as:

R(d) = J t (θ t+k , λ t ) -J t (θ t , λ t ),
(1) and where we also define:

θ t+i = U pdate θ t+i-1 , [x i j , y i j ] N j=1 x i j , y i j ∼ D l d (x) J t (θ, λ t ) = n d d=1 λ t (d) x t d ,y t d ∈ D t d l(θ, x t d , y t d ).
In these equations, N denotes the size of a batch; lr 1 is the learning rate of the sampling distribution; l(θ, x, y) is the loss of the NMT model on sample (x, y); J t (θ, λ t ) is the weighted loss aggregated over n d dev-sets corresponding to the n d domains.

To compute the reward R(d) associated to training the model with data from domain d, we simulate k training steps from the current checkpoint, using k batches sampled from D l (d) and computing the gain of the weighted dev-loss. This computation is inspired by the target prediction gain of [START_REF] Graves | Automated curriculum learning for neural networks[END_REF]. However, where [START_REF] Graves | Automated curriculum learning for neural networks[END_REF] used accumulated gains from the past as rewards, we instead predict the usefulness of each domain for improving the future performance of the system given its current state. This is achieved by simulat-ing a round of training with only the data from one domain. We also differ from these authors in the parameterization of the sampling distribution.

The work of Wang et al. (2020b) is also related: it is based on the bi-level optimization framework, which aims to find an optimal static distribution λ l that will result in the best model with respect to a given target dev set at the end of training. These authors also derive a similar form of update for ψ. However, their reward is the cosine similarity between the gradient computed with the training data from one domain and the gradient computed with the dev set. We compare this approach with ours in the experiment section.

MDAC for (multi) domain adaptation

The setting developed in previous sections is quite general and can, in principle, accommodate the variety of situations mentioned above, and many more: basic DA, multidomain adaptation with various target distributions, possibly including domains unseen in training. In our experiments, we would like to better assess the potential of MDAC in these settings and seek to study the following questions:

• is MDAC a viable alternative to fine-tuning?

In particular, does it enable to better take advantage of relevant data from other domains?

• is MDAC a viable option in multidomain adaptation scenarios?

• does MDAC enable to perform unsupervised (multi-)domain adaptation?

These questions are further explored in Section 5. We now turn to our experimental conditions.

4 Experimental settings

Data and metrics

We experiment with translation from English into French in 6 domains, corresponding to the following data sources: the UFAL Medical corpus V1.0 (M E D)2 ; the European Central Bank corpus (B A N K); the JRC-Acquis Communautaire corpus (L A W) [START_REF] Steinberger | The JRC-Acquis: A multilingual aligned parallel corpus with 20+ languages[END_REF]; documentations for KDE, Ubuntu, GNOME and PHP from the Opus collection, merged in a I T-domain; TedTalks (T A L K) [START_REF] Cettolo | Wit 3 : Web inventory of transcribed and translated talks[END_REF], and the Koran (R E L). [START_REF] Papineni | BLEU: a method for automatic evaluation of machine translation[END_REF]. 4 Significance testing is performed using bootstrap resampling [START_REF] Koehn | Statistical significance tests for machine translation evaluation[END_REF], implemented in compare-mt5 [START_REF] Neubig | compare-mt: A tool for holistic comparison of language generation systems[END_REF]. We report significant differences at the level of p = 0.05.

Additional experiments use the News Commentary

L A W B A N K I T T A L K R E L N E W S #

Baseline systems

Our baselines are standard for multidomain settings. 6 Using Transformers [START_REF] Vaswani | Attention is all you need[END_REF] implemented in OpenNMT-tf7 [START_REF] Klein | OpenNMT: Opensource toolkit for neural machine translation[END_REF], we build the following systems:

• Generic models trained with predefined mixtures of the training data taking the form:

λ α (d) = ( n d d=1 q α d ) -1 (q α d ) q d = | N s d | N s (2)
with α ∈ {0, 0.25, 0.5, 0.75, 1.0}. We denote these as Mixed-α below. Mixed-0 uses a uniform distribution, Mixed-1.0 the empirical distribution of domains.

• fine-tuned models based on Mixed-1.0, further trained on each domain for at most 50 000 iterations, with early stopping when the dev BLEU stops increasing for 5 successive iterations. The fine-tuning (FT-Full) procedure updates all the parameters of the initial model, resulting in six systems, one per domain, with no parameter sharing across domains.

• systems trained with fixed mixtures with λ l ∈ λ 0 , λ 0.25 , λ 0.5 , λ 0.75 , λ 1.0 ; these are used in the multidomain experiments of § 5.3;

CL and DDS re-implementations

We re-implement DDS in Tensorflow without any change in the choices of parameterization and hyper-parameters compared to the original code of Wang et al. (2020b). 8 We also re-implement the approach of [START_REF] Zhang | Curriculum learning for domain adaptation in neural machine translation[END_REF] according to the authors' description. For each DA experiment, we combine the training data of all other domains into one corpus then compute the cross-entropy difference score of each source sentence of this combined dataset. We then sort and split the corpus into 9 shards and execute curriculum learning with 10 shards, using the in-domain data as the first shard.

MDAC systems

The behavior of MDAC only depends on (a) the initial domain distribution at the start of training λ l t=0 , and (b) the target (dev/test) distribution λ t . We thus report these systems as MDAC(λ l t=0 , λ t ) and compare with DDS using the same settings.

In our work, we parameterize the distribution λ l as follows (with β = 2 in all experiments): 9

λ l (d; ψ) = ψ[d] β i ψ[i] β .
This parameterization avoids the "rich-get-richer" effect that we observe with λ(ψ) = softmax(ψ), which yields gradients wrt. ψ[d] that are proportional to exp(ψ[d]) (see also Figure 2). Additional settings for the hyper-parameters of our method include the number of simulation steps k = 10 and the learning rate lr data = 0.001. We update the sampling distribution via 100 gradient descent iterations for almost all experimental settings except that for adaptation with automatic clusters ( § 5.5), where we use 20 gradient descent iterations to avoid converging to degenerate distributions. We split the training into 100 short sessions that last 2000 training steps each. The choice of those hyperparameters is mostly heuristic except for the learning rate lr data which is optimized via grid search over a set of values {0.001, 0.0025, 0.005}.

The computational cost of our approach is due to the simulation step, which is conducted after every 2,000 iterations to compute the reward of each domain (eq. ( 1)). During this step, we update the temporary checkpoint with k updates for each domain, which costs as much as k training updates. Therefore, we execute k × n d updates after every 2,000 iterations. Our algorithm approximately costs 1 + k×n d 2000 times as much as a standard training.

Experimental tasks

We evaluate our method in the 5 following conditions. In the supervised domain adaptation task, given the data from 6 domains (M E D, B A N K, L A W, I T, T A L K, R E L), we aim to build expert NMT models for each domain. To challenge the flexibility of the method, we also consider a two-domain adaptation task, where given the same 6 domains, we focus on adapting to a mixture of 2 domains.

In the multidomain adaptation task, we use the same 6 domains to build one single NMT model that should perform optimally, assuming a uniform distribution of domains during the test. A fourth experiment (unseen domain adaptation), adds to the training data for 6 domains a small dev set in a new domain (NE W S): our target is a model which performs well for the unseen domain. Finally, in the unsupervised domain adaptation task, we cluster all available training data into 30 clusters using the KNN algorithm as in [START_REF] Tars | Multi-domain neural machine translation[END_REF], then learn mixture weights these clusters to one of 6 domains using the corresponding dev set. We compare MDAC to DDS for each of our 6 test sets.

Results and discussion

Domain Adaptation

In this setting, we aim to build an NMT model for one single domain: we accordingly set λ t to a deterministic distribution λ d , where the target domain d has probability 1. We consider three initializations for MDAC and DDS, using λ 0 , λ 1 and λ d . According to Table 2, MDAC achieves the overall best performance when λ t=0 = λ 0 . Doing so proves much better than initializing with λ d for small domains: T A L K, B A N K andI T. Conversely, initializing with λ d is beneficial when targeting large domains such as M E D and L A W. The same conclusion holds for DDS.

We now compare the best MDAC system (using λ t=0 = λ 0 ) to full fine-tuning. According to Table 2, fine-tuning is better for large domains such as M E D and L A W, while MDAC outperforms finetuning by approximately 1.2 BLEU for B A N K and 1.0 BLEU for R E L. This suggests that for small domains, out-of-domain data helps improve the generalization and that MDAC is able to exploit both the in-domain and the out-of-domain training data instead of edging out the out-of-domain training data as in fine-tuning. Results for DDS display similar trends but are always outperformed by MDAC. Results for CL, which does only well the large domain M E D, lag somewhat behind.

Two-domain adaptation

In these control experiments, we showcase the flexibility of dynamic sampling and adapt to (arbitrary) pairs of target domains with equal weight, contrasting MDAC with DDS in Table 3. Here, MDAC significantly outperforms DDS in two settings (M E D+I T and L A W+B A N K) out of three.

Multi-domain NMT

We now turn to a more realistic scenario and consider multidomain NMT, which aims to train one single system with optimal performance averaged over 6 domains and targets a uniform test distribution λ t = λ 0 . In this situation, CL [START_REF] Zhang | Curriculum learning for domain adaptation in neural machine translation[END_REF] does not apply: we only contrast the performance of MDAC, DDS and several fixed training data distribution λ l ∈ λ 0 , λ 0.25 , λ 0.5 , λ 0.75 , λ 1.0 , where λ α is defined according to equation (2).

We again initialize MDAC and DDS with two distribution λ 0 and λ 1 . According to Table 4, MDAC achieves the best performance with initial (uniform) λ 0 . The same conclusion holds for DDS. For this Table 2: Single domain adaptation. We report BLEU scores of each method for 6 target domains and their average: each column corresponds to a distinct system. ( * ) MDAC is significantly better than CL, fine-tuning and DDS with p < 0.05. ( * * ) MDAC is significantly better than CL and DDS with p < 0.05. ( * * * ) MDAC is significantly better than CL, fine-tuning with p < 0.05.

domain d = M E D L A W B A N K T A L K I T R E L
configuration, MDAC outperforms static training distributions including λ 0 , λ 0.75 , λ 1.0 by a significant margin, and performs slightly better than λ 0.25 , λ 0.5 . Using MDAC thus dispenses with the empirical search of an optimal training mixture.

A second observation is that MDAC again outperforms DDS by a wide margin (+1.5 BLEU on average); the only domain where DDS does better is M E D. Figure 2, which plots the evolution of the mixture weights during training, helps to understand the difference between the two methods. For DDS (Figure 2a), the sampling distribution quickly reaches a bi-modal regime in which only M E D and R E L have significant probability -hence the good performance on the former domain. In contrast, the distribution computed by MDAC evolves more smoothly; small domains such as B A N K, I T, T A L K and R E L receive a larger part of training data in the early stages; their weights then slowly decrease as larger domains such as M E D and L A W increase their share. This only happens at the end of training, when some NMT models might already be close to their peak performance for the small domains.

Unseen domain

The left part of Table 5 displays the performance on the unseen domain N E W S for systems trained with mixtures λ l ∈ λ 0 , λ 0.25 , λ 0.5 , λ 0.75 , λ 1.0 and with dynamic data selection (MDAC and DDS). These systems have insignificant differences in BLEU, suggesting that dynamic mixtures do not improve the robustness of NMT systems against unseen domains. However, the performance of MDAC and DDS remains close to the best performance, showing that they also apply in such settings.

Automatic clustering

The right part of Table 5 reports the performance of NMT systems adapted to each domain. In comparison to Section 5.1, the training data is distributed in 30 automatic clusters instead of the 6 original domains. Splitting the train data into small groups gives the learner extra degrees of freedom when selecting the best distribution. However, as these clusters are built automatically, they are noisier in nature. According to results in Table 5, this scenario is hard both for DDS and MDAC, which performs much worse than for the supervised DA setting. This again signals the importance of initialization: analyzing the clustering, we find that the data for R E L mostly correspond to one single cluster. With a uniform initialization, this cluster starts with a small weight and never succeeds in matching the good performance observed in the DA setting.

Related Work

Domain adaptation is an old problem that has been studied from many angles, both for SMT and NMT. A survey of supervised and unsupervised DA for NMT is in [START_REF] Chu | An empirical comparison of domain adaptation methods for neural machine translation[END_REF], where the authors distinguish between data-centric and model-centric DA, a view also adopted in the recent survey of [START_REF] Saunders | Domain adaptation and multidomain adaptation for neural machine translation: A survey[END_REF]. Our approach to DA in this paper falls under the former category. We refer readers interested in DA to these papers.

Multidomain NMT (MDMT) aims to develop systems that simultaneously bode well for several domains. Like for DA, techniques for supervised MDMT combine one or several ingredients: (a) the specialization of data representations [START_REF] Kobus | Domain control for neural machine translation[END_REF] or of sub-networks [START_REF] Pham | Generic and Specialized Word Embeddings for Multi-Domain Machine Translation[END_REF] to differentiate the processing of each domain; (b) the use of adversarial techniques to neutralize differences between domains [START_REF] Britz | Effective domain mixing for neural machine translation[END_REF][START_REF] Zeng | Multidomain neural machine translation with word-level domain context discrimination[END_REF]; (c) the use of automatic domain identification e.g. [START_REF] Jiang | Multi-domain neural machine translation with word-level adaptive layer-wise domain mixing[END_REF]. Unsupervised MDMT is studied in [START_REF] Farajian | Multi-domain neural machine translation through unsupervised adaptation[END_REF], as an instance of unsupervised DA.

Most approaches to adaptive/dynamic data selection take inspiration from [START_REF] Bengio | Curriculum learning[END_REF], where the notion of curriculum learning is introduced. CL relies on the notion of the "easiness" of [START_REF] Kumar | Reinforcement learning based curriculum optimization for neural machine translation[END_REF] use reinforcement learning to learn the curriculum strategy: in this work, complexity corresponds to difficulty levels which are binned using contrastive data selection. The reward is based on the increase of the devset loss that results from the current data selection strategy. This technique is applied to multilingual NMT in [START_REF] Kumar | Learning policies for multilingual training of neural machine translation systems[END_REF]. [START_REF] Zhou | Uncertainty-aware curriculum learning for neural machine translation[END_REF] propose another CL-based approach which relies on instance uncertainty as a measure of their difficulty and presents data samples starting with the easiest. Another contribution of this work is a new stopping criterium. Closest to our problems, Wang et al. (2020a) adapt CL for multidomain NMT, where an optimal instance weighting scheme is found using Bayesian optimization techniques. Each step consists of (a) weighting instances based on relevance features, (b) fine-tuning a pretrained model using the weighted training set, and is applied to train a sequence of models. The one that maximizes the devset performance is finally retained.

domain d = M E D L A W B A N K T A L K I T R E L DDS(λ0,

Conclusion and outlook

In this study, we have presented a generic framework to perform multiple adaptation tasks for machine translation, ranging from supervised domain adaptation to multidomain NMT and unseen domain adaptation. In our experiments, we have shown that the same algorithm, aimed at automatically finding an effective data sampling scheme during the course of training, can be used in all these situations. This algorithm, we believe, provides us with a more sound approach to (multi-domain) DA than existing heuristics and dispenses with the costly search of optimal meta-parameters. Another contribution of our work is an experimental comparison of recent approaches to dynamic data selection.

Our future work will continue developing this approach and improve its effectiveness. One issue that we have left unaddressed is reward normalization, which is especially important in the early stages of training [START_REF] Kumar | Reinforcement learning based curriculum optimization for neural machine translation[END_REF]. Another area where we need to progress is the unsupervised learning setting of § 5.5, where our results lag behind supervised DA. This might be due to the inability of our simplistic optimization strategy to handle situations where the number of clusters is large.

Figure 1 :

 1 Figure 1: Training and testing with distribution mismatch. We consider three domains and represent λ s and λ t in the 3-dimensional simplex. Training with weights in (a) and testing with weights in (c) is supervised multi-source domain adaptation to domain 2 (d2), while (b)-(c) is the unsupervised version, with no training data from d2; training with weights in (a) and testing with weights in (d) is multi-domain learning, also illustrated with settings (a)-(e) (training domain d1 is not seen in test), and (b)-(d) (test domain d2 is unseen in training).

•

  our implementations of dynamic sampling proposals from the literature: Curriculum Learning (CL) of Zhang et al. (2019) and Differential Data Selection (DDS) of Wang et al. (2020b) (see below); All models use embeddings and hidden layers of dimension 512. Transformer models contain 8 attention heads in each of the 6+6 layers; the inner feedforward layer contains 2048 cells. Training lasts for 200K iterations, with batches of 12,288 tokens, Adam with parameters β 1 = 0.9, β 2 = 0.98, Noam decay (warmup steps = 4000), and a dropout rate of 0.1 in all layers.

  For a given line, all the columns correspond to the same multi-domain system. ( * ) MDAC is significantly better than Mixed-α with p < 0.05. ( * * ) MDAC is significantly better than DDS with p < 0.05.

Figure 2 :

 2 Figure 2: Evolution of the sampling distribution during training.

Table 1 :

 1 Corpora statistics: number of parallel lines (×10 3 ) and proportion in the training domain mixture (exluding N E W S), number English and French tokens (×10 6 ), types and uniq types (×10 3 ): the latter are types that only appear in a given domain.

	lines	2609 (0.68) 501 (0.13) 190 (0.05) 270 (0.07) 160 (0.04) 130 (0.03)	260 (0)
	# tokens	133 / 154	17.1 / 19.6	6.3 / 7.3	3.6 / 4.6	3.6 / 4.0	3.2 / 3.4	7.8 / 9.2
	# types	771 / 720	52.7 / 63.1 92.3 / 94.7 75.8 / 91.4 61.5 / 73.3 22.4 / 10.5	-
	# uniq	700 / 640	20.2 / 23.7 42.9 / 40.1 44.7 / 55.7 20.7 / 25.6	7.1 / 2.1	-
	corpus (N E W S). Most corpora are available from
	the Opus website 3 . These corpora were dedupli-
	cated and tokenized with in-house tools; statistics
	are in Table 1. To reduce the number of types, we
	use Byte-Pair Encoding (Sennrich et al., 2016) with
	30,000 merge operations on a corpus containing
	all sentences in both languages.We randomly se-
	lect in each corpus a development and a test set
	of 1,000 lines and keep the rest for training. Val-
	idation sets are used to chose the best model ac-
	cording to the average BLEU score

Table 3 :

 3 Adapting to two domains. For a given line, non empty columns correspond to the pair of target domains. ( * ) MDAC is significantly better than DDS with p < 0.05.

		λ2)	39.5	-	-	-	50.1	-
	MDAC(λ0, λ2)	39.1	-	-	-51.8 *	-
	DDS(λ0, λ2)	-	60.8	53.3	-	-	-
	MDAC(λ0, λ2)	-	61.9 *	54.5 *	-	-	-
	DDS(λ0, λ2)	-	-	-	37.9	-	91.3
	MDAC(λ0, λ2)	-	-	-	36.9	-	90.4
	domain d =	M E D	L A W	B A N K	T A L K	I T	R E L	mean
	Mixed-0	38.6	59.3	53.7	37.3	51.0	90.4	55.1
	Mixed-0.25	38.9	59.6	53.3	37.6	50.5	90.6	55.1
	Mixed-0.5	39.0	60.2	52.5	38.5	51.9	90.3	55.4
	Mixed-0.75	39.4	59.9	51.9	38.8	50.0	87.6	54.6
	Mixed-1	40.3	59.5	49.8	36.4	49.0	80.0	52.5
	DDS(λ0, λ0)	40.1	56.9	50.7	37.4	46.8	92.0	54.0
	MDAC(λ0, λ0)	38.5 60.3 * *	54.4 *		37.3 51.3 * *	91.4 * 55.5 * *
	DDS(λ1, λ0)	40.6	55.5	48.0	36.2	46.9	60.1	47.9
	MDAC(λ1, λ0)	40.2 59.3

* * 51.0 * * 36.9 * * 48.6 * * 80.7 * * 52.8 * *

Table 4 :

 4 Multidomain adaptation.

Table 5 :

 5 Unseen sample to schedule the presentation of training data so as to start with the easiest examples and end with the hardest. Various ways to automate CL in the framework of multi-armed bandits are explored in[START_REF] Graves | Automated curriculum learning for neural networks[END_REF], which has been an inspiration for our implementation. While the initial aim was primarily to improve and speed up training, CL has also proven useful for multidomain and multilingual MT, based on alternative definitions of "easiness". For instance,[START_REF] Zhang | Curriculum learning for domain adaptation in neural machine translation[END_REF] study supervised DA and propose a curriculum approach which progressively augments the training data: early stages only use in-data, while less relevant 10 data are introduced in later stages. This is opposite to the policy of van der Wees et al. (2017), whose gradual fine-tuning progressively increases the focus on in-domain data.

	38.3	60.1	50.3	35.8 49.1	90.1	53.9

* , λ d ) domain adaptation (left) and automatic clustering adaptation (right). For a given line, each column corresponds to one distinct system. ( * ) MDAC is significantly better than DDS. a

In this paper, multi-source DA means having multiple domains to adapt from; this setting differs from multi-source translation, where several source languages are considered.

https://ufal.mff.cuni.cz/ufal_medical_ corpus. We only use the in-domain (medical) subcorpora: PATR, EMEA, CESTA, ECDC.

http://opus.nlpl.eu

We use truecasing and sacrebleu[START_REF] Post | A call for clarity in reporting BLEU scores[END_REF].

https://github.com/neulab/compare-mt

We however omit domain-specific systems trained only with the corresponding subset of the data, which are always inferior to the mix-domain strategy[START_REF] Britz | Effective domain mixing for neural machine translation[END_REF].

https://github.com/OpenNMT/OpenNMT-tf

Domain distance is computed with Lewis-Moore scores (based on the cross-entropy of in-domain LM).