Mathilde Monvoisin
email: mathilde.monvoisin@irt-jules-verne.fr

Philippe Leray
email: philippe.leray@univ-nantes.fr

Mathieu Ritou
email: mathieu.ritou@univ-nantes.fr

Bruno Castanier

Marko Cepin

David Bigaud

Christophe Berenguer

Unsupervised co-training of Bayesian networks for the diagnosis of machining spindle

Keywords: Industry 4.0, condition-based maintenance, process monitoring manufacturing, machining, tool breakage, Bayesian networks, unsupervised learning, co-training

The objective of Smart Manufacturing is to improve productivity and competitiveness in industry, based on in-process data. Indeed, failures can stop the production for a couple of days and generate costs of non-quality. Failures in industry can either damage the machine or the product being produced. In both cases, the earlier the failure is detected, the lower the impact on production. Thus, monitoring both the process and the machine condition is interesting, due to their potential interactions. Besides, the diagnosis of the nature of the incident is also important, in order to react adequately as fast as possible. It requires reliable, explainable and understandable models such as Bayesian networks for performing tasks like condition prediction. Bayesian networks can be learned with incomplete data and in a supervised or unsupervised way, which is very useful because the collect of labelled data is costly and sometimes impossible, especially in industry where problems are, moreover, very rare. In this paper, we propose a generic architecture based on two Bayesian networks and a collaborative learning strategy that improves the condition monitoring of rotating machines in unsupervised context by using information gathered from process monitoring.

Introduction

Smart manufacturing is a promising research area for the improvement of productivity and competitiveness in industry, by exploiting the manufacturing digital data [START_REF] Tao | Datadriven smart manufacturing[END_REF].

Rotating machines are widely used in many industries. This paper focuses on machine-tools and especially on the spindle, which is a very efficient and critical rotating machine. The spindle is the machine component that holds the cutting tool (in the same ways that a drill holds a drill bit). Ball bearings are used for the rotational guidance of the spindle shaft.

Condition-based maintenance of rotating machines is usually based on periodic measurements of vibration [START_REF] Ebersbach | Expert system development for vibration analysis in machine condition monitoring[END_REF][START_REF] Mathew | The condition monitoring of rolling element bearings using vibration analysis[END_REF]. Oil analysis [START_REF] Craig | Advanced condition monitoring of tapered roller bearings, part 1[END_REF] and thermal measures can also be performed. Practically, the analysis of vibration signal of rotating machines can be performed through global criteria that evaluates the vibration level, such as the root mean square velocity (V rms ). It is defined in the ISO standard 10816-3 (ISO, 2009) and widely used for condition-based maintenance. ISO standard 17243-1 (ISO, 2014) recommend the use of both velocity V rms and acceleration A rms for the monitoring of spindle condition. For more accurate monitoring, the frequency spectrum of vibration can be analyzed. Particularly, BPFO and BPFI reveal local defect on the outer race or the inner race of a bearing respectively. Besides, distributed defect can also be monitored, such as race waviness [START_REF] Vafaei | Vibration monitoring of high speed spindles using spectral analysis techniques[END_REF] or worsening of surface roughness de [START_REF] De Castelbajac | Monitoring of distributed defects on hsm spindle bearings[END_REF]. [START_REF] Randall | Rolling element bearing diagnostics-a tutorial[END_REF] present the state of the art for the diagnostic of rolling bearing by analysis of vibration signal. [START_REF] Jimenez-Cortadi | Predictive maintenance on the machining process and machine tool[END_REF] review articles from the past ten years treating about predictive maintenance. The goal is to predict the remaining useful lifetime, which is very challenging. However, this paper focuses on condition-based maintenance, since machine-tool is complex usecase. Notably, spindle condition monitoring should be performed during idle rotation of the spindle, to avoid misinterpretation from vibration data collected during cutting.

Process monitoring consists in analyzing the machining in-process signals, in order to detect tool breakage, or predict the surface quality of the machined workpiece. Most of researches about condition monitoring in machining focus on tool wear monitoring [START_REF] Teti | Advanced monitoring of machining operations[END_REF]. Other papers consider the impact of chatter on the surface quality of the machined surface [START_REF] Peigne | Impact of the cutting dynamics of small radial immersion milling operations on machined surface roughness[END_REF]. Literature deals with the prediction and monitoring of such unstable vibration, by vibration or force signals [START_REF] Munoa | Chatter suppression techniques in metal cutting[END_REF][START_REF] Kuljanic | Development of an intelligent multisensor chatter detection system in milling[END_REF]. Tool condition monitoring can rely on machine learning; generally in a supervised way, to build a predictive model [START_REF] Kim | Smart machining process using machine learning: A review and perspective on machining industry[END_REF]. For example, a Probabilistic Neural Network can classify broken tools and good ones [START_REF] Huang | A pnn self-learning tool breakage detection system in end milling operations[END_REF].

In a more general framework, Bayesian networks (BNs) are probabilistic graphical models that are suitable for maintenance or in-process monitoring in the manufacturing industry [START_REF] Weber | Overview on bayesian networks applications for dependability, risk analysis and maintenance areas[END_REF]. Indeed, their ability to model complex systems under uncertainty, to learn these models with unlabelled, incomplete, unbalanced or small datasets [START_REF] Benferhat | Belief Graphical Models for Uncertainty Representation and Reasoning[END_REF]) is well known. [START_REF] Tidriri | Bridging data-driven and model-based approaches for process fault diagnosis and health monitoring: A review of researches and future challenges[END_REF] review the main data-driven and model-based approaches for process fault diagnosis and compare their respective advantages and their drawbacks. Bayesian networks offer a great reliability to manage uncertainty. For example, [START_REF] Atoui | A single Bayesian network classifier for monitoring with unknown classes[END_REF] propose to use BNs with a Gaussian mixture model for supervised datadriven process monitoring. Supervised approaches are often proposed, whereas the collection of labelled in-process data is generally impossible in industry. Moreover, the potential interaction between the machine condition and the process condition is not studied in literature, and even less with uncertainties consideration.

The objective of this paper is to propose a generic architecture to improve the conditionbased maintenance of rotating machines in unsupervised contexts by exploiting information from in-process monitoring. In this way, the aim is to detect machine failures caused by the manufacturing process; and to explain and understand which type and intensity of incident damage the spindle, and which evolution of spindle condition monitoring criterion reveals those damages.

Our two main contributions are: (1) The proposal of a generic architecture based on BNs that models each monitoring task independently and allows to define strongly critical values for the continuous monitoring criteria provided as inputs of our architecture.

(2) Each of our models is composed of a two-level structure that is able to take into account the different possible values of the criticality threshold used in the definition of these monitoring criteria.

The next section details the monitoring criteria used, as well as the data preprocessing step. Section 3 describes the model architecture and its different components, such as the learning strategy. In section 4, experiments are performed to validate the proposed architecture. Finally, section 5 draws the conclusions and perspectives of this work.

Monitoring criteria and data preprocessing

The choice of vibration sensor, of the signal preprocessing and of the monitoring criteria was performed in previous works by [START_REF] Godreau | Continuous improvement of HSM process by data mining[END_REF].

The raw in-process monitoring criteria are A rms (t), V rms (t), Nh(t) and Ub(t). V rms (t) and A rms (t) are respectively the root mean square of the vibration velocity and acceleration, revealing the current vibration level. Nh(t) enables chatter detection, through the asynchronous vibration level, obtained by order tracking in the vibration spectrum. Lastly Ub(t) detects tool breakage, by the amplitude of the contribution at the spindle frequency that estimates the tool balancing.

These raw criteria are obtained by realtime signal processing (from 25kHz vibration signals) and collected at a frequency of 10 Hz (864000 values per day). To be compared with daily evaluation of spindle condition, the in-process monitoring criteria must be daily aggregated. Since incidents occurring during machining can be very brief, they cannot be detected by averaging over all the monitoring values measured during a day. Besides, only severe vibrations, exceeding a certain threshold, can damage the spindle. Therefore [START_REF] Godreau | Continuous improvement of HSM process by data mining[END_REF] proposed to aggregate these raw in-process criteria as described in equation 1 where τ is the criticality threshold associated to a generic criterion X(t). High τ values reveal short and violent events, when small values highlight long and moderate events.

X τ = t=864000 t=1 max(X(t) -τ, 0) (1)
In this way, in-process critical events are daily monitored, by computing the advanced criteria A τ rms , V τ rms , Nh τ , Ub τ . [START_REF] Godreau | Continuous improvement of HSM process by data mining[END_REF] proposed values for the criticality threshold of these criteria.

Regarding the spindle condition monitoring: a vibration signature is performed once a day during idle rotation of the spindle, to evaluate its condition. Faulty bearing induced vibrations are evaluated in the frequency spectrum: the amplitudes at the fundamental train frequency F T F , at the Ball Pass Frequency on respectively the Outer and Inner rings BP F O and BP F I criteria; as well as the shaft balancing (by the vibration amplitude at the spindle rotation frequency), noted 1X.

Model

Our architecture, depicted in Figure 1, proposes to model independently the spindle condition and the in-process monitoring with two Bayesian networks. The aim is to predict the state of the system from the observation of the advanced monitoring criteria defined in the previous section.

As a preamble, all the advanced monitoring criteria are discretized in a three-level scale {OK, degraded, KO} (section 3.1).

As seen in the previous section, critical thresholds τ can be defined for these monitoring criteria. We then propose to structure the in-process monitoring in a two steps procedure. The first step is a local aggregation with several critical thresholds τ for each criterion, in order to provide a local decision based on each criterion (section 3.2). The second step is a global aggregation that takes into account all the local decisions in order to predict the state of the system (section 3.3). Finally, we propose collaborative strategies to jointly learn the two models that are respectively dedicated to the spindle condition and to the in-process monitoring (section 3.4).

Discretization

The advanced monitoring criteria defined in section 2 characterize critical vibrations that could damaged the spindle. They can be classified as weakly critical (values near 0), moderately critical (medium values) or strongly critical (higher values). We thus propose to discretize in a three-level scale {OK, degraded, KO} each continuous criterion X τ (noting its discretized counterpart X τ d ). Considering the discretization step in the model itself is a usual practice when dealing with probabilistic graphical model (see for instance [START_REF] Mabrouk | Multivariate cluster-based discretization for Bayesian network structure learning[END_REF]). We then propose in a similar way to deal with this discretization step in the form of a Gaussian mixture model where the parameters are learnt in an unsupervised way by Expectation-Maximization (EM) algorithm [START_REF] Dempster | Maximum likelihood from incomplete data via the EM algorithm[END_REF].

Local aggregation

Each in-process monitoring criterion X provides us a set of (continuous) advanced monitoring criteria X τ and their discretized counterparts X τ d . The objective of this step is to provide a local decision X d for each criterion, considering several critical thresholds X τ d simultaneously.

We propose in this work to perform this aggregation in two different ways. The first one is the very naive Max operator that describes the fact that the local decision is the worst one returned by this criterion with all its possible thresholds. The second one is more evolved and consists in using a probabilistic aggregation such as the NoisyMax model [START_REF] Srinivas | A generalization of the noisyor model[END_REF]. With this kind of aggregation, the dependency between the decision and the inputs is no more deterministic, and the influence of each input on the output is no more constant.

Global aggregation

The objective of the global aggregation is to provide a global decision for one submodel based on the local decision provided by the different advanced monitoring criteria (e.g. A rms d , V rms d , Nh d and Ub d for the inprocess monitoring sub-model).

As for local decision, we propose to deal with this aggregation with a deterministic Max or a probabilistic NoisyMax function.

This global aggregation is then simplified from the initial three-level scale {OK, degraded, KO} to a binary decision {NOC, KO} for the in-process condition or the spindle condition, where NOC stands for Nominal Operating Condition {OK, degraded}.

Learning strategy

A classical approach of this condition-based maintenance issue would be to learn the spindle condition model in a supervised way, with the ground truth of daily degradation of the studied machine part. In our industrial application, the collect of labelled data is impossible, but we can access to the data collected during the in-process monitoring task.

We propose to apply a collaborative strat-egy in order to learn jointly, and in an unsupervised way, the two models respectively dedicated to the spindle condition and to the process monitoring.

In the following scenarios, the learning of the Bayesian networks is limited to the parameter learning in presence of incomplete data by using the EM algorithm.

We consider two different collaborative strategies for this learning task. The simplest naive strategy consists in learning simultaneously both models, by adding an equality constraint between the global decisions of the two models.

The second strategy is the application of the recursive unsupervised co-training strategy, proposed by [START_REF] Monvoisin | Unsupervised co-training of Bayesian networks for condition prediction[END_REF]. It is an interesting solution to jointly learn two models in an unsupervised way. The idea is to firstly learn one model in an unsupervised way (for instance the model dedicated to spindle condition monitoring) and then to learn the second one (process monitoring) enriched with the results of the first one, obtained by probabilistic inference. In the recursive strategy, this step is repeated several times until convergence.

Experiments

Experimental context

Those experiments are realized on real industrial data that contains 436 days of aggregated and unlabelled data, collected over three spindle lifetimes (a year and a half).

The results obtained by [START_REF] Godreau | Continuous improvement of HSM process by data mining[END_REF] which were afterwards confirmed by an expert are considered as the ground truth. The optimal thresholds identified in that study for each criterion (by genetic algorithm) are also used in some experiments where we will short-circuit the local aggregation step.

Incidents are rare in machining: from respectively the in-process monitoring and the spindle condition monitoring datasets only 7 (1.6%) and 5 (1.1%) events were found in that study as events that should have seriously damaged the spindle.

Another issue in this real application is the fact that the two diagnosis do not coincide perfectly: only 3 of the 5 spindle damages are explained by the 7 events detected during the in-process monitoring.

The implementation of the models and strategies described in section 3 has been performed with our library dedicated to Probabilistic Graphical Models (PILGRIM) with the help of ProBT a library. Parameters used during learning were a stopping criterion equals to 10 -4 for the EM algorithm and n step = 30 for the recursive strategy.

The performances of the spindle condition monitoring model (balanced accuracy, specificity, sensitivity) are estimated with five separate iterations of two-fold cross-validation as proposed in [START_REF] Dietterich | Approximate statistical tests for comparing supervised classification learning algorithms[END_REF]. Some preliminary experiments were realized in order to take into account the very high imbalance ratio by applying an over-sampling strategy. It did not improve the performance. So we do not use over-sampling in the following experiments.

Discretization

Fig. 2 shows an example of the discretization process proposed in section 3.1. The grey area is the histogram of the V 120 rms values. The blue and orange lines are the cut points between the three states OK, degraded and KO, a https://www.probayes.com/ determined as proposed by [START_REF] Mabrouk | Multivariate cluster-based discretization for Bayesian network structure learning[END_REF]. Weakly critical vibrations are considered for this variable when V 120 rms < 429, or strongly critical when V 120 rms > 1735.

Learning strategies

The first element we want to study is the learning strategy. As defined in section 3.4, we here compare the Naive strategy (where both models are jointly learnt with an equality constraint) with the Recursive strategy.

In this experiment, the local aggregation step is bypassed by only using the optimal critical threshold identified by [START_REF] Godreau | Continuous improvement of HSM process by data mining[END_REF] for each criterion. The global aggregation is performed by Max aggregation.

Table 1 shows the performances of both learning strategies (line (1) for the naive learning and (2) for the recursive co-training strategy). As seen in this table, the recursive co-training strategy slightly improves the balanced accuracy, with a perfect sensitivity (i.e. the proportion of KO correctly detected among the real KO), but a smaller specificity (i.e. the proportion of NOC correctly detected among the real NOC).

Global aggregation

The second choice we want to study concerns the global aggregation, by comparing the deterministic Max and the probabilistic NoisyMax functions defined in section 3.3.

In this experiment, the local aggregation step is bypassed again. Based on the previous result, the learning strategy is the recursive one.

Table 1 shows that the performances when applying respectively the Max aggregation (line (2)) and the NoisyMax one (line (3)) are almost the same, with only a very small decrease of the specifity and a non prohibitive increase of the computation time due to the increase of the number of parameters in NoisyMax model. Using the NoisyMax global aggregation can help by weighting differently the inputs of this aggregation and providing meaningful interpretation of the relevant inputs. But in this experimentation, it does not help in improving the balanced accuracy.

Local aggregation

After dealing with the learning strategy and the global aggregation, we are interested here in the choice of the local aggregation for the process monitoring model.

As defined in section 3.2, we here compare the deterministic Max and the probabilistic NoisyMax functions for both local and global aggregation. The learning strategy is the recursive one, based on the previous experiments.

Table 1 shows the performances when applying respectively the Max function in the local and global aggregation steps (line (4)), the local Max and the global NoisyMax (line (5)), and the local NoisyMax and the global Max (line ( 6)).

Let us recall here that the previous results were bypassing this local aggregation step by using the optimal inputs identified for each criterion in a previous study. We can see here that working with the whole set of inputs provides us more and more complex models in these three experiments, with more and more parameters, which explains the increase of the computational time. We can also see that we are able to get a comparable balanced accuracy when we are no more bypassing the local aggregation step. This shows us that we can learn one full model without taking into account any previous knowledge.

Conclusion

In this paper, we were interested in the unsupervised co-training of Bayesian networks for the diagnosis of machining spindle by condition monitoring, with the help of data collected during in-process monitoring. We have proposed a generic architecture that models each monitoring task independently and applies a collaborative learning strategy that takes into account both tasks. Our architecture is also able to take into account the fact Our experiments, based on a real and imbalanced dataset, were able to identify the same degradations than [START_REF] Godreau | Continuous improvement of HSM process by data mining[END_REF]. It led to surprising results where the simple aggregation schema (Max) in both local and global aggregation steps provided the best results, without considering any external knowledge (i.e. without using the Optimal local aggregation). Complexifying the aggregation with NoisyMax led to more complex models that are not useful with our limited number of data, which is conform with the principle of parsimony in Machine Learning. A paired t-test between the different models investigated, did not exhibit any significant difference between them. In another way, these more complex models could estimate which input or sensor are interesting. Based on this observation, one of our perspectives to improve this work is to add a feature selection procedure when using the NoisyMax aggregation in order to select the relevant monitoring criteria in the model, to decrease the number of parameters and to converge to similar results than the ones obtained with the Optimal models. Another perspective is related to the improvement of the unsupervised co-training procedure as described in [START_REF] Monvoisin | Unsupervised co-training of Bayesian networks for condition prediction[END_REF].

One last perspective concerns the very specific dataset we have used for our experiments. This dataset was initially unlabelled and the results obtained by [START_REF] Godreau | Continuous improvement of HSM process by data mining[END_REF], afterwards confirmed by an expert, were considered as the ground truth. This means that all the existing points labelled as KO have been confirmed as real ones in this dataset; but not that all data that are currently labelled as NOC are real ones. The sensitivity obtained by our system is 100%, which means that we are also able to detect the same KO. New and unknown KO could be discovered soon in the dataset by our improved cotraining method.

  Our generic architecture based on two Bayesian networks that models the spindle condition monitoring task (right) and the in-process monitoring one (left) independently, including a discretization step for the daily-aggregated monitoring criteria (blue nodes). Circles denote continuous random variables and (rounded) squares discrete ones. PC and MC stand respectively for Process Condition and Machine Condition.

Fig. 2 .

 2 Fig. 2. Distribution (on a logarithm scale) of the values of the continuous criterion V 120 rms and cut-points estimated after the discretization step.

Table 1 .

 1 Performances (Balanced accuracy, Specificity and Sensitivity in %) and computational time (s) obtained in different experimental setup (choice for the three parts of our architecture: local aggregation, global aggregation and learning strategy, as proposed in section 3)

		Local agg.	Global agg.	Learn. strat.	Bal. acc.	Spec.	Sens.	Time
	(1)	Optimal	Max	Naive	84.2	90.0	78.3	0.7
	(2)	Optimal	Max	Recursive	85.6	71.5	100.0	0.6
	(3)	Optimal	NoisyMax	Recursive	85.5	71.0	100.0	3
	(4)	Max	Max	Recursive	85.2	70.4	100.0	2
	(5)	Max	NoisyMax	Recursive	85.1	70.1	100.0	59.6
	(6)	NoisyMax	Max	Recursive	83.9	68.3	100.0	1945
	that several possible critical thresholds can be
	used in the daily computation of advanced
	monitoring criteria.

Acknowledgement

Thanks are addressed to the IRT Jules Verne, French Institute in Research and Technology in Advanced Manufacturing for the PhD PERFORM program. The authors would also like to thank people who laid the foundations of this project with P. Leray and M. Ritou, in particular V. Godreau, A.M. Idrissi and G. Ferrand, and the industrial partners for supplying the experimental data.


T Abadie

C Xuereb

D Legendre

J Aubin

Mixing and recirculation characteristics of gas-liquid Taylor flow in microreactors

Keywords: Gas-liquid Taylor flow, Microchannel, Microreactor, Mixing, Recirculation, CFD, VOF

. In addition, increasing the aspect ratio of the channel leads to a slight decrease in recirculating volumes but also a significant increase in recirculation times.

Introduction

Over the last decade, micro reaction technology has become of much interest to both academics and the process industries for the intensification of chemical processes. Taylor or slug flow is a commonly encountered flow regime for gas-liquid microchannel flows and has the advantage of providing high interfacial area and good liquid mixing in the liquid slug, thereby enhancing transport processes. These features of microreactors are particularly interesting for fast and highly exothermic gas-liquid reactions -amongst other applications and allow an increase in reaction performance whilst working under safe operating conditions [START_REF] Hessel | Gas-liquid and gas-liquid-solid microstructured reactors: contacting principles and applications[END_REF][START_REF] Kashid | Microstructured reactors for multiphase reactions: state of the art[END_REF]. A number of studies have focused on the understanding of hydrodynamics as well as heat and mass transfer enhancement in these flows (see reviews [START_REF] Kreutzer | Multiphase monolith reactors: chemical reaction engineering of segmented flows in microchannels[END_REF][START_REF] Shao | Flow regimes for adiabatic gas-liquid flow in microchannels[END_REF][START_REF] Gupta | Taylor flow in microchannels: a review of experimental and computational work[END_REF] but often independently. Indeed the transport efficiency appears to be closely related to the recirculation in the liquid phase, which depends on the operating conditions, fluid properties and reactor geometry.

A major feature of gas-liquid Taylor flow is the recirculation flow pattern generated in the liquid slug in the moving frame of reference as shown in Fig. 1 and detailed by [START_REF] Taylor | Deposition of a viscous fluid on the wall of a tube[END_REF]. The recirculating flow pattern is characterized by the position of the center of the circulation loop [x 0 , y 0 ], and the position of the streamline separating the recirculation zone and the liquid film at the channel wall [x 1 , y 1 ]. As the bubble velocity increases, both the loop center and the outer streamline of the recirculation zone move towards the center of the channel [START_REF] Thulasidas | Flow patterns in liquid slugs during bubble-train flow inside capillaries[END_REF]. This leads to a reduction of the volume of the recirculating zone and an increase in the volume of liquid in the film region until complete bypass flow occurs at U B ≥ U max . Studies on the recirculating flow in circular and square capillaries have been conducted previously [START_REF] Taylor | Deposition of a viscous fluid on the wall of a tube[END_REF][START_REF] Thulasidas | Flow patterns in liquid slugs during bubble-train flow inside capillaries[END_REF], however little information on the characteristics of the recirculation zone in rectangular channels is available. Recently, it has been shown theoretically [START_REF] Kececi | Recirculation time and liquid slug mass transfer in co-current upward and downward Taylor flow[END_REF] that the cross-sectional area occupied by the recirculation zone in the liquid slug is generally smaller for rectangular channels than for circular channels due to the increased film thickness in the channel corners.

The recirculation time in the liquid slug in Taylor flow is defined as the time required for an element of fluid to complete one revolution in the recirculating slug. This characteristic time is particularly relevant for transport processes Fig. 1 -Streamline pattern in the frame of reference moving with the bubble (Ca = 0.06 and Re = 16.8). The flow is directed from the left to the right. The positive part of the recirculating pattern of volume V 0 is delimited by the coordinates [x 0 , y 0 ] in the channel width and depth; the negative part of the recirculating zone of volume V 1 = V rc -V 0 is located between [x 0 , y 0 ] and [x 1 , y 1 ]; and the film region of volume V 2 = V slug -V rc is located between [x 1 , y 1 ] and the walls. occurring in the system, such as mass transfer between the bubble and slug or wall and slug, and heat transfer with the channel wall. The rate of flow recirculation can be calculated via the surface integration of the relative velocity profile across the microchannel cross-section, similarly to what is done in conventional stirred tanks to calculate circulation induced by the mechanical impeller [START_REF] Jaworski | An LDA study of the turbulent flow field in a baffled vessel agitated by an axial, down-pumping hydrofoil impeller[END_REF]. The recirculation flow rate through the microchannel can be divided into three parts (illustrated in Fig. 1):

a positive flow rate, Q 0 , in the main flow direction at central core of the microchannel occupying a volume V 0 with a cross-sectional area A 0 , a negative flow rate, Q 1 , with a volume V 1 and a crosssectional area A 1 that corresponds to the recirculating liquid in the slug, -a negative flow rate, Q 2 of volume V 2 in area A 2 , that is close to the channel wall and contributes to axial mixing between slugs instead of radial mixing within the slug.

The recirculation time is then defined as t rc = V rc /Q rc , where the recirculating volume V rc corresponds to the volume of liquid within the limit of the separating streamline and Q rc is the recirculation flow rate equal to Q 0 (and |Q 1 |). t rc can be made non dimensional by dividing it by the time taken for the bubble to travel a distance equal to the slug length:

rc = t rc /(L S /U B ).
In our previous work we have explored the effects of fluid properties, operating conditions and microchannel geometry on the size of Taylor bubbles [START_REF] Abadie | Hydrodynamics of gas-liquid Taylor flow in rectangular microchannels[END_REF] and the flow patterns in the liquid slug using micro Particle Image Velocimetry ( -PIV) [START_REF] Zaloha | Characteristics of liquid slugs in gas-liquid Taylor flow in microchannels[END_REF]. The objective of this work is to explore the effects of operating parameters (capillary, Ca, and Reynolds numbers, Re) and microchannel aspect ratio (˛= w/h = [1; 2.5; 4]) on the mixing and recirculation characteristics of the liquid slug in gas-liquid Taylor flow in microchannels. To do this, 3-dimensional VOF simulations of gas-liquid Taylor flow in microchannels have been performed. Using an approach that is analogous to the determination of circulation rate in stirred tanks, the recirculation rate in the liquid slug, as well as the size of the recirculating zone have been evaluated from the 3dimensional numerical data. An attempt has been made to relate these characteristics of the recirculating liquid slug to the enhanced transport phenomena observed in Taylor flow in microreactors. Finally, recommendations on the design and operation of microreactors employing Taylor flow are given.

Methodology

Theoretical developed velocity profile in rectangular capillaries

In the limit of infinite slug length, the analytical solution of the velocity profile in a cross-section of a rectangular capillary can be derived. The corresponding velocity profile in a rectangular capillary with a cross-section 2w × 2h is given by Eqs.

(1) and (2). This theoretical velocity profile has been used to evaluate the effects of aspect ratio and dimensionless bubble velocity on characteristic parameters of recirculation motion in gas-liquid Taylor flow: recirculating volumes and recirculation times. Under the assumption that the slug is long when compared with development lengths of velocity profile at the rear and the nose caps of the bubbles, the recirculation volume can be approximated by the cross-sectional area of the recirculating zone multiplied by the slug length:

V 0,1 ∼ A 0,1 × L S . u(x, y) = - 16c 1 w 2 3 ∞ n=1,3,...
-1 (n-1)/2 n 3 1 -cosh(ny/2w) cosh(nh/2w) cos nx 2w

(1)

U TP = - c 1 w 2 3 1 - 192 5 w h ∞ n=1,3,... 1 n 5 tanh nh 2w (2)
Numerical integration of the profiles (1) and (2) for aspect ratios ˛= [1 ; 2.5 ; 4] have been performed on fine uniform grids of 400 × 400, 1000 × 400 and 1600 × 400, respectively, to calculate the area A 0 where the velocity in the moving coordinate system is positive. Following this, an iterative integration of A 1 (cross-section of the negative part of the recirculating liquid) has been performed until the flow rates in the positive and negative parts of the recirculating zone are balanced. Integration of these profiles allows the prediction of recirculating volumes and recirculation times as a function of the aspect ratio and the dimensionless velocity W = (U B -U TP )/U B assuming long liquid slugs.

Numerical simulations

The numerical code used for this study is the JADIM code, which has been developed to simulate dispersed two-phase flows and used to simulate various multiphase flows systems [START_REF] Bonometti | An interface-capturing method for incompressible two-phase flows. validation and application to bubble dynamics[END_REF][START_REF] Dupont | Numerical simulation of static and sliding drop with contact angle hysteresis[END_REF]Sarrazin et al., 2006a,b;[START_REF] Abadie | Hydrodynamics of gas-liquid Taylor flow in rectangular microchannels[END_REF]. The interface capturing technique implemented in this code is the Volume of Fluid method (VOF), which consists of a Eulerian description of each phase on a fixed grid. Under the assumptions that (i) the fluids are Newtonian and incompressible, (ii) there is no mass transfer at the interface, (iii) the flow is isothermal and (iv) the surface tension is constant, the fluid flow can be described by the classical one fluid formulation of the Navier-Stokes equations:

∇ • U = 0 (3) ∂U ∂t + (U • ∇) U = -∇P + ∇ • + g + F ,s (4)
where is the viscous stress tensor, g is the acceleration due to gravity, F ,s = -(∇ • n) nı I is the capillary contribution whose calculation is described below, is the surface tension, n the normal to the interface, ı I is the Dirac distribution localizing the interface, and and are the local density and dynamic viscosity, respectively. The density and viscosity are deduced from the volume fraction of one phase (or color function) C by linear interpolation:

= C 1 + (1 -C) 2 , = C 1 + (1 -C) 2 (5)
where the volume fraction is C = 1 in cells filled with fluid 1, C = 0 in cells filled with fluid 2 and 0 < C < 1 in cells that are cut by the interface. Additionally, the transport equation of the color function is solved to capture the interface between the phases:

∂C ∂t + U • ∇C = 0 (6)
In many VOF methods employed to capture a fluid-fluid interface, a reconstruction technique is used to control the thickness of the interface. In the JADIM code, the location and thickness of the interface are both controlled by an accurate algorithm based on Flux-Corrected Transport schemes [START_REF] Zalesak | Fully multidimensional flux-corrected transport algorithms for fluids[END_REF][START_REF] Bonometti | An interface-capturing method for incompressible two-phase flows. validation and application to bubble dynamics[END_REF].

The equations are discretized on a staggered grid using a finite volume method and all spatial derivatives are approximated using second-order centered schemes. The time scheme used to compute the advective terms in the Navier-Stokes equations is a third-order Runge-Kutta type scheme, while the viscous stresses are solved using a semiimplicit Crank-Nicolson method. The incompressibility is ensured using a projection method, which consists in splitting the velocity field into two contributions: a rotational one, which gives a predicted velocity field calculated semiimplicitly, and a potential one, obtained from a pressure correction solution of a pseudo-Poisson equation, whose divergence is null.

In these capillary flows, care needs to be taken in the calculation of the capillary force F ,s . The numerical method used to solve the interfacial force is the Continuum Surface Force (CSF) proposed by [START_REF] Brackbill | A continuum method for modeling surface tension[END_REF]. The localization of the interface is available through a non-zero gradient of the volume fraction and the curvature is calculated from the volume fraction gradient. Thus, the surface force F ,s is transformed into a volume force F ,v by distributing its effects over grid points in the vicinity of the interface in a region that is a few cells thick:

F ,v = -∇ • (∇ C/|| ∇ C||) ∇ C.
Throughout this study, only a quarter of the channel is considered since the effects due to gravity in microchannels can be neglected (Bo ≤ 0.1). The bubble is placed in the domain and a pressure gradient between two periodic boundary conditions is imposed to induce the motion until the flow reaches a steady state. No-slip wall conditions are imposed on the planes x = w and y = h while a symmetry boundary condition is set at x = 0 and y = 0. The grids used consist of 32 × 32 × 256, 48 × 32 × 256 and 64 × 32 × 256 cells for microchannels with aspect ratios of 1, 2.5 and 4, respectively. In order to correctly resolve the thin liquid film between the bubble and the wall, these grids are refined at the walls using a geometric evolution that decreases cell size in the last 32 cells in the x and y directions. In the axial direction, a regular grid spacing is used and the length of the unit cell is 17 × h. The parameters considered in this study have been chosen so that the Laplace number La = Re/Ca = L D h / 2 is constant in all the geometries and equal to La = 280. The ratios of liquid to gas densities and viscosities have been set to L / G = 10 3 = L / G .

Results and discussion

In this section, the results concerning global observations about recirculation patterns, shapes of the recirculation zone in a cross-section of the slug, recirculating volumes, as well as recirculation times are given. As shown theoretically in circular channels [START_REF] Thulasidas | Flow patterns in liquid slugs during bubble-train flow inside capillaries[END_REF], the recirculation characteristics of infinite slugs can be related to one parameter, which is the slip velocity W = (U B -U TP )/U B . Many studies dealing with the dynamics of gas-liquid Taylor flow relate W to the capillary number [START_REF] Bretherton | The motion of long bubbles in tubes[END_REF][START_REF] Aussillous | Quick deposition of a fluid on the wall of a tube[END_REF]. In this study, the flow conditions were modified by varying the capillary number in the range 0.04 ≤ Ca B ≤ 0.5 and the corresponding Reynolds number in the range 11.2 ≤ Re B ≤ 140 whilst keeping the gas void fraction constant at 31%.

To study the effects of slug length on the recirculation characteristics, the geometry with an aspect ratio ˛= 2.5 has been considered and the capillary number has been set to Ca = 0.1 +/-0.0005 that leads to W = 0.3303 +/-0.0021. Five different cases have been simulated by varying slug length from 2.90 × h to 20.06 × h and the corresponding volume fraction of gas phase ranges from 50% to 12.5%.

Recirculation pattern and local velocity field

The recirculation patterns in a rectangular channel with aspect ratio ˛= 2.5 with varying capillary numbers are shown in It can be seen in Fig. 2 that the maximum value of the local Peclet number is found in the liquid film between the bubble rear cap and the wall. Whatever the capillary number, the maximum value of the local Peclet number in the liquid slug is found close to the wall, decreases to zero at the center of the recirculation loop and then increases whilst approaching the center of the channel. Furthermore, as the Reynolds number increases, the boundary layer close to the bubble caps (where the Peclet number is very low) decreases in size. It can also be seen that recirculation loops through both the channel width and depth are formed. However, as the capillary number increases, the recirculation pattern in the less confined direction (channel width) diminishes until the closed loop is entirely lost. These changes in the streamline patterns are mainly observable in the wake of the bubble where the flow is particularly disturbed and the streamlines do not follow the shape of the rear of the bubble in the less confined direction. These observations also depend strongly on the crosssectional shape of the channel and it can be seen in Fig. 3 that the detachment of the streamlines is enhanced when the aspect ratio increases. Indeed, although the recirculation loop in the square geometry seems qualitatively similar to that in a circular geometry, strong 3-dimensional effects of the rectangular geometry on the velocity field in the wake of the bubble are observed. In addition to the streamline pattern in the wake of the bubble, the same phenomenon is observed at the front cap of the bubble in the very wide shallow channel (aspect ratio ˛= 4). This is in agreement with previous works on semi-infinite bubbles at low Reynolds numbers [START_REF] Hazel | The steady propagation of a semi-infinite bubble into a tube of elliptical or rectangular cross-section[END_REF] in rectangular channels.

The effects of the capillary number and the aspect ratio on the dimensionless velocity at the channel centerline is illustrated in Fig. 4 and 5 respectively. It is seen that close to the bubble caps, the velocity at the channel centerline increases from the bubble velocity to the maximum velocity expected in the liquid slug that depends on the mean velocity (from Eqs. ( 1) and ( 2), U max ∼ [2.10 ; 1.92 ; 1.77] × U TP for channel cross-sections of aspect ratios ˛= [1 ; 2.5 ; 4]). It is also observed (Figs. 4 and5) that as the capillary number increases or the aspect ratio increases, there is a region in the wake of the bubble where the liquid velocity is greater than the maximum velocity of a laminar single-phase flow close to bubbles caps.

These observation suggest that the recirculation characteristics are not only governed by the dimensionless slip velocity W but also by the slug length and this will be discussed in detail in Sections 3.3 and 3.4. The influence of slug length is expected when the development lengths at the rear and front caps of the bubbles are non negligible in comparison with the slug length, as illustrated with the streamlines and local Peclet numbers in Figs. 2 and3.

Shape of the recirculation zone

The cross-sectional area occupied by the recirculation zones as approximated by the theoretical fully-developed velocity profile (see Section 2.1) has been calculated and is shown in Fig. 6 for the three aspect ratios: ˛= [1 ; 2.5 ; 4]. It can be seen that the recirculation areas follow the shape of the channel for low dimensionless velocities and reduce to circular and elliptical shapes as the velocity (or Ca) increases.

When the channel aspect ratio increases, the aspect ratio of the recirculation zone also increases. Furthermore, the center of the recirculation loop and the separating streamline both move faster towards the center of the channel in the less confined direction. As a result, the slug film thickness is greater in the width of the microchannel than in the height.

Mass transfer between the film and the recirculating regions will occur by molecular diffusion through the interface separating these regions. In a 2-dimensional representation through the channel cross-section, the perimeter that separates the recirculating and film regions can be used to understand transport processes between the recirculating and film regions. Fig. 7 shows the evolution of the perimeter that separates the recirculating and film regions in Fig. 6 as a function of the dimensionless slip velocity. To compare the three aspect ratios with circular channels, the perimeters are normalized by the perimeter of a circular channel with an equivalent hydraulic diameter D h . It is clearly seen that the perimeter that links the recirculation zone and the film region is greater in square channels than in circular channels. As the aspect ratio increases, the perimeter that separates the Fig. 4 -Dimensionless liquid velocity along the centerline of the channel between two bubbles for an aspect ratio of 2.5 and various capillary numbers. Legend: (-) U centerline /U B ; (--) U max /U B evaluated from U TP and Eqs. ( 1) and ( 2 recirculating and film regions also increases until a certain velocity W above which the opposite trend is observed. The velocity at which the transition occurs decreases as the aspect ratio increases. As a consequence, it is expected that when dealing with wall to slug mass transfer, channels with a high aspect ratio will typically provide better mass transfer than channels with square or circular cross-section. However, mass transfer from the recirculation zone to the film will also depend on the time needed for a fluid element to travel the slug distance and more generally, mixing efficiency will Fig. 5 -Dimensionless liquid velocity along the centerline of the channel between two bubbles for Ca = 0.1 and for aspect ratios 1, 2.5 and 4. Legend: (-) U centerline /U B ; (--) U max /U B evaluated from U TP and Eqs. ( 1) and ( 2

Characteristic recirculating volumes

A characteristic parameter of the recirculation motion in Taylor flow is the volume of fluid recirculating in the slug. Unlike in cylindrical tubes where the volumes of positive flow and negative flow of the recirculation region are equal, it is not the case in square and rectangular ducts. Indeed, in the channel cross-section the ratio of total recirculation area over positive recirculation area increases with the aspect ratio and tends to √ 3 for the asymptotic behavior of an infinite aspect ratio. In addition, for a fixed geometry, the ratio of recirculation areas also varies with the dimensionless velocity W.

The size of the recirculation area at the center of the slugs for channel aspect ratios [1 ; 2.5 ; 4] are shown in Fig. 8. A rc (cross-section of the recirculating liquid) and A 0 (cross-section of the positive part of the recirculating liquid) are normalized by the channel cross-section area A ch . The results integrated from the analytical velocity profile corresponding to infinite slug length and the results obtained from the numerical simulations with finite slug length are compared and very good agreement between the computational and theoretical results is observed, as long as the liquid velocity profile is fullydeveloped at the center of the slug. The computations that show a slight deviation from the theory correspond to cases where the velocity profile is not entirely fully-developed at the center of the slug. For short slugs, the laminar velocity profile is not fully-developed and the region where the flow is disturbed by the bubble is not negligible; therefore the recirculating volumes and times cannot be approximated by the theoretical velocity profile. In such cases, the recirculating volume is obtained from direct numerical simulations by integrating small slices of recirculating zones along the slug.

Fig. 9 shows the influence of the slug length on the recirculating volume V rc and the positive part of the recirculating volume V 0 normalized by the slug volume V slug . In every simulation, care has been taken that the bubble is long enough and the dimensionless bubble velocity is constant whatever the slug and bubble lengths. From this data, it is clear that the regions close to bubble caps do not affect the ratio of recirculating volume over slug volume significantly. Fig. 10 shows the ratio of the recirculating volumes V 0 and V rc over total slug volume V slug as a function of the capillary number for the three rectangular microchannel

V 0 /V slug ; (open symbols) V rc /V slug .
geometries. The recirculating volume for the infinite aspect ratio and axisymmetrical cases are also given using the relation W = f(Ca) from [START_REF] Aussillous | Quick deposition of a fluid on the wall of a tube[END_REF] for the tubes and [START_REF] Abadie | Hydrodynamics of gas-liquid Taylor flow in rectangular microchannels[END_REF] for 2D cases; both assume long liquid slug in order to neglect the changes in the flow close to the bubble caps. It should also be kept in mind that axisymmetrical and planar cases are given as an indication only since the relations W = f(Ca) for these cases are valid for negligible inertia, which is not the case in the present work since Re ≥ 10. Experimental data obtained by -PIV in [START_REF] Zaloha | Characteristics of liquid slugs in gas-liquid Taylor flow in microchannels[END_REF] (obtained for a Laplace number La = Re/Ca ∼ 1000) are also presented. As expected, the recirculating volumes are close to the slug volume at low capillary numbers (Ca B ∼ 10 -3 ) due to the fact that the bypass flow is negligible compared Fig. 10 -Dimensionless recirculating volume as a function of the capillary number. Legend: (-) V 0 /V slug , infinite slug length assumption; (--) V rc /V slug , infinite slug length assumption; (gray) circular channels; (light gray) 2D channels; (filled symbols) V 0 /V slug ; (open symbols) V rc /V slug ; (⊳ , black)˛= 2.5, experimental [START_REF] Zaloha | Characteristics of liquid slugs in gas-liquid Taylor flow in microchannels[END_REF]; ( , red)˛= 1, numerical; (⊲ , black)˛= 2.5, numerical; (♦ , blue)˛= 4, numerical. (For interpretation of the references to color in this figure legend, the reader is referred to the web version of the article.) with the total flow rate. In other words, almost all the liquid contributes to the recirculation motion. Then, as the capillary number increases, the recirculating volume significantly decreases in all the geometries. Although the volume of the positive recirculating flow is hardly dependent on the geometry, the total recirculating volume decreases as the aspect ratio increases.

3.4.

Characteristic recirculation times

The characteristic recirculation times 0 for the positive recirculation flow, 1 for the negative flow and rc for the total recirculation as a function of the dimensionless velocity are plotted in Fig. 11 in comparison with the infinite slug length assumption. The effects of slug length on positive ( 0 ) and total ( rc ) recirculation times are represented in Fig. 12 and the recirculation times as a function of the capillary number are then plotted in Fig. 13.

Although the recirculation areas show very good agreement with the theoretical predictions, the recirculation times rc are greater than those determined using the infinite slug assumption (Fig. 11), and mainly for aspect ratios 1 and 2.5 at low capillary numbers. Indeed, close to bubble caps where the flow is disturbed, the axial flow rate is reduced and the recirculation time is increased. However, as the capillary number or the aspect ratio increase, the velocity in the wake of the bubble increases (see Figs. 4 and5) and balances the lower velocities in the vicinity of the interface. It can be seen in Fig. 11 that the recirculation time for the channel with aspect ratio of 4 (represented with diamonds in blue) is lower in our calculations compared with the infinite slug assumption. The trends observed for all geometries are qualitatively similar and an increase in the aspect ratio leads to an increase in both the positive and negative recirculation times and thus in the total recirculation time until divergence when complete by-pass flow is attained.

Although the recirculating volume only depends on the dimensionless velocity W, it is shown in Fig. 12(a) that the dimensionless recirculating times vary with the slug length. Indeed, whilst the disturbed flow around the bubble does not change the recirculating volume, it does cause an increase in the recirculation time when the slug length decreases, as mentioned in the previous paragraph. As expected, the recirculation times tend towards the value evaluated from the theoretical laminar single-phase flow representing the case of an infinite slug length. For the shortest slug considered here, the total recirculation time is 20% longer than the asymptotic value. [START_REF] Leung | Heat transfer in well-characterised Taylor flow[END_REF] studied the effect of slug length on the Nusselt number with a constant heat flux on the wall of a circular tube and found an increase in the Nusselt number as the slug length decreases. Indeed, fluid elements will refresh the interface separating the recirculating and film regions more frequently in short slugs. The evolution of the mean time for a fluid element to complete a cycle can be represented in a dimensional form to illustrate this decrease in recirculation time as the slug length decreases. It is seen in Fig. 12(b) that although the recirculation time is slightly underestimated with the infinite slug length assumption due to the increased time taken for a fluid element to travel the development length close to bubble caps, the recirculation time decreases as the slug length decreases.

In the cases presenting a recirculation motion in the liquid slug, the bubble velocity is lower than the liquid velocity at the centerline of the channel. Thus, around the bubble caps, there is a region where the velocity on the channel centerline decreases to the bubble velocity, as shown in Figs. 4 and5. This explains the increase in recirculation time when the slug length decreases since the ratio of development length over slug length increases.

Here, it has been observed that in the development region close to the bubble caps, the centerline velocity can be greater than the expected velocity on the axis (e.g. U max = 1.92 × U TP for the 2.5 aspect ratio channel, Fig. 4 for Ca ≥ 0.1) before stabilizing if the slug is long enough. This local increase in the velocity can lead to a decrease in the recirculation time as the slug length decreases. This phenomenon is accentuated as the aspect ratio increases and it can be seen in Fig. 11 that recirculation times evaluated in the numerical simulations are lower than that for an infinite slug length with the aspect ratio 4 and W > 0.3. Fig. 13 reports the evolution of the recirculation times as a function of the capillary number. It is shown that an increase in the capillary number also leads to an increase in the recirculation time for all the geometries considered. Thus, when the capillary number increases, the volume involved in the recirculation motion is reduced but the dimensionless time needed for a fluid element to travel the slug length is increased. It is also observed that the 2D theory is not a good approximation for high aspect ratios. In addition, it can be seen that even if the recirculating volume shows a weak dependency on the channel geometry, the effect observed on the recirculation time is more significant. Indeed, at fixed operating conditions (Ca, Re and L S ), an increase in the aspect ratio leads to significantly longer recirculation times.

Concluding remarks

Direct numerical simulations of gas-liquid Taylor flow in microchannels have been performed. The effects of channel cross-section and operating conditions on mixing and recirculation in the liquid slug have been investigated. The results of the numerical simulations have been compared with theoretical approximations for the case of infinite slug length in terms of cross-sectional recirculation area, recirculating volume and recirculation times. It has been shown that as long as the slug is long, the recirculation area in the slug can be well predicted using the fully-developed velocity profile. Whatever the geometry, increasing the capillary number leads to a decrease in recirculating volume and an increase in recirculation time, which means that transport processes will be hindered. Radial mixing is then reduced while axial dispersion is enhanced with increasing the capillary number. It has been shown that when the aspect ratio increases for a given capillary number, the recirculating volume decreases slightly and the recirculation time increases while the exchange surface between film and recirculating regions is enhanced at low slip velocities, which is the case for capillary numbers lower than approximately Ca ∼ 0.1. Indeed, high aspect ratio channels may appear attractive for heat transfer in plate microreactors, however the decrease in the recirculation rate and effective recirculation volume is disadvantageous for both heat and mass transfer operations. It is therefore expected that some intermediate aspect ratio geometry would be most effective. This would require a systematic study that investigates the effects of process parameters (aspect ratio, void fraction, capillary number, Reynolds number) on the efficiency of the gas-liquid system that could be expressed by a cost-function integrating global parameters (pressure drop, circulation volumes and times, mass transfer, heat transfer).

The disturbed flow close to the bubble caps has shown to generally increase recirculation time. It has also been shown that for certain cases in rectangular channels, the bubble wake causes a local increase in the velocity at the center of the channel and this phenomenon is enhanced in high aspect ratio channels and high capillary number flows. Although this increase in the axial velocity can lead to a slight decrease in the recirculation time, these wake effects do not balance geometrical effects and the recirculation times in wide shallow channels remain greater than the ones obtained in channels of lower aspect ratio.

  Fig. 2. The local velocity field is shown in terms of local liquid Peclet number Pe l = Re × Sc = |U l |D h /D L where |U l | is the local velocity magnitude in the frame of reference moving with the bubble, D L is the mass diffusivity in the liquid and Sc is the Schmidt number and is set to 1.

Fig. 2 -

 2 Fig. 2 -Local liquid Peclet numbers and streamlines through the channel width and depth for channel aspect ratio of 2.5. (a) Ca = 0.06, Re = 16.8; (b) Ca = 0.1, Re = 28; (c) Ca = 0.2, Re = 56; (d) Ca = 0.4, Re = 112. The bold line represents a Peclet number of 2 for (a)-(c) and 10 for (d) where Ca = 0.4 and Re = 112.

Fig. 3 -

 3 Fig. 3 -Local liquid Peclet numbers and streamlines through the channel width and depth for Ca = 0.1 and Re = 28. (a) ˛= 1; (b) ˛= 2.5; (c) ˛= 4. The bold line represents a Peclet number of 2.

  Fig. 4 -Dimensionless liquid velocity along the centerline of the channel between two bubbles for an aspect ratio of 2.5 and various capillary numbers. Legend: (-) U centerline /U B ; (--) U max /U B evaluated from U TP and Eqs. (1) and (2); (black) Ca = 0.06; (gray) Ca = 0.1; (light gray) Ca = 0.2.

  Fig. 5 -Dimensionless liquid velocity along the centerline of the channel between two bubbles for Ca = 0.1 and for aspect ratios 1, 2.5 and 4. Legend: (-) U centerline /U B ; (--) U max /U B evaluated from U TP and Eqs. (1) and (2); (black) ˛= 1; (gray) ˛= 2.5; (light gray) ˛= 4.

Fig. 6 -

 6 Fig. 6 -Cross-sectional shape of the recirculating zone for aspect ratios ˛= [1 ; 2.5 ; 4] for different dimensionless velocities W. Legend: (-) center of the recirculating loop; (--) streamline separating the circulation loop and the film. From the walls to the center (green; blue; red; black): W = [0.1 ; 0.2 ; 0.3 ; 0.42]. (For interpretation of the references to color in this figure legend, the reader is referred to the web version of the article.) depend on the total recirculation time, as well as the size of the recirculating volume.

Fig. 7 -

 7 Fig. 7 -Perimeter of the recirculation zone through the channel cross-section (as depicted by dashed lines in Fig. 6) normalized by the perimeter of a circular channel with an equivalent hydraulic diameter D h in a cross-section of the slug as a function of the dimensionless velocity W. Legend: (gray) circular channels; (red) ˛= 1; (black) ˛= 2.5; (blue) ˛= 4. (For interpretation of the references to color in this figure legend, the reader is referred to the web version of the article.)

Fig. 9 -

 9 Fig. 9 -Dimensionless recirculation volumes as a function of dimensionless slug length. Legend: (filled symbols) V 0 /V slug ; (open symbols) V rc /V slug .

Fig. 11 -

 11 Fig. 11 -(a) Dimensionless recirculation times 0 and 1 versus the dimensionless bubble velocity. Legend: (-) 0 , infinite slug length assumption; (-.-) 1 , infinite slug length assumption; (light gray) 2D channels; (filled symbols) 0 , numerical ( , red)˛= 1; (⊲ , black)˛= 2.5; (♦ , blue)˛= 4. (b) Dimensionless total recirculation time rc versus the dimensionless bubble velocity. Legend: (--) rc , infinite slug length assumption; (light gray) 2D channels; (open symbols) rc , numerical; ( , red)˛= 1; (⊲ , black)˛= 2.5; (♦ , blue)˛= 4. (For interpretation of the references to color in this figure legend, the reader is referred to the web version of the article.)

Fig. 12 -

 12 Fig. 12 -(a) Dimensionless recirculation times as a function of dimensionless slug length. Legend: (--) asymptotic values 0 and rc ; (filled symbols) 0 ; (open symbols) rc . (b) Dimensional recirculation times t 0 and t rc (with D h = 571 m and U B = 0.5 m/s) as a function of the time taken for a bubble to travel the slug length L S /U B . Legend: (--) asymptotic values t 0 and t rc ; (filled symbols) t 0 ; (open symbols) t rc .

Fig. 13 -

 13 Fig. 13 -Dimensionless recirculation times versus the capillary number. Legend: (-) 0 , infinite slug length assumption; (--) rc , infinite slug length assumption; (gray) circular channels; (light gray) 2D channels; (filled symbols) 0 ; (open symbols) rc ; (⊳ , black)˛= 2.5, experimental (Zaloha et al., 2012); ( , red)˛= 1; (⊲ , black)˛= 2.5; (♦ , blue)˛= 4, numerical. (For interpretation of the references to color in this figure legend, the reader is referred to the web version of the article.)

Acknowledgements

This work was financed by the French "Agence Nationale de la Recherche" in the framework of the project MIGALI no. ANR-09-BLAN-0381-01. We also acknowledge the support from this project from the CNRS research federation FERMaT, such as the CALMIP, IDRIS and CINES projects for providing computational resources.


Fan Zhang

Audrey Fluck

Stephane Baudron

Wais Mir

Hosseini

Strapping a benzaldehyde-appended 2

come

Strapping a benzaldehyde-appended 2,2'-bis-dipyrrin Zn(II) double-stranded helicate using imine bond formation Fan Zhang, Audrey Fluck, Stéphane A. Baudron* and Mir Wais Hosseini* Laboratoire de Tectonique Moléculaire, UMR CNRS 7140, icFRC, Université de Strasbourg, F-67000, Strasbourg, France. Fax: 33 368851325; Tel: 33 368851323; E-mail: sbaudron@unistra.fr, hosseini@unistra.fr Abstract: A novel linear binuclear Zn(II)/2,2'-bis-dipyrrin double-stranded helicate bearing peripheral p-benzaldehyde moieties has been strapped by formation of imine bonds upon reaction with m-xylylenediamine and characterized by 1 H-and 13 C-NMR spectroscopy as well as by UV-Visible spectroscopy and by high-resolution mass spectrometry.

Since their discovery, metallahelicates have been the object of intense interest owing to their intrinsic beauty and their potential optical properties, for example. 1 This class of compounds has also been employed as templates for the construction of more complex architectures such as molecular knots. 2 This approach relies on the use of functionalized ligands for the construction of helicates bearing peripheral reactive groups like alkene, aldehyde, amine or alcohol. The latter moieties may be reacted with a strapping fragment. Whereas bis-pyridine based ligands have been extensively investigated for the elaboration of helicates, 1 the openchain tetrapyrrolic derivatives of the bis-dipyrrin type have been also recently studied. [3][4][5] In particular, the 2,2'-bis-dipyrrin ligands have been shown to form linear as well as circular helicates, upon assembly with Zn(II) cations (Scheme 1). 3 We have recently demonstrated that formation of such helical arrangements can be envisioned as a template offering the preorganization of both the pyrrolic strands and the peripheral aromatic units allowing the targeted complexation of Ag(I) cations by Ag-π interactions. 4 In light of this observation, it appeared interesting to take advantage of this effect beyond the generation of such interaction poles but also for the pre-organization of reactive sites at the periphery of the helicate. To the best of our knowledge, only one report has addressed this issue. 5 The group of Maeda has described a series of 2,2'-bis-dipyrrin Zn(II) helicates bearing peripheral olefin groups leading to strapped systems upon olefin metathesis (Scheme 1). 5 Scheme 1. Schematic representation of a 2,2'-bis-dipyrrin ligand (left), a linear binuclear helicate (middle) and of a strapped helicate (right).

Over the past few years, the imine bond formation and reduction have been successfully employed for the elaboration of complex (supra)molecular systems. 6 In order to take advantage of this reactivity and versatility, we have designed and prepared helicate 1 bearing peripheral p-benzaldehyde moieties (Scheme 2) as a template platform for the preparation of strapped helicates upon reaction with an appropriately chosen diamine. In this contribution, we report on the synthesis and characterization of 1 as well as its reactivity with m-xylylenediamine affording a strapped tetra-imine analogue.

The novel helicate 1 has been prepared following the Ni-templated strategy developed by Scott and coworkers (Scheme 2). 7 The homoleptic Ni(II) complex 2, prepared as described for the Zn(II) analogue, 8 has been reacted with DDQ in refluxing toluene to afford compound 3 in 80 % yield. The latter was demetallated by treatment with HCl and then NaHCO3 leading to the free ligand 4 in 99 % yield. Upon reaction of this tetrapyrrolic derivative with Zn(OAc)2, helicate 1 was obtained in 94 % yield. These compounds have been characterized in solution by 1 H-, 13 C-NMR and UV-visible spectroscopies, by mass spectrometry and, in the solid state, by single-crystal X-Ray diffraction (see ESI for details).

As observed for other Ni(dipyrrinato)2/Ni(2,2'-bis-dipyrrinato) systems, 3h,4,5,7 a flattening of the coordination sphere around the tetra-coordinated Ni(II) center is observed in the crystal structures of the compounds (Fig. 1). Whereas the two bipyrrolic chelates form an angle of 60.88° in 2, this angle decreases to 19.15° in 3. Also, the Ni-Ni bond distances are in the rather narrow 1.8835(12)-1.8930(12) Å range in the structure of 2, while they vary from 1.8584(17) to 1.9027(18) Å for 3, owing to the formation of the 2,2'-bis-dipyrrin backbone. For the free ligand 4, the two bipyrrolic chelates are oriented in a trans fashion and coplanar as described for the benzonitrile analogue. 4a Furthermore, an intramolecular hydrogen bond is observed 1). In the structure of the binuclear helicate 1, each Zn(II) cation is in a distorted tetrahedral environment coordinated to two dipyrrin chelates belonging to the two different strands forming an angle of 53.6°. The metal centers are separated by 3.280(8) Å, as observed for other helicates (Fig. 1). 3- Scheme 2. Synthetic pathway for the preparation of helicate 1 bearing p-benzaldehyde units and of strapped helicates 5 and its reduced analogue 6.

within each dipyrrin moiety (dN•••N = 2.753(3) Å; αN-H•••N = 122.45°) (Fig.
Aiming at strapping the helicate on each side by reaction of the peripheral aldehyde groups with an appropriately chosen reactant, we considered m-xylylenediamine, since the distance between the two nitrogen atoms within this compound should be compatible with the one between the carbonyl oxygen atoms of 5.797 Å observed in the structure of 1 (Fig. 1). The reaction of two equivalents of m-xylylenediamine with one equivalent of 1 in the presence of a catalytic amount of TFA (0.01 eq.) was thus investigated. A first attempt consisting in the direct mixing of the reactants in a 1.10 -4 mol.L -1 CHCl3 solution of 1 led, after three days, to the complete disappearance of the CHO signal in the 1 H-NMR spectra suggesting the formation of the desired product 5. However, the purity of the compound was not satisfactory and a low yield was observed given the concomitant formation of insoluble material in large amount (about 50 %). The latter could not be further characterized but may correspond to non-strapped oligo-/poly-meric assemblies. Lowering the initial concentration to 1.10 -5 mol.L -1 led also to product 5, although insoluble polymeric material was again obtained. Pure product 5 could be isolated in almost quantitative yield by the controlled slow addition of the diamine over 24 hours. The 1 H-NMR spectrum of the product in CDCl3 confirms the disappearance of the carbonyl proton signal along with the appearance of a signal at 8.50 ppm corresponding to the imine (Fig. 2). Furthermore, it should be noted that the two protons of the methylene bridges of the m-xylylenediimine straps appear as two doublets at 4.90 and 5.06 ppm, indicating that they are in a different environment (Fig. 2). For both 1 and 5, the aromatic protons of the meso groups are present as four distinct doublets, as described for RCM strapped analogues 5 and for a recently reported amide substituted helicate. 9 Interestingly, the chemical shift of the α-pyrrolic proton is substantially modified, as it is observed at 6.98 ppm for 1 and at 7.37 ppm for 5 (Fig. 2 and ESI for assignment of the protons and carbon atoms by 2D NMR). Such a phenomenon has been rationalized by Maeda and coworkers 5 as the complex adopting different helical forms with varying Zn-Zn distances and hence different organizations of the pyrroles within a strand and around the metal centers. With a short or constraining strap, the complex is in a so-called long form with the metal cations further apart (dZn-Zn of ca 5.1 Å), in comparison with unstrapped derivatives or analogues strapped by a longer chain thus adopting a short form (dZn-Zn of ca 3.50 Å). 5 This suggests here that the helicate adopts a short mode in 1, as confirmed by X-ray diffraction (Fig. 1), and a long mode upon strapping by imine bond formation, in 5.

Unfortunately, no crystals of sufficient quality for structure determination could be obtained for the strapped helicate. The nature of 5 was also confirmed by high-resolution mass spectrometry showing the presence of the desired compound as the main peak with the appropriate isotope pattern (Fig. ESI 14). The UV-Visible spectra of the starting ligand 4 and helicate 1 as well as of strapped compound 5 were recorded in CH2Cl2, as presented in Fig. 3. As for other reported helicates, 3-5 a red-shift of the absorption band and the appearance of a sharp intense band around 430 nm and a broader one around 650 nm are observed for 1 as a result of the coordination of 4 to Zn(II) cations.

Upon strapping of the helicate, a blue shift is observed with the broad low energy band being slightly better resolved into two bands around 550 and 650 nm for 5, and a relative lowering of the absorption coefficient of the band around 430 nm. A similar behaviour has been described with species strapped upon olefin metathesis. 5 Once again, it suggests that strapping affects the relative orientation of the dipyrrin chromophores within a tetrapyrrolic strand and around metal cations with the latter being further apart from one another. This is consistent with the abovementioned observation by 1 H-NMR of these compounds. We have attempted to convert the tetra-imine based helicate 5 into the reduced tetra-amine analogue 6, by reaction with NaBH3CN in CHCl3. 10 In spite of numerous attempts, we have not been unable to obtain 6 in pure form. However, 1 H-NMR spectroscopy in CD2Cl2 and highresolution mass spectrometry confirmed the formation of compound 6. Indeed, while the imine proton signal disappeared, the methylene bridge protons, formed upon reduction, appeared at 3.90 ppm (Fig. ESI 15). The signal associated with the protons of the methylene bridge of the m-xylylenediamine strap appeared to be shifted from around 5.0 to 4.0 ppm upon reduction of the imine bonds. In this case, the signal corresponding to the α-pyrrolic proton appears at 7.06 ppm suggesting that the helicate is in an intermediate form between long and short forms. This is consistent with the reduction to amines releasing some of strain brought by the two di-imine based junctions.

In conclusion, a novel 2,2'-bis-dipyrrin ligand bearing peripheral p-benzaldehyde units has been designed and prepared. Its assembly with Zn(II) cations leads to the formation of the linear binuclear double-stranded helicate 1, featuring the reactive aldehyde groups on each side of the complex at a relative distance of ca 5.8 Å according to the structure investigation by singlecrystal X-Ray diffraction. This pre-organization of the functions at the periphery resulting from the template effect of the Zn(II) bnding has been taken advantage of for the strapping of the helicate upon reaction with m-xylylenediamine. Slow addition of the latter to a dilute solution of the helicate in CHCl3 in the presence of catalytic amount of TFA afforded quantitatively the strapped tetra-imine helicate 5 as demonstrated by 1 H-and 13 C-NMR spectroscopy as well as by high-resolution mass spectrometry. Furthermore, these characterizations combined with UV-Visible spectroscopy suggest that the strapped helicate 5 adopts a long form with Zn cations further apart than in the starting complex 1, in agreement with recently reported theoretical calculations. 5 Reduction by NaBH3CN has been explored and indicates the formation of the tetra-amine analogue 6. Future work will focus on the use of other functionalized diamines towards the construction of diverse assemblies as well as on the use of other reactive groups at different positions on the bis-dipyrrin backbone. Furthermore, the photophysical properties of the strapped helicates, in particular their emission, will be thoroughly investigated. Indeed, as demonstrated herein in the absorption studies, strapping impacts the organization of the tetrapyrrolic system but it may also affect the rotation of the peripheral meso groups known as a favourable non-radiative deactivation pathway in such systems. 3g,11 Hence, it may result in modified emission wavelength and quantum yield. These investigations are currently under way in our laboratory.

Fig. 1 .

 1 Fig. 1. Crystal structures of Ni complexes 2 and 3, of 2,2'-bis-dipyrrin ligand 4 and of helicate 1. Solvent molecules have been omitted for clarity.

Fig. 2 .

 2 Fig. 2. 1 H-NMR spectra (CDCl3, 500 MHz, 293 K) of helicates 1 (top) and 5 (bottom).

Fig. 3 .

 3 Fig. 3. Absorption spectra of ligand 4 and helicates 1 and 5 in CH2Cl2 solution at room temperature.

In the crystal (Tetragonal, P-4n2), both P and M helicates are present.

Acknowledgments

We thank the Université de Strasbourg, the C.N.R.S, the International centre for Frontier Research in Chemistry (icFRC), the Labex CSC (Ph. D. fellowship to F. Z.; ANR-10-LABX-0026 CSC) within the Investissement d'Avenir program ANR-10-IDEX-0002-02, the Institut Universitaire de France, the Ministère de l'Enseignement Supérieur et de la Recherche for financial support. Dr. Antoine Béziau, Dr. Hervé Ruffin and Fiona Schmitt are acknowledged for preliminary investigation of the synthesis of 2 and 3; Nathalie Kyritsakas is thanked for her assistance in the structure determination of 1.

Conflicts of interest

There are no conflicts to declare.


Jeremy Riviere
email: jeremy.riviere@imag.fr

Carole Adam
email: carole.adam@imag.fr

Sylvie Pesty
email: sylvie.pesty@imag.fr

Catherine Pelachaud
email: catherine.pelachaud@telecom-paristech.fr

Nadine Guiraud
email: nadine.guiraud@irit.fr

Dominique Longin
email: dominique.longin@irit.fr

Emiliano Lorini
email: emiliano.lorini@irit.fr

Expressive Multimodal Conversational Acts for SAIBA Agents

Keywords: ECA, Emotions, Interaction, Dialogue

We discuss here the need to define what we call an agent conversational language, a language for Embodied Conversational Agents (ECA) to have conversations with a human. We propose a set of Expressive Multimodal Conversation Acts (EMCA), which is based on the Expressive Speech Acts that we introduced in a previous work, enriched with the multimodal expression of the emotions linked to these acts. We have then implemented these EMCA for SAIBA-compliant ECA, and specifically for Greta. We were able to use Greta in experiments aimed at assessing the benefits of our language in terms of perceived sincerity and believability of an ECA using it to interact with a human user.

Introduction

Nowadays, human-computer interactions and agent-agent interactions are omnipresent, making dialogue a major research problem. Therefore, some recent works in the field of human-computer interaction have turned to dialogue annotation [START_REF] Core | Coding dialogs with the DAMSL annotation scheme[END_REF][START_REF] Bunt | The dit++ taxonomy for functional dialogue markup[END_REF] for understanding dialogue, recognizing its structure and generating dialogues. In parallel, some other works have turned to natural language communication. With the advance of virtual agents, multimodal aspects have become essential in human-agent dialogue. In particular these agents are now able to express (their) emotions, but mainly in a non-verbal manner. We believe that it is important to closely link verbal and non-verbal aspects in order to improve the conversational capabilities of these agents. This is particularly true for what we call complex emotions (involving mental representations as cause, self, action ...) that are (in humans) mainly conveyed through language, and are therefore often neglected in current research. In this paper we thus focus on complex emotions such as guilt or reproach and propose a set of Expressive Multimodal Conversational Act (EMCA) that allows an ECA to express them in a multimodal manner: verbal and non-verbal.

We ground on the set of Expressive Speech Acts that we proposed in [START_REF] Guiraud | The face of emotions: a logical formalization of expressive speech acts[END_REF] and enrich them with the multimodal expression of their intrinsic emotion, thus allowing the automatic generation of the expression of the emotions linked to the act uttered by the agent. We believe that an agent able to not only verbally express sentence but also multimodally express the underlying emotions will be more sincere and believable to the user. To check our claim, we implemented our MCL in the standard SAIBA architecture [START_REF] Vilhjálmsson | The behavior markup language: Recent developments and challenges[END_REF], in particular in the Greta agent [START_REF] Poggi | Greta: A believable embodied conversational agent[END_REF], and then conduced experiments with this agent.

Speech Acts to Express complex emotions

Starting from the formalization of counterfactual emotions proposed by [START_REF] Lorini | A logic for reasoning about counterfactual emotions[END_REF], in [START_REF] Guiraud | The face of emotions: a logical formalization of expressive speech acts[END_REF], we defined complex emotions as based on one hand upon the agent counterfactual reasoning and on the other hand upon reasoning about responsibility, skills and social norms. We formalized eight complex emotions (rejoicing, gratitude, regret, disappointment, moral satisfaction, admiration, guilt and reproach) in terms of three types of logical operators representing the agent's mental state: beliefs (Bel i ϕ), goals (Goal i ϕ) and ideals (Ideal i ϕ); and one operator for the attribution of responsibility (Resp i ϕ). We then used the operator Exp i,j,H ϕ, representing what an agent i expresses to an agent j in front of group H, to define eight Expressive Speech Acts, each of which expresses one complex emotion. Each emotion is said to be "intrinsically" attached to each act. The figure 1 shows the set of Expressive Speech Acts according to the complex emotions they are expressing.

Goal iϕ ∧ Bel iResp i ϕ déf = Exp i,j,H (Rejoicing i ϕ) déf = Rejoice i,j ϕ Goal iϕ ∧ Bel iResp j ϕ déf = Exp i,j,H (Gratitude i,j ϕ) déf = Thanki,j ϕ Goal i¬ϕ ∧ BeliResp i ϕ déf = Exp i,j,H (Regret i ϕ) déf = Regret i,j ϕ Goal i¬ϕ ∧ BeliResp j ϕ déf = Exp i,j,H (Disappointment i,j ϕ) déf = Complain i,j ϕ Idealiϕ ∧ BeliResp i ϕ déf = Exp i,j,H (MoralSatisfaction i ϕ) déf = IsMorallySatisfied i,j ϕ Idealiϕ ∧ BeliResp j ϕ déf = Exp i,j,H (Admirationi,j ϕ) déf = Compliment i,j ϕ Ideali¬ϕ ∧ Bel iResp i ϕ déf = Exp i,j,H (Guiltiϕ) déf = FeelGuilty i,j ϕ Ideali¬ϕ ∧ Bel iResp j ϕ déf = Exp i,j,H (Reproach i,j ϕ) déf = Reproach i,j ϕ
Thus, this Expressive Language links sharply complex emotions and speech acts, meaning you can't compliment somebody without expressing the complex emotion of admiration. Indeed, according to the speech act theory, an ECA has to express the intrinsic emotion of the expressive acts to "successfully" perform this act (conditions of success). We affirm that the multimodal expression of the EMCA's intrinsic emotions expressed by an ECA increases its perceived sincerity and believability, under the human-ECA dialogue context and the cultural context.

Proposition: EMCA for SAIBA Agents

Our proposal consists in adding to the Expressive Speech Acts we defined the multimodal expression of each act's intrinsic emotion. Expressive Speech Acts thus become Expressive Multimodal Conversational Acts (EMCA), containing the emotions intrinsically attached to each act, as well as their multimodal expression. We have implemented these EMCA into the SAIBA architecture to allow ECA SAIBA-compliant to express their complex emotions. We are only interested here in the expressive category of acts, whose aim is specifically to express emotions.

The SAIBA framework consists of three step process of intention planning, behavior planning and the virtual agent behavior realization. The FML (Function Markup Language [START_REF] Mancini | The fml-apml language[END_REF]) and BML (Behavior Markup Language [START_REF] Vilhjálmsson | The behavior markup language: Recent developments and challenges[END_REF]) languages link the three modules. An xml file (the lexicon file) is used by the behavior planning module to convey patterns in multimodal signals before encoding into BML. Such patterns can be emotions, some performatives (like inform, offer or suggest) and other behavior. The lexicon file contains the list of these patterns as well as their multimodal expressions, which are mainly defined from psychological theories and video corpus analysis. Thanks to this lexicon file the behavior planning thus links the present patterns in the agent's communicative intention and their multimodal signals before encoding into BML and sending to the behavior realization module.

Our work consisted in offering a library of EMCA to SAIBA-compliant ECA. These ECA can thus choose the speech act matching with its communicative intention, thanks to the appropriate dialogue processing or reasoning module (which should have a set of rules based on our formalization). We have modified the lexicon file, where we have implemented our library of EMCA. Every act in the library is associated with its intrinsic emotion (see Section 2): concretely, this emotion is manifested by various multimodal signals that are used in the behavior planning module, and it will thus be expressed by the ECA alongside the conversational act.

For example, let us consider the next turn in a dialogue with the user: "I am disappointed in your behavior, It's not very nice of you.". Here, the ECA's intention is to express that the ECA had a certain goal and that the human is responsible for not reaching this goal. In this case, the EMCA complain, whose intrinsic emotion is the complex emotion of disappointment, fits best the ECA's communicative intention. The first module translates this communicative intention in FML-APML 1 , and completes it with a BML specification of the agent's utterance, as shown in the example below: <fml-apml> <bml> <speech id="s1" language="english" voice="openmary" text="I am disappointed in your behavior, It's not very nice of you."> </bml> <fml> <performative id="p1" type="complain" > </fml> </fml-apml>

In the FML-APML specification, the performative tag represents the EMCA that is linked with the utterance. The multimodal expression of the emotion associated with each EMCA, which was identified in FML by the performative tag, is now detailed in the lexicon file. For example, here is the EMCA complain in the EMCA specification general format: <behaviorset name ="complain"> <signal id="1" name="faceexp=disappointment" modality="face"/> Here, the complex emotion of disappointment associated to the EMCA complain is multimodally expressed by its facial expression (signal 1). The behavior realization Module finally receives all the signals translated in BML, and it adapts them to animation rules and constraints: it checks coherence and absence of conflicts, before animating the verbal and non-verbal behavior of the ECA. So, the ECA has to express the emotion of disappointment to "successfully" complain (i.e. to seem sincere and believable to the user in the dialogue and cultural context).

Unlike current works [START_REF] Lee | Nonverbal behavior generator for embodied conversational agents[END_REF][START_REF] Gebhard | Adding the emotional dimension to scripting character dialogues[END_REF], we established a link between acts and emotions from both the speech act theory [START_REF] Searle | Speech acts: an essay in the philosophy of language[END_REF] and our complex emotion definition. This approach enable a consistency of these acts and should improve the virtual agent's sincerity and credibility. To confirm this hypothesis, we assessed the Greta agent's sincerity and credibility when using our EMCA (see next section).

EMCA Evaluation: Application to Greta

The hypothesis that we want to test is the following: adding non-verbal signals of the expression of the complex emotions coherent with the verbal signals expressed by an ECA increases its perceived sincerity and believability. This hypothesis is to be considered through the human-ECA dialogue context, and the cultural context. We ran a first experiment with a subset of our EMCA and Greta, a SAIBA-compliant ECA ( [START_REF] Poggi | Greta: A believable embodied conversational agent[END_REF]). We have implemented our library of EMCA in the lexicon file of Greta (as described in 3), thus endowing this character with the ability to use our EMCA.

Test Scenario and Protocol

The protocol that we set up allowed us to test our hypothesis for two EMCA that express complex emotions: to rejoice, that expresses rejoicing, and to apologise, that expresses regret 2 . In the setting of a dialogue between a user and a ECA (acted out by Greta) whose aim is to recommend movies, these two EMCA are expressed by the ECA in three conditions: a congruent condition, a neutral condition and a non-congruent condition. In the congruent condition, the ECA expresses the EMCA as it is defined in our language, i.e. with the intrinsic emotion matching the corresponding act (the regret emotion with an Apologise utterance, and the rejoicing emotion with a Rejoice utterance). In the neutral condition, the ECA only expresses the utterance without expressing its intrinsic emotion of the act. Finally in the non-congruent condition, the ECA expresses the utterance with the emotion that is the opposite (see table 1) of the intrinsic emotion of the act (an Apologise utterance is expressed with the Rejoicing emotion, while a Rejoice utterance is expressed with a Regret emotion).

In the first stage of the protocol, videos of two scenarios are proposed, each containing the expression by the ECA of one of the two tested EMCA. For each scenario, three videos are submitted to the subject (one per condition, i.e altogether six videos) to be evaluated. Both scenarios imply a dialogue with Greta after the user supposedly comes back from the cinema where they saw a movie that had previously been recommended by Greta. In the first scenario, they loved the movie and has just thanked Greta for her good advice. The user is then asked to watch the three videos (in a predefined order that differs for each user), and to evaluate how Greta rejoices. Similarly in the second scenario, the user hated the movie and has just reproached Greta for her bad advice (cf. figure 1). The user is then asked to watch the three videos of the EMCA to apologize (also in a predefined order) and to evaluate how Greta apologises. The videos are assessed on two criteria: sincerity ("Does Greta seem to express what she thinks?") and believability ("Does the scene played by Greta look plausible?") of the ECA. These criteria can each take four possible values on the same qualitative scale (Not at all -Rather not -Rather -Totally).

The second part of the protocol consists in a questionnaire whose aim is to collect the user's impressions, mainly regarding the impact of the agent, and their subjective feelings. During this part, users can express themselves with their own words.

Results

23 users selected within students in science between 18 and 26 participated in the first experiment. To perform the statistical tests, we associated a score with each qualitative value of the sincerity and believability criteria: 1 for "Not at all", 2 for "Rather not", 3 for "Rather", and 4 for "Totally". Since we evaluate 2 characteristics of 2 EMCA, we have 4 dependent variables: Sinc Apologize, Believ Apologize, Sinc Rejoice, Believ Rejoice. The only independent variable is the condition (congruent, neutral, non-congruent).

The first step was to check that the dependent variables follow a normal distribution. Secondly, the analysis of the results told us that the EMCA expressed in the congruent condition (C) have the best scores in both sincerity and believability. When Greta apologises under this condition, 65% of participants found Greta "rather" sincere (39%) or "totally" sincere (26%) and 70% of participants found Greta "rather" believable (52%) or "totally" believable (18%). Similarly when Greta rejoices under this condition, 74% of participants found Greta "rather" sincere (52%) or "totally" sincere (22%) and 78% of participants found Greta "rather" believable (52%) or "totally" believable (26%). This tends to confirm our hypothesis; we can thus assume that expressing the coherent emotion (C) with the act makes the ECA more sincere and believable to the user than expressing nothing (N) and expressing the opposite emotion (NC).

We then ran a one-way Analysis of Variance (ANOVA) test on each of our dependent variables with respect to the independent variable. The hypothesis H0 we checked through this ANOVA test is the following: "There is no difference between the sincerity means (resp. the believability means) under the condition of congruence for the EMCA Apologize and Rejoice". H1 thus says that "there is a difference between the sincerity means (resp. the believability means) under the condition of congruence for the EMCA Apologize and Rejoice". This ANOVA test allowed us to reject H0 significantly in both sincerity and believability: for example, it showed a significant effect of the condition of congruence on Greta's sincerity (F(2.66)=22.80 p<0.05) for the EMCA Apologize. A Tukey range test showed us that there is a significant difference between the sincerity (resp. believability) mean in the congruent condition (C) and the sincerity (resp. believability) means in the other conditions for both EMCA Apologize and Rejoice, while there is no significant difference between the non-congruent condition (NC) and the neutral condition (N) (see for ex. table 2). Table 2. Tukey range test of the Sinc Apologize variable: the difference between the congruent condition (C) and the neutral(N) and the non-congruent (NC) conditions is significant Groups Difference Statistic Probability C -NC 1.478 q = 9.290 0.0000 C -N 1.478 q = 9.290 0.0000 NC -N -0.435 q = 2.732 0.1376 Finally, when asked "What do you think about Greta?" during the qualitative interview, a lot of users assessed the Greta's personality. Despite our efforts to avoid this evaluation, Greta's personality Prudence [START_REF] Mcrorie | A model of personality and emotional traits[END_REF] has been described as "austere" and "severe". We do not know how it has influenced the results, but it appears that several users has mentioned the trust notion; Greta's personality and aspect seems to be as important as its expression of coherent emotion to the user, as shown in [START_REF] Gong | How social is social responses to computers? the function of the degree of anthropomorphism in computer representations[END_REF].

Conclusions and Perspectives

In this paper we have presented a library of Expressive Multimodal Conversational Acts, that make a link between the expression of emotions through language and through other modalities. We formally defined and implemented this set of eight Expressive MCA in the SAIBA architecture. Two of them were evaluated through an implementation in Greta, during a first experiment that showed that they contribute to a better perceived sincerity and believability of ECA. A second phase of evaluation will concern all our library in longer dialogue scenarios involving several turns of speech of an ECA and a human user. In order for this library to be useful to other SAIBA-compliant ECA, we intend to complete it with MCA from the other categories: assertive (inform, tell ...), directive (ask, offer ...), commissive (promise, assure ...).

Fig. 1 .

 1 Fig. 1. Apologise EMCA in the congruent condition: Greta apologises in reply to user's reproach

Table 1 .

 1 Expressive Speech Acts expressing complex emotions

	Exp i,j,H (...)	Expressive Speech Acts

Since FML is still under development, we will work here on the FML-APML specification proposed in[START_REF] Mancini | The fml-apml language[END_REF].

To apologise is to regret something which is against the human's goal.

Acknowledgments. This work is supported by the French National Research Agency (ANR), project CECIL (www.irit.fr/CECIL), grant number ANR-08-CORD-005. We thank the statistician engineer and psychologist Nadine Mandran of the LIG Marvelig unit which is specialized in experimentation protocol.