Matthieu Lemerre

Emmanuel Ohayon

Damien Chabrol

Mathieu Jan

Marie-Benedicte Jacques

Marie -Bénédicte

LIST, Jacques Cea

Method and Tools for Mixed-Criticality Real-Time Applications within PharOS

come

I. INTRODUCTION

Over the last decades, the number of electronic equipments embedded in complex systems such as automotive vehicles or airplanes has increased in drastic proportions. At first, system architects used to embed one computing unit per functional system, which resulted in large distributed and heterogeneous real-time systems. The current number of equipments makes this approach no more sustainable, mainly due to its high costs. As a result the trend is to integrate several (possibly unrelated) applications on fewer computing units, resulting in what is called a mixed-criticality system.

This integration raises however a number of issues, such as:

• efficiently and securely sharing resources (and especially CPU time) between these applications; • providing the same degree of isolation than achieved when using separate computing units; • combining different real-time paradigms and APIs;

• handling multicore computing units, that are becoming the usual solution of the chip designers to increase performances.

This paper presents the formal and/or technical answers we provided to each of these questions. Each proposed solution is illustrated through PharOS, the system we have built following these methods and that addresses all the problems hereabove exposed.

PharOS features a rich model of tasks that accurately describes the task constraints, and admits an optimal dynamic scheduling algorithm. Moreover, it allows exact feasibility analysis, which avoids wasting CPU time by oversizing the time resources for the tasks. Its tool suite allows for automatic fine-grained memory partitioning, if the underlying platform comes with a Memory Protection Unit (MPU). At last, some additional extensions allow compliance with other real-time paradigms, or paravirtualization of other existing RTOSes and APIs, while taking advantage of the additional power provided by multicore computing units.

The article is divided as follows: Section II presents the PharOS formal task model. Section III proposes implementation mechanisms for temporal and spatial isolation, fault containment and recovery. Section IV focuses on cohabitation with the event-triggered paradigm and with another Operating System. In particular, we present here a paravirtualization usecase of Trampoline [START_REF] Bechennec | Trampoline An Open Source Implementation of the OSEK/VDX RTOS Specification[END_REF], a RTOS with an OSEK API, over PharOS. Section V presents related work, and Section VI concludes and suggest directions for future research.

II. TIME-CONSTRAINED AUTOMATA (TCA): THE TASK MODEL FOR PharOS

At the heart of PharOS is a formal model of task, the time-constrained automata model [START_REF] Lemerre | An introduction to time-constrained automata[END_REF], which is based on time. This section gives a short overview of this model and its implementation, which is necessary to understand the design of the temporal isolation in PharOS.

A. Short introduction to the TCA model

A PharOS application is composed of a set of tasks. The TCA model is a formal specification of all the possible temporal behaviors for each task. This specification is expressed using a directed graph, where arcs represent the successive jobs to be executed (one at a time), and the nodes bear the temporal constraints on the task. There are only four kinds of nodes:

• an after(d) (written d ) node constrains the execution of the following jobs to start after date d; • a before(d) (written d ) node constrains the execution of the preceding jobs to end before date d; • an advance(d) (written d ) node is a combination of an d and a d nodes;

• a no-constraint (written ) node imposes no temporal constraints.

All constraint dates are expressed relatively as an increment to the previous d or d node (i.e. d nodes do not change the "reference date"). Nodes with multiple outgoing arcs represent conditional execution (e.g. if statements).

Execution of an application can be seen as choosing one (possibly infinite) walk in the graph of each task, and choosing for each encountered job the intervals of time where it is actually executed. The execution is correct if, and only if:

• successive jobs of a task execute in order;

• the execution intervals for each job fulfill the timing constraints of all nodes.

B. An example using ΨC

The graph for each task need not be written manually: it can be automatically extracted from the control flow graph of an application. The language used to implement PharOS applications is ΨC, a variant of the C language with extensions for timing constraints. Figure 1 gives an example of a code excerpt for a task with the corresponding time-constrained automaton. Timing constraints statements are transformed into timing constraints nodes. The beginning of a possible correct execution for this example is represented on Fig. 2. The code executes the while loop twice, and chooses the "else" part the first time, the "then" part the second time. The upper part of the figure represents a timescale, with the timing constraints and their date (expressed absolutely), and the blocks of code that are executed between these blocks. The lower part represents execution as a Gantt diagram. Execution is correct because blocks execute in order, and all timing constraints are fulfilled.

C. Implementation in PharOS

PharOS implements the previously described TCA formal model. 1) ΨC compilation and scheduling: The automaton for each task is extracted from the ΨC code, and used for different timing analysis phases 1 . Eventually, ΨC code is compiled as regular C code, where timing constraints statement have automatically been transformed into system calls to the scheduler. Each system call indicates to the kernel the node that has been reached. The node contains the new start time, deadline and timing budget of the task, which the kernel uses to update the task parameters. These informations are used to implement the scheduling algorithm, which is a variant of EDF. Details can be found in [START_REF] Lemerre | An introduction to time-constrained automata[END_REF].

2) Communication mechanisms: In addition to temporal constraints, ΨC provides different communication primitives. These primitives take advantage of the timing constraints to provide properties such as determinism, but do not have any impact on the TCA model (for instance, as they are nonblocking, they do not impact TCA scheduling).

In general, the actual PharOS implementation strictly follows the formal model, which does not impede it from addressing real-world applications.

III. ISOLATION MECHANISMS IN PharOS

Isolation is the set of mechanisms that ensures that a misbehaving task cannot alter the execution of other tasks. In this section, we introduce the spatial and temporal isolation mechanisms for a safety-oriented implementation of the TCA model described in section II.

A. Determinism

An important principle for fault isolation is determinism. It means that the application behavior is independent of the scheduling or the execution time of the tasks (as long as timing constraints are met). It implies that variations in execution time of the tasks cannot impact the behavior of other tasks, and thus contributes to temporal faults isolation. In PharOS, memory segmentation imposes that all communications happen only through dedicated mechanisms that enforce determinism. This makes impossible for an application developer to write a nondeterministic application. Whatever the application, its behavior is always fully reproducible, even if spatial or temporal faults do happen.

1) Principles: Nondeterminism must be prevented in two ways:

• by complying with the observability principle: a message is formally sent only when the d node following the sending is reached, and can be formally received only when the receiver has passed a d with d > d. This condition ensures timing consistency, i.e. that a message can be received only after it was formally sent. • by ensuring that the order in which messages are received is independent of the actual date when the messages were sent (a shared FIFO being a typical counter-example). This property can only be reached with a careful implementation of the communication primitives. 2) Design and implementation: PharOS offers several communication primitives, such as temporal variables (i.e. periodic data-flow) or messages. More examples are given in Section IV-A3. Communications between tasks in PharOS follow a strict observability principle: only temporally visible data can be used in a code block, and already visible data can not be modified.

B. Spatial isolation and fault containment

1) Principles: Spatial isolation protects tasks against unauthorized accesses to code or data, either by the other tasks, by the system (i.e. the kernel) or even by the task itself. In our approach, an application is made of a set of tasks that are statically described. Tasks are protection units whose external communications are statically and entirely defined. Based on this static knowledge, and using the available hardware for memory protection (Memory Protection Unit, MPU, or Memory Management Unit, MMU) and privilege protection (CPU protection levels), the general principles to achieve spatial isolation are:

1) Separate the application into a set of memory contexts, following the separation of privilege principle [START_REF] Saltzer | The protection of information in computer systems[END_REF]. It is more secure to have as many different memory contexts as possible. In particular, it is better to have more than one memory context per task, that changes depending on the actions it performs. 2) Generate for each memory context a constant table of the memory rights. Off-line generation of memory tables makes the kernel much simpler, and thus more secure.

The generation is performed during the linking stage. Moreover tables should have as few rights as possible, following the least privilege principle [START_REF] Saltzer | The protection of information in computer systems[END_REF]. 3) Use access control checks and memory protection to control communication between tasks. Only tasks explicitly allowed to communicate should be able to do so, following the closed-system [START_REF] Denning | Fault tolerant operating systems[END_REF] and complete mediation [START_REF] Saltzer | The protection of information in computer systems[END_REF] principles. Also, applications should be able to communicate only through the dedicated deterministic mechanisms (recall Section III-A). 4) Switch the memory rights when memory context changes (task switch, or switching of a memory context inside of a task). Memory rights switching should be performed by the kernel and the kernel only -this can be ensured with CPU protection levels. Note that memory segmentation serves two separate purposes: spatial isolation and fault detection. For instance, forbidding write access to the stack of another task is spatial isolation and fault detection, while forbidding read access to another task is only fault detection. When coping with hardware limitations (e.g. limited number of MPU descriptors), early fault detection may be sacrificed to help providing full spatial isolation, as we are about to see in the next section.

2) Design and implementation: a) Segments and layers: Our TCA-based system is made of three layers. The user layer hosts the functional code of the tasks. The system layer implements the communication mechanisms between tasks, as well as the management of state transitions within each task. Finally, the microkernel manages the time and the scheduling of tasks. Put together, the system layer and the microkernel make the kernel.

For these different layers, the memory is separated into segments with necessary and sufficient access rights depending on the execution context. That is for each physical memory segment, the memory access rights are defined for the following 4 contexts:

• the user layer;

• the system layer of the task that "owns" the segment, i.e.

that may be allowed to write to it (there is at most one writer per memory segment); • the system layer of the other tasks (multiple readers can be authorized); • the microkernel. This partitioning of each execution context is performed off-line by the toolchain of our system.

In the most general case, 2n + 1 execution contexts are defined in an application made of n tasks2 , each context defining the access rights for several memory segments. This allows detection and confinement of operational anomalies (e.g. referencing fault) between the tasks, between a task and the kernel, inside the kernel and inside the microkernel.

Note that a memory address translation mechanism is not necessary to implement the principles exposed here. In other words, a MPU is sufficient, as long as it allows defining enough memory segments to partition a memory context. Indeed, the grain of the final segmentation depends ultimately on the flexibility of the memory protection hardware. When the MPU device defines only few segment descriptors (typically 8, on our PharOS experimentation platform), different memory segments can only be protected by common descriptors. For instance, code segments of different tasks are protected by a single descriptor, with therefore the same read-only memory rights for all tasks. In every case however, we never give write access to read-only memory segments. Thus isolation is fully achieved, and only fault detection is weakened.

C. Temporal isolation 1) Principles: Temporal isolation protects tasks against timing problems. In the case of a real-time system, this protection ensures that a task cannot affect the fulfilling of timing constraints (start times and deadlines) of other tasks. For instance, in a classical priority-based real-time system without temporal isolation, an infinite loop in a high priority task can deny CPU to less priority tasks, which will miss their deadlines.

The general principles to achieve temporal isolation are: 1) Give a timing budget to each task or group of tasks, that is replenished at some predefined points; 2) Ensure that given the timing budgets, all critical tasks will meet their deadlines in every possible case (ensure this for non-critical tasks if possible); 3) Ensure that a task or group of tasks cannot execute for longer than defined by its timing budget; 4) Ensure that scheduling cannot be outside the allowed cases. Static scheduling as used in [START_REF] Arinc | Avionics application software standard interfarce. arinc specification 653[END_REF], is the simplest example of technique that ensures isolation. But it is also feasible to combine temporal isolation with dynamic scheduling.

2) Design and implementation: a) Definition of timing budgets: We attribute a timing budget to every block of code between any two successive temporal constraints. Replenish occurs when the code executes one of the system call corresponding to a node. For instance on Fig. 1, there would be a budget for b; h, one for b; c, one for f ; a, etc.

b) Feasibility analysis: We can perform an off-line feasibility analysis step using the automata of the tasks and their timing budgets. First, a synchronous product of the automata is computed, shows every possible configuration of simultaneously executing blocks of code among all tasks. From this information, a linear programming problem is extracted. Solving it tells whether every task can reach its deadline with the given initial timing budgets. Details on feasibility analysis in PharOS can be found in [START_REF] Aussaguès | A method and a technique to model and ensure timeliness in safety critical real-time systems[END_REF], [START_REF] Lemerre | Equivalence between schedule representations: Theory and applications[END_REF].

Note that feasibility analysis "only" provides the assurance that all deadlines of the systems are met, i.e. contributes to the security of the global application. Tasks can still lack some time to terminate execution within the given timing budgets.

To avoid this situation, timing budgets should be set to an upper bound of the worst-case execution time (WCET) of the corresponding block of sequential code. But WCET analysis is outside of the scope of this paper.

c) Optimal scheduling: The microkernel is considered as a trusted piece of software, whose main role is to execute the scheduler. The scheduling algorithm it implements is a variant of EDF for time-constrained automata. We proved this algorithm to be optimal [START_REF] Lemerre | An introduction to time-constrained automata[END_REF]. Optimality means that the proof of feasibility of a set of tasks ensures that EDF will always schedule tasks so that they meet their deadline (provided that they execute for no longer than their allocated timing budgets). We also have prototyped optimal scheduling algorithms for multiprocessor [START_REF] Lemerre | Equivalence between schedule representations: Theory and applications[END_REF].

d) Timing budgets and deadlines monitoring: The microkernel also monitors timing budgets online using a hardware timer: whenever a task switch occurs, a hardware timer is set to the value of the remaining budget. Upon preemption, the budget is decreased by the amount of time spent executing the task, and the hardware timer is set to monitor the new task. Thus, the timer is triggered only when a task exhausts its timing budget.

A similar mechanism is used to monitor deadline overruns. This monitoring is not strictly necessary, since timing budgets are monitored, and they should have been checked against feasibility analysis. However, it is used as a defensive programming measure.

e) Node transitions monitoring: When a system call is made because a node has been reached, the system layer checks if the transition between the previous node and the current one is allowed. This design prevents a task from following a forbidden path on its execution graph. Indeed, should such an erroneous case go undetected, the set of currently executing blocks for all the tasks would be different from the ones checked by the synchronous product. This design also protects the microkernel from system calls flooding.

Furthermore, the communication design with statically limited sending/reception rate prevents tasks from flooding each other. The set of defensive programming mechanisms implemented by the kernel ensures that undetected denial of service attack can never occur within applications.

f) Non-blocking design: A task can affect scheduling only in two ways:

• by the path choosed on conditional nodes;

• by the execution time used to reach a temporal node. There is no other way for a task to affect scheduling, for instance by yielding the CPU between nodes by blocking on a semaphore 3 . Combined with node transition and timing budget monitoring, this ensures that scheduling always happens as foreseen by the synchronous product.

g) Summary: how temporal isolation principles are applied: PharOS is a strict implementation of the temporal isolation principles based on timing budgets that we gave previously:

• Step III-C2b and III-C2c ensure that all critical tasks will meet their deadlines if, and only if, scheduling follows the cases anticipated in the feasibility analysis and timing budget is sufficient; • Step III-C2d, III-C2e and III-C2f ensure that scheduling always happens as expected by the feasibility analysis, and that tasks cannot execute for longer than specified by their timing budgets; • Ensuring that timing budgets for critical tasks are sufficient has to be proved separately, and is outside the scope of temporal isolation.

D. Group of tasks and error recovery

PharOS incorporates the concept of "groups of tasks", which is a first step towards mixed-criticality.

Tasks are partitioned into several groups. Local failure of one task (be it timing budget overrun, segmentation fault, illegal instruction, etc.) makes every task in its group to be stopped and recovered by restarting at some predefined point later in time. Failure of a task has no impact on the behavior of other tasks belonging to other groups. The time needed for recovery can be modeled as special transitions in each task graph, and can be accounted for in the feasibility analysis.

Closely-related tasks should be put together in a group. A typical example are tasks that concur together to provide a single high-level functionality, and for which there is no advantage to allow them to fail individually. All tasks in a group are guaranteed to fail and restart consistently (e.g. without time shifts).

Groups of tasks and error recovery, on top of the previous isolation mechanisms, allow for a first level of mixedcriticality. Indeed, if non-critical tasks (i.e. tasks likely to fail) and critical tasks are separated in different groups, the failure of a "non-critical group" would have no incidence on a "critical group".

IV. COHABITATION WITH OTHER EXECUTION MODELS A. Cohabitation of time-driven and event-driven applications

In embedded environment such as automotive domain, many activities (e.g. signal capture) require very tight timing constraints. Hardware performance does not allow addressing such activities with time-constrained automata. Therefore, PharOS allows coexistence of time-and event-driven tasks.

In such systems, two execution paradigms are active at the same time:

• the time-triggered paradigm where the behavior of the tasks (called TT tasks) is specified according to timeconstrained automata; • the event-triggered paradigm where the behavior of the tasks (called ET tasks) and their activation depends on the occurrence of the associated events (I/O interrupts). In a system that provides both paradigms, an important point is to ensure noninterference among them while preserving the characteristics of the previously described model: safety, dependability, temporal determinism, spatial and temporal isolation.

PharOS addresses both time-triggered and event-triggered paradigms on dual core architectures. To reach the previously described objectives PharOS provides some isolation mechanisms that ensure safe communications between ET and TT tasks, with accurate temporal and spatial isolation. Most of these properties are automatically extracted from software design and automatically generated for the runtime.

1) Spatial partitioning and isolation: The partitioning mechanisms previously described are reinforced in order to maintain protection between the ET and TT domains. Figure 3 gives a graphical description of the memory protection mechanism.

Moreover, as PharOS can run on dual-core hardware architectures, a high degree of spatial and temporal isolation between the TT and ET domains is achieved by running TT tasks on one core and ET tasks on the other core. Nevertheless, mixed execution of ET and TT on a same core is currently under study. In addition, multi-cores architectures are also supported for TT tasks and TT Tasks are dynamically assigned to available cores at runtime. TT and ET tasks could therefore also be active at the same time on multi-cores architectures (with more than 2 cores), even though this has not yet been performed.

2) Temporal partitioning and isolation: We have seen in section III-C that all PharOS TT tasks are carefully monitored during execution. For the event-triggered tasks, the temporal control is performed by an interrupt configuration monitoring set by the following parameters:

• a time interval ∆ 0 ;

• a maximum number N of authorized interrupt occurrences during the interval ∆ 0 . When the number N is reached for a given source of interrupt, any new attempt to handle the same interrupt is detected by the system. Then, the corresponding interrupt source is disabled. As a consequence, the whole PharOS application cannot be disturbed by a temporal overflow of the ET domain. Note that N and ∆ 0 are specification parameters of the behavior of event-triggered tasks, i.e an ET task is considered to be correct if and only if its number of occurrence during ∆ 0 is equal or less than N . Therefore this mechanism protects the software against external or hardware errors, and cause no harm when the system behaves normally.

3) Communications: ET tasks do not communicate with each other, but can send and receive messages from the TT tasks.

The ET→TT interface allows point-to-point communication from one ET task (the sender) to one TT task (the recipient). Data is explicitly pushed by the sender and is made visible for the recipient according to the temporal observability principles (recall section III-A). Figure 4 shows an example of ET→TT communications. An ET task pushes some data when processing a2. This data will be available for the TT task after its first activation following the "push", i.e. at T 2. Determinism of the TT tasks is preserved because a) FIFO do not receive information from multiple ET tasks, and b) the temporal observability principle ensure that behavior is not impacted by the concurrent execution of the ET and TT tasks.

The TT→ET interface allows point-to-point communication from one TT task (the producer) to one ET task (the consumer). Data are regularly updated by the system according Fig. 4. Flow from ET to TT tasks to temporal behavior of the producer. On event occurrence, the consumer accesses the latest value. Figure 5 shows such an example of TT→ET communications. A TT task produces a temporal flow X. When the event e1 occurs, the ET task observes that the value of X is 0. However, when e3 and e4 events occur, the ET task observes that the value of X is 4 in both cases.

B. Virtualization support

The main goal of virtualization is to run safely heterogeneous applications on the same platform, with reduced development costs (e.g. applications using different APIs, mixed-criticality applications, etc.). We present here the main technical challenges that OS virtualization implies on embedded platforms, then we briefly describe a prototype of a virtualized OS running under the supervision of PharOS.

A lot of techniques exist to run multiple OSes on a same CPU. We will only focus on how to run one virtualized OS (or guest) under the supervision of a host OS that ensures the execution of real-time critical tasks, such as PharOS. We will only consider the case where the guest executes noncritical tasks, or at least less critical than the ones executed by the host 4 . In this context, we use the paravirtualization technique, which means we modify the guest OS to make it runnable by the host. Obviously, the host OS also require some modifications to execute and schedule the guest: these modifications form a minimal hypervisor.

We enumerate here the key-points that ought to be carefully considered when implementing such a paravirtualization solution, then we illustrate them with a practical case based on PharOS.

1) CPU execution mode: if the CPU supports protected execution, the guest OS should not be able to run in privileged mode. That way, hardware protection ensures that the guest is not able to execute sensitive instructions that would compromise the host OS functionalities (interrupts masking, CPU halting, etc.). Instead, as a paravirtualized OS, the guest "asks" for the corresponding service to the host OS using dedicated system calls, or hypercalls. The handlers of these hypercalls are implemented as part of the hypervisor within the host kernel.

It is formally impossible to ensure any isolation between the guest and the host without the support of a protected execution mode. Likewise, it is impossible to ensure isolation between a regular OS and its tasks.

At last, note that some CPUs provide a "harware-assisted" virtualization feature [START_REF] Fisher-Ogden | Hardware support for efficient virtualization[END_REF] -especially the Intel x86 architecture. It basically consists in an additional protection level that ultimately allows the kernel and the tasks of the guest system to be isolated from each other. This separation ensures that a failing guest task can not crash the whole guest system -assuming that the guest kernel implements such protection mechanisms. Hardware-assisted virtualization also greatly helps in virtualizing un-modified OSes (full virtualization).

2) Guest scheduling -CPU resource sharing: we propose three scheduling policies for executing a guest OS on a realtime host OS.

The first policy consists in using the "idle-time": the guest is scheduled only once all real-time tasks of the host have been executed and completed. The guest can run as long as no real-time task reaches its earliest starting date. This policy dedicates all the "idle-time" of the CPU to the execution of the guest, thus providing an optimal utilization of the CPU resource. However, this is also a drawback: the CPU time allocated to the guest depends on the real-time behavior of the host application, and can potentially be very irregular. This policy is suited for executing "best-effort" guest OSes, e.g. Embedded Linux.

Another solution consists in executing the guest OS as a special real-time task of the host OS. This task is executed with a period p within a time window of duration ∆ ≤ p, and is given a CPU time budget b × ∆, where b ∈ [0, 1] is a percentage of the temporal window width (see figure 6). The guest task is scheduled as any other real-time task, but when its CPU budget is entirely consumed, the task is considered to be completed, and must wait for the next temporal window to be schedulable again. It is also possible to allow the guest OS to release the CPU before its timing budget is consumed, using a specific hypercall.

This policy gives the developer a precise control of the CPU share given to the virtual machine, equals to b•∆ p . As it turns out, it is very well suited for virtualizing a real-time guest OS on a real-time host OS. In that case, p can be wisely chosen to simulate a regular clock interrupt to the guest OS.

A last solution consists in mixing the two policies above to handle mixed-criticality among the tasks of the guest. For instance, guest tasks with soft real-time requirements can be executed within a dedicated host real-time task, but can also be given the remaining idle CPU time that will help reaching their deadlines.

3) Guest interrupts handling: when the guest is running, it may handle its own interrupts. However, when it is idle, the host must handle them on behalf of the guest, while still ensuring safeness and determinism. Once again, two policies are possible:

• Interrupts channels of the guest are active only when the guest is actually being executed. Interrupts received while the host is executing any critical task are lost. This solution is simple and safe, but the loss of interrupts is a drawback. • Interrupts for the guest are received by the host even if the guest is "asleep", but under certain conditions. A mechanism similar to the one described in section IV-A should be implemented to prevent the host from an "Interrupt Storm". This would also guarantee a bounded overhead due to guest interrupts handling. The interrupts received should be buffered and delivered to the guest when it is executed again. This policy theoretically avoids interrupts lost for the guest, however its implementation is a bit more complex than in the first policy. Note also that this solution can simulate a periodic interrupt (e.g. clock interrupt) to the guest. A dark spot though: waking-up the guest is delicate. Indeed, it will probably not endure well to be flooded by all the interrupts missed while it was asleep, and a "temporization" policy would be more appropriate. This question deserves further study and experimentation and is part of our future work. 4) Memory isolation: this feature cannot be implemented without a hardware support, i.e. without a MMU or at least a MPU device. When such a hardware memory protection unit is available, the host OS must ensure that it remains the only one able to modify memory access rights. If the guest needs to modify those rights, it should do so through a dedicated hypercall, that will check the sanity of its request. Moreover, the guest OS should only be able to access to a dedicated part of the physical memory, but not to the memory of the host or of its critical tasks.

5) Communication between the host and the guest: the ultimate integration step of the guest OS is to allow him to communicate with the host. We voluntarily will not go into the details here, because a lot of solutions can be considered: shared memory buffers, dedicated hypercalls, virtual network devices, etc. In addition, these solutions depend on the requirements set on data communications between the guest OS and the host as well as the performance of the underlying hardware.

6) Use-case -Trampoline over PharOS: A first prototype of mixed systems has been realized with both PharOS and Trampoline PharOS also provides hypercall services to Trampoline to manage its basic functionalities. Communications between host and guest OSes are not implemented yet in this prototype. However, a solution based on a specific ET task interface, as presented in section IV-A, is currently being studied.

The prototype is designed to execute safely two representative automotive applications. The first one is executed by PharOS: six real-time tasks put together a subset of a body controller. That controller includes some system output commands, communication mechanisms through a CAN-bus and sensors signal measurements. The second application, hosted by Trampoline, is basically a diagnostic function (aliveness monitoring) made of two real-time tasks based on the OSEK API.

The Trampoline task has been designed as a single PharOS task "group" (made of only one task), thus able to recover from an error and restart (see section III-D). Various error scenarii are used to stress PharOS safety mechanisms. The same tests were transposed to the Trampoline application, all resulting in successful detection and isolation by PharOS. The two main tests are:

• Illegal memory access: an attempt to read or write to PharOS memory space from Trampoline is detected and Trampoline is restarted. • Infinite loop: an erroneous Trampoline task entering an infinite loop does not interfere with PharOS tasks scheduling. It is the responsibility of the Trampoline kernel however to detect this fault, and possibly restart the task. As a conclusion, let us stress that we have focused on one paravirtualization technique through a practical use-case, but a lot of alternatives can be considered. The choice of a solution should depend on the type of the virtualized OS, and mostly on the underlying architecture. Typically, multi-core architectures offer many different ways to solve the CPU-sharing problemthe simplest but clearly suboptimal solution being to dedicate one core for each guest.

V. RELATED WORK Spatial and temporal isolation mechanisms have been introduced to support the virtualization of the execution of multiple operating systems on a single platform [START_REF] Meyer | A virtual machine time-sharing system[END_REF], [START_REF] Borden | Multiple operating systems on one processor complex[END_REF], a long time ago. Since then, protection of information within operating systems has been the subject of numerous work, such as [START_REF] Denning | Fault tolerant operating systems[END_REF], [START_REF] Saltzer | The protection of information in computer systems[END_REF] to name a few. Lately, the avionic domain has repopularized and defined the principles of spatial and temporal partitioning [START_REF] Rushby | Partitioning in avionics architectures: Requirements[END_REF], while Xen [START_REF] Barham | Xen and the art of virtualization[END_REF] has also re-popularized research work in the area of virtualization. However, fewer work have focused on issues of virtualization within real-time constraints [START_REF] Henning Schild | Faithful virtualization on a realtime operating system[END_REF], especially in highly-constrainted memory environments.

In [START_REF] Rushby | Partitioning in avionics architectures: Requirements[END_REF], underlying mechanisms and issues to implement spatial and temporal partitioning have been shown. Classically, spatial partitioning relies on hardware memory protection units. However, a software approach has also been proposed called Software Fault Isolation [START_REF] Wahbe | Efficient software-based fault isolation[END_REF] (SFI). In this approach, memory references are statically checked and, when it is possible, additional code is embedded in order to check accesses during execution. Proof-carrying code (PCC) [START_REF] Necula | Safe kernel extensions without run-time checking[END_REF] can be seen as a generalization of SFI, in which a safety policy is expressed using first-order logic that is certified and verified by a kernel prior to execution. But the problem of these static-checking methods is that they do not provide any protection against hardware bugs, nor do they allow execution of unchecked third-party code, which PharOS does.

Currently, several operating systems implements spatial and temporal partitioning such as DECOS [START_REF] Obermaisser | Decos: an integrated time-triggered architecture[END_REF] and successors or LynxOS-178 [START_REF] Lynuxworks | LynxOS-178. Certifiable RTOS for safety-critical computing[END_REF], an ARINC 653 compliant OS [START_REF] Arinc | Avionics application software standard interfarce. arinc specification 653[END_REF]. However, in both cases, temporal partitioning is ensured by static scheduling of partition time slots.

VI. CONCLUSION

The PharOS model, method and tools provide a strong basis for real-time systems isolation. We have explained how we combined temporal isolation with optimal scheduling and dynamic task behavior, using offline compilation and analysis, and online monitoring techniques. We showed how we achieved spatial isolation by automatic generation of finegrained memory tables at compile time, and control of communication using a trusted "system layer" in a separate address space. Using these isolation principles as a basis, we have presented many different ways of extending PharOS to support several levels of mixed-criticality:

• error confinement and recovery with groups of tasks;

• interrupts processing with appropriate monitoring and physical core separation; • paravirtualization of external OSes Future work on PharOS will mainly focus on improving resource-sharing techniques for paravirtualization (CPU cores, interrupts channels, peripherals), and on new scheduling algorithms on multi-core embedded platforms, based on our previous work [START_REF] Lemerre | Equivalence between schedule representations: Theory and applications[END_REF]. In addition, we plan to improve the safeness of our approach by making the system layer untrusted. The node transitions checking would be handled by the microkernel (as critical activities), and the system layer would ultimately be allowed to fail without impacting other tasks. Such a solution would also allow much more flexibility and facilitate the virtualization of communicating tasks.

Fig. 1 .

 1 Fig. 1. A ψC code excerpt and corresponding automaton. An arc of label l represents the block of code to which label l belongs.

Fig. 2 .

 2 Fig. 2. Example of a correct execution for the automaton in Fig. 1.

Fig. 3 .

 3 Fig. 3. Spatial isolation between ET and TT domains

Fig. 5 .

 5 Fig. 5. Flow from TT to ET tasks

Fig. 6 .

 6 Fig. 6. Temporal parameters of the execution of the guest OS. The light gray time zone is the temporal execution window; the black zones are actual execution times. On this example, the CPU time consumed is x × ∆, with x ≤ b.

[ 1 ]

 1 OSes, used respectively for high and low critical functions on a Freescale S12XE architecture. Trampoline is an open source Real-Time Operating System, that implements the OSEK/VDX API, a standard in automotive embedded software. This prototype is therefore a rather straightforward solution to feature OSEK/VDX tasks on a PharOS driven platform.The paravirtualization of Trampoline is ensured by the PharOS microkernel. In accordance with the principles exposed previously, Trampoline is encapsulated in a timetriggered (TT) task and executed in protected mode. This guest TT task has a period of 10 ms, which corresponds to the elementary time quantum of the Trampoline real-time kernel. Using the notations from figure6, the parameters of the guest task are: p = 10ms, ∆ = 10ms, b = 10%. In other words, Trampoline is given 10% of the total CPU resource. Every time the Trampoline task is awaken by the scheduler, the mini-hypervisor implemented in the PharOS kernel simulates a clock interrupt for Trampoline. This regular "tick" indicates that 10ms have elapsed in the real world.

In the actual version of PharOS, there are no d nodes, only d and d nodes are allowed.

Per task, one context for its user layer and one context for its system layer, plus one dedicated context for the microkernel.

Actually, the only provided mean of synchronization in the TCA model is time: PharOS provides no semaphores or mutexes whatsoever.

This model complies to a intuitive and simple rule: the OS that is able to execute the most critical tasks should be the one in full control of the hardware.

ACKNOWLEDGMENT

The authors would like to thank Vincent DAVID (CEA LIST), for its pioneer works on OASIS [6], that lead to the current implementation of PharOS.

This research on PharOS is supported in part by Agence Nationale de la Recherche (ANR), through the SCARLET project. The paravirtualization prototype of Trampoline over PharOS was especially developed within this project.


J Kuzmanis
email: yanisku@gmail.com

A simple solvability criterion for the negative Pell equation

Keywords: Positive/negative Pell equation, quadratic irrationalities. AMS Classification: Primary 11D09, Secondary 11R29

Provided solution of positive Pell equation is known, a simple approach to Mollins/Srinivasan criterion for the solution existence of the corresponding negative Pell equation is proposed.

Introduction

The Pell equation is fascinating mathematicians by its seeming simplicity and many deep connections to various theoretical problems. Therefore almost all introductory number theory courses include chapters on Pell equation; we recommend [START_REF] Nagell | Introduction to Number Theory[END_REF].

While the positive Pell equation

a 2 -Db 2 = 1 (1)
always has infinitely many non-trivial solutions (a, b) for every non-square natural number D, this is not the case for the negative Pell equation

x 2 -Dy 2 = -1. (2)
Theory states the following crucial moments concerning positive and negative Pell equations.

(1) Equation ( 1) has a least positive solution, where variables a > 0 and b > 0 have their least positive integer values; it is often named as fundamental solution (a 0 , b 0 ). The same can be told about equation (2) -when it is solvable; here we denote its fundamental solution as (x 0 , y 0 ).

(2) All (infinitely many) solutions (a, b) of (1) for one fixed D value are obtained from its fundamental solution (a 0 , b 0 ) as rational and irrational parts of the expression

a + b √ D = (a 0 + b 0 √ D) k ,
where k = 1, 2, etc.

(3) All (infinitely many) solutions (x, y) of ( 2) are obtained analogously from its fundamental solution (x 0 , y 0 ), but only odd powers of k are permitted.

(4) Provided equation ( 2) is solvable for given particular D and its fundamental solution is (x 0 , y 0 ), the fundamental solution (a 0 , b 0 ) of equation ( 1) for the same D can be obtained as

a 0 + b 0 √ D = (x 0 + y 0 √ D) 2 . (3)
The intriguing point is that there are many non-square D values, for which equation ( 2) is insolvable. One well-known limitation is D ≡ 1, 2 (mod 4) with no 4n + 3 factors (see [START_REF] Nagell | Introduction to Number Theory[END_REF], page 201), but this is not sufficient. Sequence A 031398 from oeis.org includes square-free D values with no 4n + 3 factors, for which negative Pell equation is insolvable. There are known several methods of determining the solvability of (2) and we can mention [START_REF] Perron | Die Lehre von den Kettenbrüchen[END_REF], [START_REF] Matthews | Primitive Pythagorean triples and the negative Pell equation[END_REF] and [START_REF] Grytchuk | The Negative Pell Equation and Pythagorean triples[END_REF], but they include careful analysis of continued fraction expansions of √ D. On 2010 an excellent elementary criterion was published by R.A. Mollin and A. Srinivasan [START_REF] Mollin | A Note on the Negative Pell Equation[END_REF], but the way to this number theory gem also deals with not so simple treatment of continued fractions [START_REF] Mollin | Pell equations: non-principal Lagrange criteria and central norms[END_REF].

Proposed approach

In the following we present a short and elementary approach to solvability criterion for negative Pell equation [START_REF] Perron | Die Lehre von den Kettenbrüchen[END_REF]. Given (a 0 , b 0 ) is obtained as fundamental solution of equation ( 1) and taking into account (3), we can seek for unknowns (x 0 , y 0 ) by standard algebraic procedure. From (3):

a 0 + b 0 √ D = (x 0 + y 0 √ D) 2 = x 2 0 + Dy 2 0 + 2x 0 y 0 √ D.
We equalize rational and irrational parts: Inserting this in the first equation of the system, we get

   a 0 = x 2 0 +
4a 0 x 2 0 = 4x 4 0 + Db 2 0 or 4x 4 0 -4a 0 x 2 0 + Db 2 0 = 0, or (2x 2 0 ) 2 -2a 0 (2x 2 0 ) + Db 2 0 = 0. Thus 2x 2 0 = 2a 0 ± 4a 2 0 -4Db 2 0 2 = a 0 ± a 2 0 -Db 2 0 = a 0 ± 1.
As we are seeking for fundamental solution, we take minimal positive value, so

x 0 = a 0 -1 2 and y 0 = b 0 2x 0 must be natural numbers.
Are these conditions sufficient?

x 2 0 -Dy 2 0 = a 0 -1 2 -D b 2 0 4x 2 0 = a 0 -1 2 -D b 2 0 2(a 0 -1) = = a 2 0 -2a 0 + 1 -Db 2 0 2(a 0 -1) = 1 -2a 0 + 1 2(a 0 -1) = -1
So for the non-square natural D ≡ 1, 2 (mod 4) we get a solution of negative Pell equation (2) with roots x 0 = a 0 -1 2 and y 0 = b 0 2x 0 as natural numbers; here (a 0 , b 0 ) is the fundamental solution of corresponding positive Pell equation [START_REF] Nagell | Introduction to Number Theory[END_REF]. Proposed approach easily gives Mollins/Srinivasan criterion a 0 ≡ -1 (mod 2D) from [START_REF] Mollin | A Note on the Negative Pell Equation[END_REF], because this criterion easily follows from a 0 = x 2 0 + Dy 2 0 and (2). In both cases the fundamental solution (a 0 , b 0 ) of corresponding positive Pell equation ( 1) is necessary.

At the end three illustrative examples. Negative Pell equation ( 2) is solvable for D = 2, 5, 10, 13, because we can get natural number values for (x 0 , y 0 ).

D a 0 + b 0 √ D x 0 + y 0 √ D 2 3 + 2 √ 2 1 + √ 2 5 9 + 4 √ 5

  The same is for D values with 4n + 3 factors.

	2 + 3 + √ √ 10 5 13 18 + 5 10 √ √ 13 649 + 180 10 19 + 6 √ 13
	For D values from A 031398 irrationality in the third column looks different.
	D		a 0 + b 0	√ D	x 0 + y 0	√ D
	34		35 + 6	√	34	√ 17 (17 + 3 17	√ 34)
	146	145 + 12	√ 146	√ 2 (12 + 2	√	146)
	178 1601 + 120	√ 178	√ 2 (40 + 3 2	√	178)
	194	195 + 14	√ 194	√ 97 (97 + 7 97	√	194)
	205 39689 + 2772	√ 205	√ 41 (902 + 63 41	√	205)
	221 1665 + 112	√ 221	√ 13 (104 + 7 13	√ 221)
	D a 0 + b 0	√	D	x 0 + y 0	√ D
	3	2 +	√	3	√ 2 (1 + 2	√ 3)
	7	8 + 3	√ 7	√ 14 (7 + 3 14	√	7)
	15	4 +	√ 15	√ 6 (3 + 6	√	15)
	18 17 + 4	√ 18	√ 2 (4 + 2	√	18)


Peipei Zhou

Dominique Tarlet
email: dominique.tarlet@univ-nantes.fr

Min Wei

Yilin Fan

Lingai Luo

Novel multi-scale parallel mini-channel contactor for monodisperse water-in-oil emulsification

Keywords:

This paper presents the development of a novel parallel minichannel contactor for monodisperse water-in-oil emulsification. The novel contactor having 16 Y-junctions in parallel equipped with multi-scale tree-like fluid distributor and collector was realized and experimentally tested. Oils with different physical properties and water were used as continuous phase and dispersed phase, respectively. Influences of different factors including the parallelization, the water/oil flow-rate ratio and two-phase physical properties on the two-phase flow characteristics and on the coefficient of variation were studied and discussed. Results

showed that for water-pure oil flow, the novel parallel minichannel contactor is capable of producing monodisperse water-in-oil emulsions with high throughput and low energy consumption under our tested conditions. The Ohnesorge number (Oh) has a significant influence on the emulsion monodispersity and the capillary number (Ca) appears as a proper indicator for liquid-liquid break-up regimes. A predictive correlation for the Oh number was proposed by fitting the experimental data obtained in the study.

Introduction

Microfluidic devices have recently received great attentions in the field of emulsification process because of the ability to generate uniform-sized droplets and plugs, also known as monodisperse emulsion. As an important step to produce micromaterials (e.g., microparticles and microcapsules), monodisperse emulsions are intensively required in numerous applications such as microreaction synthesis [START_REF] Meng | Plug-n-play microfluidic systems from flexible assembly of glass-based flow-control modules[END_REF], drug delivery [START_REF] Duncanson | Microfluidic synthesis of advanced microparticles for encapsulation and controlled release[END_REF] and the production of colloid particles [START_REF] Jeong | Controlled generation of submicron emulsion droplets via highly stable tip-streaming mode in microfluidic devices[END_REF]. Many researchers studied systematically the emulsification phenomena in a single micro-channel, including the formation of droplets [START_REF] Link | Geometrically mediated breakup of drops in microfluidic devices[END_REF], the size distribution [START_REF] Gunther | Multiphase microfluidics: from flow characteristics to chemical and materials synthesis[END_REF] and the influence of physical properties on emulsion dispersity [START_REF] Nisisako | Microfluidic large-scale integration on a chip for mass production of monodisperse droplets and particles[END_REF].

Nomenclature

Ca

Capillary number (-) D w average diameter (mm) E energy dissipation rate (J kg -1 ) I range for average value (-) J o oil superficial velocity (m s -1 ) J w water superficial velocity (m s -1 ) n number of droplets (-) Oh

Ohnesorge number (-) P injection pressure (Pa) P w coefficient of variation (%) Q volume flow rate (mL min -1 ) (m 3 s -1 ) Re

Reynolds number (-) S section area (m 2 ) u (-) normal law (-) W,Z channel width (m) the emulsion throughput, the rate of increase depending on the ratio of parallelization. However, the fluid flow maldistribution among parallel channels might significantly reduce the emulsion monodispersity. [START_REF] Li | Simultaneous generation of droplets with different dimensions in parallel integrated microfluidic droplet generators[END_REF] confirmed that the size of droplets varied with the flow-rate ratio of two liquids.

Greek
Even though the uniform droplets can be obtained in each single channel, the mixture of emulsions produced from different channels will cause the reduced emulsion monodispersity.

Various attempts were then made to tackle this problem. For example, [START_REF] Kobayashi | Formulation of monodisperse emulsions using submicron-channel arrays[END_REF] proposed a tapered shape before the parallel channels to increase the flow-rate uniformity of each single phase. The width of total inlet is decreased from 22.5 mm (initial inlet) to 14 mm (the cross-section of the inlets of parallel channels). 1500 channels are ranged in such tiny plate, with identical width of 3.2 m for each channel.

To achieve the uniform flow-rate distribution among channels, the size of plate cannot be too large, which in turn also limits the productivity. [START_REF] Nisisako | Microfluidic large-scale integration on a chip for mass production of monodisperse droplets and particles[END_REF] reported a new microfluidic module for mass production of monodisperse droplets and particles. The uniformity of pressure head was realized by the parallelization of 128 cross-junctions circularly arranged on a chip. However, the dimensions of holes (or tubes) were not optimized so that extra pressure drop might occur, leading to higher pumping power consumption.

Recently, a new multi-scale tree-like fluid distributor has been proposed by [START_REF] Tondeur | Design and scaling laws of ramified fluid distributors by the constructal approach[END_REF]. It is composed of multi-scale bifurcated internal channels which are optimized under certain constraints, such as minimal viscous dissipation and constant void volume. Numerical and experimental investigations (Luo et al., 2007;[START_REF] Fan | Flow distribution property of the constructal distributor and heat transfer intensification in a mini heat exchanger[END_REF] confirmed that the multi-scale tree-like distributor could provide better flow distribution uniformity preferentially for laminar flow, which usually fits the working condition for emulsification process. This promising concept has been widely applied in thermal, fluidic and reactive devices (Luo et al., 2007) such as liquid-liquid mixers [START_REF] Fan | Evaluation of the performance of a construcal mixer with the iodide-iodate reaction system[END_REF][START_REF] Tondeur | Uniform flows in rectangular lattice networks[END_REF], heat exchangers [START_REF] Fan | Flow distribution property of the constructal distributor and heat transfer intensification in a mini heat exchanger[END_REF]Luo et al., 2007[START_REF] Luo | Experimental study of constructal distributor for flow equidistribution in a mini crossflow heat exchanger (mche)[END_REF], gas-liquid contactors [START_REF] Yue | Flow distribution and mass transfer in a parallel microchannel contactor integrated with constructal distributor[END_REF], multifunctional reactors [START_REF] Guo | Mixing performance assessment of a multi-channel mini heat exchanger reactor with arborescent distributor and collector[END_REF](Guo et al., , 2014), etc. ), etc. In this study, we try to apply the multi-scale tree-like design for the water-oil two-phase emulsification, in order to obtain monodisperse emulsion with high productivity and low energy consumption. For this purpose, a novel parallel mini-channel contactor having 16 Y-junctions equipped with multi-scale tree-like distributor/collector was developed and experimentally tested. To further investigate the local twophase flow patterns at the junction, a single Y-junction replica was also fabricated and tested for comparison.

In the rest of this paper, we shall first introduce the experimental set-up used in this study. Then the multi-scale network with 16 parallel mini-channels and the single Y-junction replica will be described in detail. After that, the flow patterns in the single Y-junction mini-channel replica will be discussed. Then the water-oil emulsification in the multi-scale contactor will be systemically studied, including the influence of different factors (parallelization, the physical properties and flow rates) on the emulsification and its coefficient of variation. The energy dissipation rate will also be calculated and compared to that of other emulsification devices reported in the literature. A predictive correlation based on non-dimensional groups will be proposed by fitting the experimental data. Finally, main conclusions and future perspectives will be summarized.

Methodologies

Experimental test-rig

As schematically shown in Fig. 1, the experimental setup consists of test channels, oil and water circuits and the data acquisition system. Water was fed in the channel by a FRESE-NIUS Injectomat Agilia syringe pusher, with 3% uncertainty on the water flow-rate. Pressurized air (1.6-3.0 bar) regulator and needle valve were used to control the flow rate of oil (less than 5% variation within 5 minutes). To visualize and record the two-phase flow patterns in minichannels, a fast CCD camera LAVision HighSpeed Star 6 equipped with the macro objective was used. The focus was made in the channel depth. Background illumination was provided by a 50W halogen spot located in the line of sight, right above the transparent test section.

Test section

Fig. 2a shows the designed parallel minichannel contactor which is composed of three parts: multi-scale tree-like fluid distributor for oil and for water respectively, 16 Y-junctions for water-oil contact and multi-scale tree-like collector for emulsion collection. The tree-like fluid distributor has four bifurcating generations (2 4 = 16), connecting single inlet port to 16 end-ports. Its detailed geometric dimensions are listed in Table 2. The big difference in viscosities enhance mixing. When Ca is low, interfacial tension dominates to produce sharp plugs. [START_REF] Garstecki | Formation of droplets and bubbles in a microfluidic t-junction-scaling and mechanism of break-up[END_REF] Water-Oil Single T-junction 100 m × 33 m

Re 1, Ca < 0.01, Eq. ( 5) can be used to predict the plug length. [START_REF] Steegmans | Characterization of emulsification at flat microchannel y junctions[END_REF] Milli

-Q-Water Ethanol-Water Single Y-junction 23 m × 5 m
Droplet size in Y-junction is mainly determined by Ca of the continuous phase. [START_REF] Link | Geometrically mediated breakup of drops in microfluidic devices[END_REF] Plug emulsions A series of three T-junctions satisfying l 0 w 0 > 1 (l 0 : plug length, w 0 : channel width.) Sequential passive break-up to make up smaller droplets with low coefficient of variation (<3%). [START_REF] Li | Simultaneous generation of droplets with different dimensions in parallel integrated microfluidic droplet generators[END_REF] Deionized Water-Oil Four parallel flowfocusing devices

The device produces monodisperse emulsions (<4%) and the droplet size increases with the wider orifices. [START_REF] Kobayashi | Formulation of monodisperse emulsions using submicron-channel arrays[END_REF] Oil-Water Array 1500 channels 22.5 mm × 8 mm Minimum coefficient of variation (CV) 9%; droplet production rates is high, reaching up to 2.4 × 10 4 droplets per second. [START_REF] Nisisako | Microfluidic large-scale integration on a chip for mass production of monodisperse droplets and particles[END_REF] Water-Oil A glass chip with 128 cross-junctions or co-flow geometries (42 mm × 42 mm)

The cross-junction model and co-flowing model produce 320 ml/h emulsions with low coefficient of variation of 1.3% and 3.3%, respectively. [START_REF] Kashid | Numbering-up and mass transfer studies of liquid-liquid two-phase microstructured reactors[END_REF] Kerosene-Acetic Acid-Water

Numbering-up of six 120 o Y-shape mixing elements.

Influence of interfacial tension on the mass transfer performance. It obeys the scaling relation that accounts for both the fluid equidistribution among the end-ports and the minimized pressure drop [START_REF] Tondeur | Constructal optimization of arborescent structures with flow singularities[END_REF]. Through their proper fluid distributor, the introduced oil and water flows are divided into 16 sub-streams each, which are then guided into the Yjunctions to produce emulsions. The 16 Y-junctions are spaced 4 mm between the axis of one junction and another. In more detail, each Y-junction has a rectangular cross-section (Fig. 3a) of 1.0 mm in width and 0.6 mm in depth, with a contact angle of 60 • (shown in Fig. 3b). The endpoints of the 16 Yjunctions are connected to the fluid collector which has the same geometry as the distributor. Using this split-contactrecombine technique, the parallel minichannel contactor is expected to produce high throughput of water-in-oil emulsion (compared to that in single minichannel) with low coefficient of variation.

The smallest channel dimensions of the present device are 600 m in depth and 1 mm in width, resulting in a hydraulic diameter of 750 m. Considering the larger channel sizes in the distributor/collector (Table 2), they are called minichannels rather than microchannels [START_REF] Kandlikar | Evolution of microchannel flow passages -thermohydraulic performance and fabrication technology[END_REF][START_REF] Zhao | Two-phase microfluidic flows[END_REF]. Generally, two-phase flows are confined in small channels whose hydraulic diameter is smaller than the capillary length [START_REF] Aarts | Capillary length in a fluid-fluid demixed colloid-polymer mixture[END_REF] so that interfacial tension dominates [START_REF] De Gennes | Gouttes, bulles, perles et ondes[END_REF].

The integral contactor prototype (Fig. 2b) was fabricated with transparent resin (WaterClear Ultra 10122) using laser stereolithography. In this process, photopolymer is cured selectively at the local focus of ultraviolet laser. With post processing, the quality of final product (surface planarity) is guaranteed.

A single Y-junction replica (shown in Fig. 3) is also constructed to investigate the local two-phase flow patterns at the junction, where the oil and water meet. It was fabricated in-house, using the milling technique in two PMMA rectangular blocks: one with milled Y-junction, the other one being a cover. They were stuck against each other using screws, and tested without leakage. The measurements of droplet shape/size and emulsion coefficient of variation in this single Y-junction replica under different tested conditions will provide a reference case (ideal condition) for the comparison with those in the parallel minichannel contactor.

Wettability properties are tested by measuring the contact angle of a sessile drop on both PMMA and photopolymer surfaces. They are found to be similar. The average roughness of the PMMA surface mechanically etched is less than 1 m [START_REF] Belkadi | Water-in-oil emulsification in a microfluidic impinging flow at high capillary numbers[END_REF] whereas the maximum roughness of the photopolymer surface is estimated to be about 10 m [START_REF] Campbell | Surface roughness visualisation for rapid prototyping models[END_REF].

Working fluids

Tap water (Nantes métropole) and filtered sunflower oil (Carrefour) were used as working fluids. In order to modify the interfacial tension and viscosity of the two-phase fluid, additives of Butanol (Sigma-Aldrich) and SPAN83 (sorbitan sesquioleate, Fisher Scientific) were used. The chosen SPAN83 mass percentage is up to 0.4 wt.% since its critical micellar concentration (CMC) is from 0.8 to 1.0 wt.% [START_REF] Belkadi | Water-in-oil emulsification in a microfluidic impinging flow at high capillary numbers[END_REF]. Using these products and formulation, the obtained water-oil dispersion was observed to be stable for more than 50 minutes with a small diameter (5-10 m) [START_REF] Tarlet | Distribution of thermal energy of child-droplets issued from an optimal micro-explosion[END_REF]. Stability remains significant with more than about 10 minutes in the present experiments, at a larger mean diameter around 1 mm. Thus, it can still be considered as an emulsion in the present work.

Parameters measurement and uncertainty analysis

Interfacial tension w/o and viscosity were measured using a KRÜSS K-12 tensiometer and a LAMY RM-200 rheometer with their uncertainty of 5% and 1%, respectively. Water flow-rate Q w was controlled by the syringe pump (0-2 ml/min), with 3% uncertainty of the full scale. Oil flow-rate Q o was measured by means of weighing scales of Sartorius-MSE (1 Hz sampling, accuracy of 0.001 g). Concerning oil flow rate, the uncertainty is 5% within a time scale of 5 min for each measurement. Inlet pressures of water and oil were measured by using compact pressure transmitters of Gems 3100 series (0-2.5 bar with an uncertainty of 0.25% of the full scale). The pressure at the outlet of the test section was atmosphere pressure. K-type thermocouples were used to measure the fluid temperature at different locations (see Fig. 1). Details about the instruments and their measuring ranges and precisions are presented in Table 4.

Different parameters and dimensionless numbers are then introduced to describe the characteristics of water-in-oil emulsions produced both in single Y-junction and in parallel minichannel contactor. The coefficient of variation P w of the obtained emulsions is defined in Eq. ( 1). It is also called the polydispersity of emulsions [START_REF] Link | Geometrically mediated breakup of drops in microfluidic devices[END_REF]. P w is the standard deviation w divided by their mean diameter D w of more than 15 water droplets in single Y-junction, or more than 100 water droplets in the outcoming area of parallel minichannel contactor (as shown in red frame in Fig. 12).

P w = w D w (1)
The sampling size (15 and 100) is determined based on Eq. ( 2), where w is the standard deviation (of water droplets diameters), D w is the mean diameter and I is the reliable range for the mean diameter (I = 2.5% in this study). With a low level of assumed risk ˛r = 1%, the normal law gives u(0.995) = 2.58. As a result, the minimum sampling size n min calculated is always smaller than 100. Thus, the present average diameters D w are based on samples of 100 droplets as it was done in [START_REF] Li | Simultaneous generation of droplets with different dimensions in parallel integrated microfluidic droplet generators[END_REF]. In this study, the mean diameter D w of water droplets were obtained statistically by using image processing software Image-Pro Plus 7.0. The image resolution rate was set as 30 pixels per millimeter, which gives an estimated uncertainty of 6.7% considering 2 pixels of maximum uncertainty. For some non-spherical droplets, the equivalent diameter was measured that corresponds to the diameter of a sphere with the same projected area as the droplet.

c = w D w ⇒ n min = c u(1 -˛r/2) I (2)
Dimensionless numbers such as capillary number (Ca) and Ohnesorge number (Oh) are also introduced for describing the liquid-liquid two-phase flow. Ca number is widely known as the criterion to discern among different break-up regimes in the liquid-liquid flows [START_REF] Zhao | Two-phase microfluidic flows[END_REF]. It compares the viscous shear stress within the continuous phase to the interfacial tension, as defined in Eq. (3) where J o is the superficial velocity of the continuous phase.

Ca = J o o w/o (3)
Meanwhile, the Oh number relates the viscous force to the inertial and interfacial tension force of the dispersed phase [START_REF] Belkadi | Water-in-oil emulsification in a microfluidic impinging flow at high capillary numbers[END_REF], as defined in Eq. ( 4).

Oh = o ( w D w w/o ) 0.5 (4)
For convenience of comparison, all the tests were carried out between 21.3 • C and 22.9 • C to avoid big variation of fluid physical properties with temperature. Viscosity is measured at 22 • C to better fit the temperature range of the experiments. The values of density and interfacial tension w/o are measured at 25 • C since they do not significantly vary over this range of temperature.

For each measurement, the two-phase flow was stable over more than 5 minutes and good reproducibility was observed. The capillary number (Ca) varies between 0.004 and 0.08, the Reynolds number (Re) between 0.0033 and 0.45, and the Ohnesorge number (Oh) between 0.3 and 0.7 under the tested conditions (cf. Table 5).

In the present work, the uncertainty analysis for the derived quantities was performed following the method of [START_REF] Moffat | Describing the uncertainties in experimental results[END_REF]. The calculation indicates the maximum uncertainty of 11% for superficial velocities J o and J w , 12% for the capillary number (Ca), 4.3% for the Ohnesorge number (Oh) and 20% for the total energy dissipation rate (E). It should be claimed here that these estimations are conservative because most uncertainties between channels generated by imperfections and irregularities in fabrication were considered to be 50 m for channel size in the parallel minichannel contactor.

Flow patterns in a single Y-junction minichannel replica

The hydrodynamic characteristics of the water-oil flow have been investigated in the single Y-junction minichannel replica. Typical liquid-liquid flow patterns such as droplet flow or plug flow are generally reported in the existing literature. Fig. 4 shows the observed flow patterns in the present Y-junction. Characterization of the observed flow patterns is discussed as follows.

In minichannel, droplet flow is characterized by water droplets dispersed in the continuous oil phase, of which the shape is close to spherical with a diameter smaller than or close to the width of the channel. Plug flow is characterized by elongated water plug with axial length clearly larger than the channel diameter and separated by the oil plugs in the axial direction, and usually the lengths of water and oil plugs are constant. The pictures of Expts. 11, 25 and 31 in Fig. 4 show plugs with their length more than 20% longer than the channel width. Consequently, they belong to the characteristic plug flow regime. Photography of Expt. 33 in Fig. 4 shows the droplets having their size +/-15% close to the channel width. However it is not strictly droplet flow since the water droplets are not really spherical because of the rectangular cross-section of the mini-channel. We call them "apparent droplet flow" as being a transitional pattern between the characterized droplet flow (not observed in the current study) and the characteristic plug flow. The characteristic droplet flow pattern as reported in the literature is not observed under the tested conditions of this study. Fig. 5 shows the two-phase flow regime map for the single Y-junction replica as a function of oil and water superficial velocities. Two groups can be seen: experiments 1 to 14 use pure oil, and experiments 15 to 35 use oil with 5% Butanol. The transition line concerns the cases with 5% Butanol. The cases with pure oil are different due to a significant difference in viscosity and interfacial tension (see Table 3). The transition line from droplet flow to plug flow proposed by [START_REF] Kim | Liquid-liquid segmented flows in polycarbonate microchannels with cross-sectional expansions[END_REF], and the transition from apparent droplet flow to characterized plug flow are both indicated in Fig. 5 for comparison. It may be observed that most of the present liquid-liquid flows belong to the plug flow regime, according to [START_REF] Kim | Liquid-liquid segmented flows in polycarbonate microchannels with cross-sectional expansions[END_REF]. However, this criterion can be discussed mainly because of two reasons. Firstly, the break-up of the dispersed phase is based on a Y-junction here, instead of a T-junction by [START_REF] Kim | Liquid-liquid segmented flows in polycarbonate microchannels with cross-sectional expansions[END_REF]. Secondly, the dimensions of the channel with rectangular cross-section and the two-phase flow conditions are not exactly the same.

It is further observed that there is a critical value of Ca number that distinguishes different break-up mechanisms of droplet/slug formation. This critical value is typically reported equal to 0.01 in the literature [START_REF] Zhao | Two-phase microfluidic flows[END_REF]. When Ca ≤ 0.01, the interfacial tension force is stronger than the viscous shear stress. It is known as the squeezing regime, as explained by [START_REF] Garstecki | Formation of droplets and bubbles in a microfluidic t-junction-scaling and mechanism of break-up[END_REF] that proposed the mechanism of squeeze stage and growth stage for the formation of final plug in a T-junction minichannel. At first, the dispersed phase enters and fills the main channel. In this stage, the plug length is L 1 = εz, where z is the width of the main channel and ε is usually equal to 1. Then the plug neck with characteristic width d decreases at a rate equal to 2006)) to predict the length of plugs is:

L z = ı Q w Q o + ε (5)
where ı is a dimensionless parameter in relation to the minichannel geometry and dimensions.

In this study, we found that this mechanism proposed for T-junction [START_REF] Garstecki | Formation of droplets and bubbles in a microfluidic t-junction-scaling and mechanism of break-up[END_REF] is logically applicable to Yjunction minichannel. Figure 6 shows the experimental data obtained in this study and the fitting curves using the correlation of [START_REF] Garstecki | Formation of droplets and bubbles in a microfluidic t-junction-scaling and mechanism of break-up[END_REF] (Eq. ( 5)). The experimental data of [START_REF] Garstecki | Formation of droplets and bubbles in a microfluidic t-junction-scaling and mechanism of break-up[END_REF] and [START_REF] Xu | Correlations of droplet formation in t-junction microfluidic devices: from squeezing to dripping[END_REF] for rectangular T-junction microdevices are also introduced for comparison. It can be observed that the predicting results (ı = 2, ε = 1 for the current study) fit well the experimental data, with 9.7% of mean absolute percentage error (MAPE) concerning the squeezing regime (Ca < 0.01). MAPE is a parameter defined in Eq. ( 6) which is often used to estimate the deviation between the measured (M i ) and predicted values (P i ) [START_REF] Kim | Liquid-liquid segmented flows in polycarbonate microchannels with cross-sectional expansions[END_REF].

MAPE(%) = 1 N n i=1 100 × M i -P i M i (6)
Fig. 6 also shows experimental results for 0.08 > Ca > 0.01. It may be observed that the droplet lengths follow a linear regression with water/oil flow-rate ratio, implying that the flow regime can also be well described by Garstecki model which was proposed for Ca number typically smaller than 0.01 [START_REF] Garstecki | Formation of droplets and bubbles in a microfluidic t-junction-scaling and mechanism of break-up[END_REF]. But comparing to the datapoints for Ca < 0.01 of the current study, the slope is reduced by a factor of 2 while the y-intercept remains as the same (ı = 1 ; = 1). Since the channel geometry/dimension is the same for both working ranges, the difference in slope (ı value) is mainly due to the influence of physical properties. The datapoints for Ca > 0.01 were obtained by addition of 5% Butanol, so the viscosity o and interfacial tension w/o were decreased by 24% and 54% respectively compared to those of pure oil.

In the meantime, we found that the droplet lengths obtained for 0.08 > Ca > 0.01 were related to the capillary number Ca of the continuous phase, as proposed by [START_REF] Steegmans | Characterization of emulsification at flat microchannel y junctions[END_REF]:

L z = ˛+ ˇCa -0.5 (7)
where ˛is the intercept, which in most cases turns into zero and ˇis the fitting parameter in relation to mini-channel geometry/dimensions. Fig. 7 presents the lengths of plugs measured experimentally in this study for 0.08 > Ca > 0.01 and the fitting curve using the correlation of [START_REF] Steegmans | Characterization of emulsification at flat microchannel y junctions[END_REF] (Eq. ( 7)). It can be observed that the predicting results (˛= 0, ˇ= 0.32) and the experimental data are generally in good agreement (MAPE = 11.7%), which coincides with the constants proposed by the Steegmans model (˛= 0, ˇ= 0.334). On the contrary, the datapoints for Ca < 0.01 cannot be described by the shear-driven model of [START_REF] Steegmans | Characterization of emulsification at flat microchannel y junctions[END_REF]. Based on the above discussion, one may conclude that the droplet lengths measured in this study for Ca < 0.01 can only be predicted by Garstecki model indicating the interfacial tension force dominant squeezing regime. For 0.08 > Ca > 0.01, the droplets lengths can be well predicted by both Garstecki model and Steegmans model, implying the balance between interfacial tension and viscous shear stress for the droplet break-up. Finally, the values of coefficient of variation P w-s (Eq. ( 1)) for emulsions produced in the single Y-junction minichannel replica are all smaller than 2.51% under both pure oil and 5% Butanol addition conditions. Nevertheless, bigger values of P w-s (up to 9.16%) are observed with the addition of SPAN83 or 13% Butanol, indicating that emulsions become less homogeneous when the interfacial tension w/o decreases. These results are listed in Table 5 and will be further discussed as reference case in comparison with emulsion coefficient of variation obtained in the parallel minichannel contactor.

Water-oil emulsification in the parallel minichannel contactor

Experiments on water-oil emulsification were also performed in the parallel minichannel contactor (Fig. 2) having 16 parallel identical Y-junctions. It is shown working in the video provided as supplementary material. Especially, the influences of two-phase flow properties, such as the interfacial tension o/w , the water and oil flow rates Q w and Q o , and the oil viscosity o , on the coefficient of variation of emulsion P w were investigated. Fig. 8 shows three images as examples of high, intermediate and low coefficient of variation P w measured downstream of the 16 parallel Y-junctions. Since they are in a much larger channel section, the droplets appear as nearly spherical but distorted by shear stress. This is due to the high viscosity of the oil phase and the flow convergence at the Tjunction. Their visible black tip is actually a bright spot from the halogen light source, as the images have been inverted from their original format that can be seen, for instance in Fig. 12. The obtained experimental results are reported in Table . 5.

The capillary number Ca of the continuous phase is calculated at the outcome of an individual Y-junction (cross-section 0.6 × 1 mm, see Fig. 3a). Thus, it describes the two-phase flow patterns precisely where emulsification happens, as it does for the Y-junction replica in the previous section.

Influence of parallelization on the coefficient of variation: water-pure oil emulsification

For the case of water-pure oil emulsion, the values of P w are between 3.0% and 5.8% under different water/oil flow-rate ratios (0.1 < Q w /Q o < 0.6). This is compared to the P w-s values (rightmost column in Table 5) between 1.3%-2.51% for the single Y-junction replica under the same working conditions. The noticeable augmentation in P w values is mainly due to two reasons. Firstly, the multi-scale tree-like distributor, although strictly symmetrical in geometry, could have intrinsic flow distribution non-uniformity because of the existence of elbows and the inertial forces of fluids [START_REF] Fan | Flow distribution property of the constructal distributor and heat transfer intensification in a mini heat exchanger[END_REF][START_REF] Guo | Multi-channel heat exchanger-reactor using arborescent distributors: a characterization study of fluid distribution, heat exchange performance and exothermic reaction[END_REF]. Regarding two-phase flow, the effect of distribution nonuniformity on the coefficient of variation of emulsion could be multiplied due to the proportion variation of fluids among the multiple channels. Secondly, each individual channel may not have the same surface roughness due to the limitation of fabrication accuracy during the manufacturing process. A small difference in minichannel geometric dimension may result in non-negligible influence on the coefficient of variation of emulsion after parallelization.

Recall that monodisperse emulsion is considered when P w ≤ 5% [START_REF] Xu | Correlations of droplet formation in t-junction microfluidic devices: from squeezing to dripping[END_REF][START_REF] Nisisako | Microfluidic large-scale integration on a chip for mass production of monodisperse droplets and particles[END_REF][START_REF] Nisisako | High-volume production of single and compound emulsions in a microfluidic parallelization arrangement coupled with coaxial annular world-to-chip interfaces[END_REF]. Hence the parallel minichannel contactor is capable of producing homogeneous water-pure oil emulsions since the values of P w are below or near the threshold, with the throughput multiplied by a factor of 16 owing to the parallelization.

Influences of water and oil flow rates and two-phase flow properties on the coefficient of variation

The effect of interfacial tension w/o on the coefficient of variation of emulsion P w is investigated by adding SPAN 83 ( = 989 kg m -3 ) in the oil phase, while keeping the viscosity of oil constant ( o = 58 mPa s). Fig. 9 reports the coefficient of variation of emulsion as a function of water/oil flow rate ratio for three different interfacial tension cases (water-pure oil:

w/o = 27.6 mN m -1 , water-0.3% SPAN oil: w/o = 11.7 mN m -1 , water-0.4% SPAN oil: w/o = 10.4 mN m -1 ).

It can be observed that the values of P w increase with the decreasing interfacial tension. When 0.3% or 0.4% wt. SPAN83 is added in the oil phase, the value of P w steps from 3-5% to 7-10%. Although the main diameter D w is smaller as can be seen in Table 5, values of P w ≥ 5% indicates that the emulsion produced in the parallel minichannel contactor can no longer be considered as monodisperse.

Note that the water/oil flow-rate ratio has no significant influence on the coefficient of variation under constant interfacial tension condition, as shown in Fig. 9. In microfluidics, the effect of inertia can generally be ignored due to Reynolds numbers smaller than unity [START_REF] Zhao | Two-phase microfluidic flows[END_REF]. Indeed, Reynolds number is smaller than 2.5 in the present experiments.

The negative influence of Ca number on the coefficient of variation can also be noticed, out of the squeezing regime previously identified in the Y-junction replica for Ca ≤ 0.01. Shear stress plays a significant role for Ca > 0.01. It is known as sheardriven break-up [START_REF] Thorsen | Dynamic pattern formation in a vesicle-generating microfluidic device[END_REF]. As can be seen in the corner of Fig. 9 around the value of Ca 0.02, the P w value of emulsions ranges from 3% to 9%. Globally, a monodisperse emulsification is obtained in the present parallel device for a capillary number Ca significantly below 0.02 to maintain the squeezing regime within the Y-junctions.

In order to further investigate the interplay between shear stress and interfacial tension, the dimensionless Ohnesorge number Oh which relates the viscous force to the inertial and interfacial tension force as defined in Eq. ( 4) is introduced. Fig. 10 plots the coefficient of variation of emulsion as a function of Oh for various emulsions produced in the parallel minichannel contactor (Fig. 10a) and in the single Y-junction replica (Fig. 10b). Concerning the single channel, Oh cannot be calculated in the case of a plug flow because of their noncircular shape. Thus, only 5 experiments are plotted because they produce an apparent droplet flow downstream of the Yjunction. One of them is near 5% of coefficient of variation, the others are below 2.5%. In terms of monodispersity, the cost of parallelization can be seen: the gap between Fig. 10a and b is nearly doubled concerning low-viscosity, low-interfacial tension case with 13% butanol.

A frontier line between monodisperse and polydisperse emulsions may be clearly observed around Oh = 0.45 for parallel minichannels. This means that in order to produce a monodisperse emulsion, the viscous shear stress exerted by the surrounding flow must not exceed too much the inertia of the water droplets. As noted above, it must avoid the regime of shear-driven break-up [START_REF] Thorsen | Dynamic pattern formation in a vesicle-generating microfluidic device[END_REF] to achieve monodisperse emulsion. It should be noted that the critical Oh number should be identical for both single minichannel and for parallel channels under ideal conditions. The departure comes from the small differences/errors accumulated during the parallelization process, as already explained above.

Influence of the water flow-rate on the number of water droplets in emulsions

We also counted the number n of water droplets presented in the principal branch of the collector (red area in the photo on Fig. 12) so as to quantify the concentration of water droplets in the emulsion. Results are listed in Table 5 and depicted in Fig. 11. It may be observed that the number n rises with increasing water flow-rate Q w . For example in water-5 wt.% Butanol oil emulsions, when Q w increases from 0.4, 0.7, 1.0 to 1.35 mL min -1 with almost constant oil flow-rate Q o , the number n augments from 22, 36, 45 to 59. As a result, the water flow-rate is the primary factor influencing the number of water droplets in the produced emulsions.

It is also found that the number of droplets in two types of emulsions: water-pure oil and water-5 wt.% Butanol oil is approximately equal at low water flow-rate (Q w = 0.4 mL min -1 ). However, larger difference may be observed when Q w increases, as shown in Fig. 11 the oil help to produce monodisperse emulsions with a high concentration of water droplets.

Energy dissipation rate

The energy dissipation rate E is defined as the pumping power consumption per kg of emulsion produced (Eq. ( 8)).

E = Q o p o + Q w p w o Q o + w Q w -→ [J kg -1 ] (8)
For all the experiments performed in this study, the maximum value is 43.7 J kg -1 , with the maximum throughput. Similar estimations for E are also performed using reported data in the literature, as shown in Table 6. It can be seen that the proposed device is capable of producing a quite high throughput of monodisperse emulsion at a relatively low energy dissipation rate, owing to the parallelization of mini-channels. However, the cited studies in Table 6 obtained droplets having 1 or 2 orders of magnitude finer sizes, implying an inevitable energy cost. Furthermore, different applications could be aimed: [START_REF] Belkadi | Water-in-oil emulsification in a microfluidic impinging flow at high capillary numbers[END_REF] delivered high flow rates of emulsion (25 m) for fuel applications, whereas the present device would better suit for producing solid spherical particles with a diameter about 1 mm.

4.5.

Fitting a predictive correlation [START_REF] Kiss | Formation of O / W emulsions by static mixers for pharmaceutical applications[END_REF] uses small dimensionless groups to build dimensional correlation to predict the droplet sizes. Dimensional analysis [START_REF] Zlokarnik | Scale-up in chemical engineering[END_REF] is applied to obtain a predictive formula of Oh pre fitting with the measured one.

In the field of liquid-liquid emulsification, the viscosity ratio o w is a crucial parameter as expressed by [START_REF] Grace | Dispersion phenomena in high viscosity immiscible fluid systems and application of static mixers as dispersion devices in such systems[END_REF]. Then, the capillary number Ca [START_REF] Zhao | Two-phase microfluidic flows[END_REF] and the water-to-oil flow rates ratio Qw Qo are also considered as great influencing factors. Finally, Eq. ( 9) is obtained by minimizing its gap with the present experimental results.

To evaluate the intensity of emulsification process which is affected by the liquid physical properties, the Ohnesorge number Oh compares viscous shear stress to inertia of the produced water droplets. Eq. ( 9) results from an optimization following the Nelder-Mead simplex method of [START_REF] Lagarias | Convergence properties of the Nelder-mead simplex method in low dimensions[END_REF] giving the 1.64 factor and the exponents: The results show that the capillary number Ca has a dominant influence on the emulsification process, with the highest exponent (0.38). As explained above, the experiments with addition of Butanol and SPAN are in the "shear-driven breakup" [START_REF] Thorsen | Dynamic pattern formation in a vesicle-generating microfluidic device[END_REF] regime within the Y-junctions. Thus, shear stress plays a significant role in the emulsification.

The exponent of the viscosity ratio ( o w ) and the flow rates ratio ( Qw Qo ) are close, 0.11 and 0.12, respectively. As the value of the viscosity ratio is approximate 20 times larger than the value of the flow rate ratio (0.2 0.6), the viscosity ratio is very influent to the emulsification. In this extent, [START_REF] Grace | Dispersion phenomena in high viscosity immiscible fluid systems and application of static mixers as dispersion devices in such systems[END_REF] explained that the criterion for liquid-liquid break-up is governed by the viscosity ratio. The flow rate ratio has less influence. This is interpreted as the low influence of inertia generally stated in microfluidics [START_REF] Zhao | Two-phase microfluidic flows[END_REF]. The predicted value (Oh pre ) provides a good fit with experimental values Oh expts with a MAPE of 11%, as shown in Fig. 12.

Conclusion

Water-in-oil emulsions are produced in a single Y-junction replica and in a parallel minichannel contactor. The flow regimes and the influences of the two-phase flow physical properties on the coefficient of variation of emulsion have been experimentally investigated. Based on the results and discussions presented above, the following conclusions can be drawn:

(1) For water-in-oil flows through the single Y-junction replica, an apparent droplet flow and a characterized plug flow are observed. The length of droplets/plugs can be well predicted by the correlations of [START_REF] Garstecki | Formation of droplets and bubbles in a microfluidic t-junction-scaling and mechanism of break-up[END_REF][START_REF] Steegmans | Characterization of emulsification at flat microchannel y junctions[END_REF] in the case of Ca > 0.01, and by both correlations of [START_REF] Garstecki | Formation of droplets and bubbles in a microfluidic t-junction-scaling and mechanism of break-up[END_REF] and [START_REF] Steegmans | Characterization of emulsification at flat microchannel y junctions[END_REF] in the case of 0.08 > Ca > 0.01. (2) For water-pure oil flow, both the single Y-junction replica and the parallel minichannel contactor are capable of producing monodisperse emulsion (P w-s < 2.51%, P w < 5.8%) under our tested conditions. (3) The liquid-liquid break-up regimes and the resulting monodispersity are governed by the balance between the interfacial tension and the viscous shear stress in the continuous phase of oil. Concerning coefficient of variation P w , the Ohnesorge number Oh which relates the viscous force to the inertial and interfacial tension force is a proper indicator. For our cases tested in this study, monodisperse emulsions could be produced in the parallel minichannel contactor when Oh < 0.45. (4) Compared to other microfluidic emulsification devices reported in the literature, the novel parallel minichannel contactor integrated with multi-scale tree-like distributors and collector has a comparable mono-disperse emulsification performance (droplet size about 1 mm) at a relatively low energy dissipation rate (47 J kg -1 ), with a relatively high throughput (1.5-6 mL min -1 ). (5) A predictive correlation for the Ohnesorge number has been derived for the parallel minichannel contactor. Its difference with the experimental data is 11%.

Future investigations should enable a finer emulsification, to obtain a smaller droplet diameter. The liquid-liquid breakup has to be further investigated at the Y-junction, to minimize the pressure drop and further optimize the emulsification.

  symbols ˛r assumed risk (%) ı dimensionless parameter (-) dimensionless parameter (-) superficial tension (mN m -1 ) dynamic viscosity (Pa s) density (kg m -3 )

Fig. 1 -

 1 Fig. 1 -Schematic of the experimental apparatus for the measurement of two-phase flow pattern in a single minichannel and in a parallel minichannel contactor.

Fig. 2 -

 2 Fig. 2 -Parallel minichannel contactor for water-in-oil emulsification. (a) side view; (b) photography (front viewoil and water inlets are at the top) next to a 1 D coin.

Fig. 3 -

 3 Fig. 3 -Single Y-junction replica. (a) Schematic view; (b) Photography next to a 1 D coin.

Fig. 4 -

 4 Fig. 4 -Photographies of the oil-water two-phase flow in the single Y-junction minichannel replica.

Fig. 5 -

 5 Fig. 5 -Emulsion flow patterns in the single Y-type micro-channel replica, as a function of oil and water superficial velocities. Photographies of Expts. 11, 25, 31 and 33 (green spots) are shown in Fig. 4.

Fig. 6 -

 6 Fig. 6 -Experimental data and fitting curves for the length of plugs. Current study: Y-junction (60 • ), 1000 m × 600 m; Garstecki et al., 2006: T-junction, 100 m × 33 m; Xu et al., 2008: T-junction, 200 m × 150 m.

Fig. 7 -

 7 Fig. 7 -Experimental data and fitting curve for the length of plugs for Ca > 0.01.

Fig. 8 -

 8 Fig. 8 -Three example images of the collected emulsion downstream of the 16 parallel Y-junctions in the minichannel device: (Top) High P w = 9.99% with 0.4% SPAN, (Middle) Intermediate P w = 5.81% with pure oil, (Bottom) Low P w = 3.04% with pure oil.

Fig. 9 -

 9 Fig. 9 -coefficient of variation P w as a function of water/oil flow rate ratio Q w /Q o for various emulsions produced in the parallel minichannel contactor -Corner figure: coefficient of variation P w as a function of capillary number Ca.

  Fig. 10 -coefficient of variation P w as a function of Ohnesorge number Oh for various emulsions produced in the parallel minichannel contactor (a) and in the single Y-junction replica (b).

Fig. 11 -

 11 Fig. 11 -Number of droplets n in the downstream area of the collector, as a function of the water flow rate Q w for different emulsions.

  Fig. 12 -Predicted Ohnesorge number Oh pre by Eq. (9) compared with experimental Oh expts data in the parallel minichannel contactor -Corner photography: The red area of diameter measurement downstream the multi-scale emulsion collector.

  Oh pre = 1.64 Ca0.38

Table 1 -Selected literature on droplet formation and emulsion coefficient of variation within single or parallel micro-channels.

 1

	Authors	Liquid-liquid	Experimental devices	Key findings
	Tice et al. (2004)	Non-viscous	Single T-junction
			50 m × 50 m

Table 2 -Dimensions of the multi-scale tree-like fluid distributor / collector. Bifurcating Channel length Channel width Channel depth generation

 2

		(mm)	(mm)	(mm)
	1	16	2.1	1
	2	8	1.35	1
	3	4	1	1
	4	2	0.6	1

  Table 3 lists the physical properties of working fluids measured between 22 • C and 25 • C.

Table 3 -Physical properties of the used fluids.

 3

	Liquids -temperature	Water	Pure oil		Butanol (wt.%)	SPAN (wt.%)
				5	13	0.3	0.4
	[mN/m] -25 • C	73.5	33.7	-	-	-	-
	w/o [mN/m] -25 • C	-	27.6	12.8	5.4	11.3	10.4
	[mPa s] -22 • C	0.96	58	44	34	58	58
	[g/L] -25 • C	998	865	862	857	865	865
	Supplier	Tap water	Carrefour		Sigma-Aldrich	Fisher Scientific

Table 4 -Range and precision of experimental instruments.

 4

	Parameter	Device	Range	Unit	Precision
	Water flow rate	Syringe pump	[0,2]	ml/min	1%
	Oil flow rate	Sartorius-MSE 2203	[0,1000]	g/s	5%
	Working temperature	K-type thermocouple	[-75,250]	• C	0.1 • C
	Pressure transmitters	Gems 3100 Series	[0,2.5]	bar	0.25%
	Viscosity	LAMY RM-200 rheometer	-	Pa S	1%
	Surface tension	KRÜSS tensiometer	-	mN/m	5%
	Camera	LAVision HighSpeed Star 6	-	-	-
	P w = w Dw	Image-Pro Plus 7.0	-	-	-

Table 5 -Results of emulsification in the parallel minichannel contactor. n is the number of droplets in the red appointed area in Fig. 12. As an indicator, P w-s is the coefficient of variation of emulsion measured in the individual Y-junction replica under the same working conditions.

 5

	Q w	Q o	P water	P o	n	Ca	Oh	D w	P w	P w-s
	mL. min -1		bar	bar				mm	%	%
					Pure oil, o = 58 mPa. s, w/o = 27.6 mN. m -1
	1	2.05	1.24	1.19	35	0.0074	0.33	1.07	4.48	1.36
	0.7	2.62	1.15	1.17	28	0.0096	0.36	0.97	4.42	1.35
	0.7	2.9	1.2	1.2	30	0.011	0.36	0.99	3.04	1.25
	1.35	2.56	1.3	1.2	43	0.0093	0.34	1.04	4.53	2.14
	1	2.15	1.24	1.2	37	0.0079	0.34	1.05	5.39	2.51
	0.7	2.47	1.2	1.19	32	0.009	0.37	0.89	4.32	1.53
	0.4	1.1	1.08	1.03	23	0.004	0.33	1.09	4.59	2.18
	1.35	2.88	1.24	1.28	42	0.011	0.37	0.89	4.29	1.51
	0.4	2.84	1.16	1.21	20	0.01	0.38	0.86	5.81	2
				0.3% SPAN, o = 58 mPa. s, w/o = 11.7 mN. m -1
	1	2.84	1.2	1.25	40	0.024	0.59	0.844	8.02	2.76
	1	2.29	1.21	1.19	42	0.02	0.58	0.86	7.65	2.31
	1	3.35	1.23	1.27	40	0.03	0.58	0.864	8.4	3.1
	1	4.02	1.24	1.33	41	0.034	0.59	0.83	8.9	2.56
	0.7	2.18	1.2	1.15	36	0.019	0.54	0.96	8.65	3.1
				0.4% SPAN, o = 58 mPa. s, w/o = 10.4 mN. m -1
	1	3.45	1.21	1.29	40	0.033	0.62	0.83	8.52	3.92
	1	5.04	1.32	1.38	48	0.049	0.66	0.75	9.99	5.88
	1	2.35	1.18	1.2	42	0.023	0.62	0.85	7.94	4.2
				5% Butanol, o = 44 mPa. s, w/o = 12.8 mN. m -1
	0.4	2.95	1.15	1.12	20	0.018	0.4	0.91	3.41	1.9
	0.7	2.76	1.17	1.15	36	0.016	0.38	1.05	4.95	1.57
	0.7	3.1	1.17	1.14	35	0.019	0.39	1.01	3.81	1.2
	1	2.91	1.2	1.18	47	0.017	0.38	1.04	4.25	1.32
	1	2.3	1.17	1.14	45	0.014	0.36	1.16	4.56	1.24
	1	4.28	1.18	1.3	60	0.025	0.39	0.96	4.49	2.26
	1.35	3.41	1.21	1.28	75	0.021	0.38	1.04	6.44	2.12
	0.7	4.44	1.17	1.29	35	0.026	0.41	0.87	3.05	1
	0.4	3.74	1.18	1.24	22	0.022	0.44	0.81	5	2
	1	2.68	1.18	1.15	47	0.016	0.36	1.14	3.96	1.05
	0.4	2.7	1.12	1.13	22	0.016	0.38	1.07	4.4	2.02
	1.35	3.03	1.23	1.25	49	0.019	0.38	1.06	5.54	1.95
				13% Butanol, o = 34 mPa. s, w/o = 5.4 mN. m -1
	1	3.26	1.13	1.1	53	0.036	0.51	0.81	11.98	9.16
	1	2.94	1.1	1.1	50	0.032	0.51	0.817	11.4	5.27

Table 6 -Comparison of the energy dissipation rates E (Eq. (8)) for different liquid-liquid microfluidic studies.

 6

	Authors	Flow rate of emulsion E (Eq. (8))	D w
		mL min -1	J kg -1
	Sugiura et al. (2001)	<10 -2	3.5	20 m
	The present study	1.5 -6.0	≤35	1 mm
	Nisisako et al. (2002)	2.5 × 10 -2	170	200 m
	Thorsen et al. (2001)	8 × 10 -4	176	25 m
	Belkadi et al. (2015)	80	681	25 m

Acknowledgments

The authors gratefully acknowledge the support of China Scholarship Council (CSC).

Appendix A. Supplementary Data

Supplementary data associated with this article can be found, in the online version, at http://dx.doi.org/10.1016/j.cherd. 2017.03.010.