Marieke Huisman
email: marieke.huisman@inria.fr

Cl L Ément Hurlin
email: clement.hurlin@inria.fr

PERMISSION SPECIFICATIONS FOR COMMON MULTITHREADED PROGRAMMING PATTERNS

Multithreading is the next challenge for program verification. To support modular verification of multithreaded programs, one should know when data might be accessed or updated by the different threads in the system. We propose a permission-based annotation system that is designed to do exactly this, i.e. it specifies when a thread can read or write a variable. The annotation system ensures that threads have exclusive access to a variable whenever they have the possibility to write it, thus avoiding data races. Moreover, the annotation system allows to change permissions dynamically throughout the execution. The information from the permission annotations can be used for further verification of the program. This paper shows how the annotation system can be used to specify variable access in several typical multithreaded programming patterns.

INTRODUCTION

Because of increasing requirements on performance and reactivity, multithreading nowadays is unavoidable for programmers. However, multithreaded programs are notoriously difficult to write correctly. Therefore, it is important to have the means to establish a formal correctness statement for the application (it behaves as specified, does not contain security bugs etc.). However, just as writing multithreaded applications is more difficult than writing sequential ones, verifying multithreaded applications is also more complicated than verifying sequential ones.

A promising approach to verify multithreaded programs is to abstract the behaviour of threads by specifying what they cannot do. This can be considered as the worst-case behaviour of a thread. Then, to verify a particular thread, one can reduce all other threads to their abstractions, which makes it unnecessary to know the precise state each thread is in. Our annotation system for permissions does exactly this: it specifies whether a thread can potentially change or read a variable, without specifying what it reads or writes. In addition, our system also ensures that (i) if a thread can write to a location no other thread may write or read to this location simultaneously, and (ii) if a thread can read to a location no other thread can write to this location simultaneously. The annotation system has been designed to be as general as possible, and to handle a large class of programming patterns. It can be used to show the correctness of lock-free algorithms, algorithms based on locking, and mixtures of both.

The annotation system declares for each object one or more permission states, listing all permissions on an object that may exist at a particular point in time. If a thread can be shown to hold all permissions on a particular object, i.e. it has the exclusive access to the object, then the object's permission state can be changed. Method's preconditions specify which permissions are needed to call a method, while postconditions specify which permissions are given back to the caller. As mentioned above, the annotations should ensure that at each point in time (i) there exists at most one permission to write a variable, and (ii) if a thread has the permission to write a variable, then no other thread has the permission to read this variable.

Our annotation system is strongly based on Boyland's fractional permission system [START_REF] Boyland | Checking interference with fractional permissions[END_REF], but whereas Boyland's system is defined for a simple language with parallel composition, we extended it to a Java-like language with synchronisation, thread creation and joins 1 . Moreover, our annotation system supports proper encapsulation of the object state by introducing names for the different permissions.

This paper gives a quick overview of the syntax of our annotation system, described as an extension of JML [START_REF] Leavens | JML Reference Manual[END_REF]. The remainder of the paper shows how the annotations can be used to specify typical multithreaded programming patterns. In particular, we discuss the following patterns: Immutability, Lock Splitting, Confinement across methods, Worker threads, Fork/Join algorithms, and concurrent Subject/Observers. These patterns and their implementations are all taken or inspired by Lea's book on concurrency in Java [START_REF] Lea | Concurrent Programming in Java: Design Principles and Patterns (Second Edition)[END_REF] or the Design Pattern book [START_REF] Gamma | Design Patterns: Elements of Reusable Object-Oriented Software[END_REF].

SHORT OVERVIEW OF THE ANNOTATIONS

A class can contain several permission and lock clauses. Figure 1 gives their grammar (where we use regular expression notation: X ? for 0 or 1 X, X * for 0 or more X, and X + for 1 or more X). A permission consists of a list of names and for each name, a (possibly empty) list of permission declarations. A permission declaration can specify a permission for a field of primitive type, a field of reference type, or it can be a special lock permission. A permission declaration for a primitive field contains the field name and a base permission. A permission declaration for a reference field, in addition, specifies the encapsulated permissions of the object referenced to. Finally, references can be declared to be special lock permissions. A base permission is of the form \split(n) (where n is a natural number). Its intuitive meaning is that there are at most 2 n threads accessing the field at the same time. Permissions can be split, using the equivalence \split(n) ≡ \split(n+1)*\split(n+1), where the * operator comes from separation logic [START_REF] Reynolds | Separation logic: A logic for shared mutable data structures[END_REF]. For standard JML expressions, * and & are equivalent, but their meaning is slightly different for permission formulas. Intuitively, if P and Q are permission formulas, "P * Q is valid w.r.t. to a heap h" means that h can be split in two disjoint heaps h 1 and h 2 such that "P is valid w.r.t. h 1 " and "Q is valid w.r.t. h 2 ". For this paper, it is sufficient to understand * as &. The equivalence doubles the number of threads that can access the field. Notice that permission splitting is unbounded.

For convenience we introduce W to abbreviate \split(0): since 2 0 = 1, this means that the thread has exclusive access to the field, and thus can be allowed to write to it. Any permission \split(n), where n > 0, allows only to read, because there might be other threads accessing the field. We use R to abbreviate \split(1), a basic read permission.

Encapsulated permissions are of the form \split(p,n), meaning that the reference permission contains (1/2 n ) th of the permission named p, from the object that the reference points to. We use the permission name p to abbreviate \split(p,0), i.e. exclusive hold of the permission. As base permissions, encapsulated permissions can be split using the equivalence \split(p,n) ≡ \split(p,n+1)*\split(p,n+1).

For each lock permission declaration, the class contains a special lock clause. This specifies which permissions are obtained when the lock is acquired. Intuitively, at every point in the execution, a thread has certain permissions on an object. If a thread acquires a lock, it obtains the permissions that are held by that lock, and when it releases the lock, it returns these permissions to the lock. Fields that are declared as locks are required to be final.

To ensure that the annotations avoid data races, each permission should contain at most one write permission per field, and if a permission contains a write permission, it cannot contain read permissions on the same field. Notice that this only has to be guaranteed at the current class level: it is implicitly guaranteed to hold for all objects that are referenced to by fields in the class.

Method pre-and postconditions can be extended with permission predicates, see Figure 2 for their grammar (where a path denote 0 or more indirections beginning with this or a variable name, i.e., a sequence of the form this.f1.....fn or v.f1.....fn). Intuitively, a precondition specifies which permissions are necessary to execute a method, a postcondition specifies which permissions are returned to the caller2 .

class Fraction{ //@ permission rd = num : R, den : R; protected final long num, den;

//@ ensures rd; public Fraction(long n, long d) { // ... normalize num = ...; den = ... } //@ requires \part(rd) * \part(f.rd); //@ ensures \part(rd) * \part(f.rd) * \result.rd; public Fraction plus(Fraction f) { return new Fraction(num * f.den + f.num * den, den * f.den);}}

Figure 3: Fragment of immutable class Fraction

Permission predicates can be base permissions on the fields of the object, or on the encapsulated fields of the object (named by a qualified expression, containing the name of the surrounding permission). Additionally, permission predicates can express whether a thread owns a lock or not: permission \has(path,i,j) means that the considered thread has acquired the lock on the object pointed to by path j times, and in addition at most i other threads have permission to acquire this lock. Further, we introduce to useful abbreviations: \part(p) and \hasnot(path), defined in terms of JML/permission expressions. Expression \part(p) abbreviates (\exists int n. \split(p,n) & n >= 0), i.e. it denotes a read permission on the fields contained in p. Expression \hasnot(path) indicates that a thread does not hold the lock pointed to by path: it abbreviates (\exists int i.\has(path,i,0) & i >= 0). Other useful abbreviations could be imagined, e.g. to express that exactly the same permission is returned. More experience with writing permission annotations will show which abbreviations are useful.

To allow modular verification, standard JML uses modifies clauses. However, because permissions can be considered as modifies clauses (i.e., a write permission on a field in a precondition is similar to a modifies clause mentioning this field), we do not need to specify modifies clause in this paper.

ANNOTATIONS OF DIFFERENT CONCURRENCY PATTERNS

This section shows how the annotation system presented above is particularly suited to specify several common concurrent programming patterns. The patterns that we present here are derived from Lea's pattern collection [START_REF] Lea | Concurrent Programming in Java: Design Principles and Patterns (Second Edition)[END_REF] and the Design Pattern book [START_REF] Gamma | Design Patterns: Elements of Reusable Object-Oriented Software[END_REF].

3.1. Immutability. The first example is the immutability pattern. An object is said to be immutable if after initialisation it can never change its internal state (see also [START_REF] Haack | Immutable objects for a Java-like language[END_REF]). This is useful for multithreaded programming, because accesses to (initialised) immutable objects do not have to be protected by locks. Typical applications of immutable objects are abstract data types, value containers and shared state representations. Figure 3 presents an example of an immutable object, representing the abstract data type of fractions (taken from [6, §2.1.1]). The annotations specify that the class has only one permission rd, in which both numerator (num) and denominator (den) can only be read. The constructor returns the full permission rd to its caller. These permissions can then freely be split, to give read access to the fraction object to any thread that requires this. To be able to compute with fractions (as does method plus), one needs to have some part of the rd permission. After plus has finished, a (possibly different) part of the rd permission is returned to the caller. In addition, the caller obtains the complete rd permission on the newly created object. Notice that our permission system is also suitable to express partial immutability, where only some fields are immutable, or an object only is immutable during a part of the execution (e.g., before or after calling a certain method). To support the first case, a class can declare different permissions: one containing read permission declarations for all immutable fields, and one (or more) containing write permission declarations for the mutable fields. To support the second case, different permissions can be specified, so that the object can change its state from one permission to another. Methods that need the object to be immutable, require that the object is in the appropriate permission state, and thus that the appropriate read permissions are held.

Lock Splitting.

To increase performance, it often is useful to associate locks only with certain functionalities of a class. Thus, different (groups of) fields are protected by different locks. Naturally, this makes it even more important to clearly specify which fields are protected by which lock (and to ensure that the application respects this).

Figure 4 gives a typical example of this lock splitting pattern [6, §2.4.2]. The class Person contains fields age and salary. It has a single permission d, which declares two lock permissions. All methods in the class require and ensure a fraction of this permission. Further the class contains two lock clauses: the first specifies that if a thread acquires lock l1, it obtains the permission to write the field age; similarly acquiring lock l2 gives the permission to write salary. Notice that the method specifications only mention the permission d. However, the lock clauses ensure that any access to age or salary is protected by the appropriate locks. Note that the synchronized(..) statements in class Person can be reentrant acquirements of l1 or l2. A lock acquirement is said to be reentrant if the thread acquiring the lock already owns it. In this case, Java's semantics is such that the thread continues normally. To allow reentrant acquirement of locks, our system proceeds as follows: whenever a thread acquires a lock, the thread must show it does not already have the lock (i.e., the thread can gain the permissions inside the lock) or the thread must show it already has the permissions inside the lock. Thus, after acquiring a lock, one is sure that the permissions inside the lock have been transfered to (or were already owned by) the thread considered. For example, after executing the synchronized(l1) statement in getAge, a thread gains the permissions inside l1 (i.e., write access to the field age) and can execute return age;.

Confinement Across Methods.

Thread locality, i.e. the case where an object is only accessible via a single thread, is another means to guarantee non-interference of other threads. A generalisation of this is the case where there is at most one thread at the time having a reference to an object. After passing an object to another method (which might cause it to be accessed by different threads), the calling method guarantees that it does not access it anymore. Thus, possible changes to the object do not change the correct behaviour of the method. This pattern is called confinement across methods [6, §2.3.1].

Figure 5 shows a typical example of this pattern. The basic Point class contains a single permission set that allows to write both its fields. In class Plotter, the method showNextPoint creates a new point class, and properly initialises it. After initialisation, it gives the newly created point to the display method (that will typically use special threads for doing the graphics). The fact that the object is given away, is made explicit by the method's postcondition: the caller of showNextPoint does not get back any permissions on the newly created point. Also the display method requires the full permission on the point, and does not return any permission.

Notice that other confinement patterns (confinement within a thread, confinement within an object) can also be expressed in a natural way with our annotation system.

Worker Threads.

The next pattern that we discuss is the worker thread pattern; an important pattern for many industrial applications. A main thread prepares different tasks, and sends them off to a worker thread for computation. This ensures that the main thread is never blocked for a long time, and thus that the application can stay reactive. It is important that access to the task class MainThread extends Thread{ //@ permission before = wk : R (p), t : R (wr, regain); //@ permission after = wk : R (p), t : R (wr); //@ permission wr = ...; //@ permission regain; //@ requires wr; //@ ensures wr; public void prepare(){} //@ requires wr; //@ ensures wr; public void doTheJob(){ ... } //@ requires wr; public void setIsDone(){ ... } //@ requires regain; //@ ensures \result ==> wr && !\result ==> regain; public boolean isDone(){ ... }} Figure 6 shows an implementation of a main thread, that prepares jobs for a worker thread. In our simple example, it has read-only references to a single worker thread and a single task object. It has two different permissions: before and after. The only difference between these two permissions is whether the server thread has a special regain permission on task t. The role of this regain permission will be explained below. The permissions on the worker thread and task objects are encapsulated in the permission state for the ServerThread object. Notice that since both before and after contain the full encapsulated permissions p and wr (explained below), the main thread cannot hold both the before and the after permission at the same time. Notice further that ServerThread's constructor requires that its parameter has an appropriate permission on the worker thread. In contrast, since the task is a newly created object, no further requirements are necessary.

final
Figure 7 shows the implementation of a worker thread (inspired by [6, §4.1.4]). The worker thread has a single permission, which ensures that access to the task pool is protected by a lock. This avoids possible race conditions.

Figure 8 gives a possible implementation of tasks. It has a permission wr giving write access to the fields of the task (not further detailed here), and a special regain permission, mentioned above. This permission is only a marker, it does not give access to any of the fields. The holder of the regain permission is allowed to regain the wr permission, once the task is done. The constructor of Task returns this permission to the thread that constructed the object. This (or some other) thread can inspect whether the task is finished, and if so, it regains the write permission wr on the fields of the task, and the regain permission is destroyed.

The run method of the ServerThread first prepares the task, and then it puts it in the task pool of the worker thread. Notice that the specification of addTask in WorkerThread class Fib extends FJTask{ int number; //@ permission p = number : W; //@ ensures p; Fib(int n) { number = n; } //@ requires p; //@ ensures p; public void run() { int n = number; if (n = 1) number = 1; else {Fib f1 = new Fib(n -1);Fib f2 = new Fib(n -2); coInvoke(f1, f2); number = f1.number + f2.number;}}} requires the wr permission on the task object, and does not return it, i.e. the server thread looses its permission to do something with the task. When the WorkerThread takes the task from the task pool (using getTask) it acquires the permission wr on the task, and thus it is able to do the job, and then set a volatile field 3 to signal that it is done. By executing the setIsDone method, it looses its permission to access the task. In the mean time, the ServerThread can continue with other things (e.g. react on other requests). When it needs the result of the task, it can inspect whether the task is finished, using the isDone method. For this it needs the special regain permission. If the isDone method returns true, the regain permission is destroyed, and instead the wr permission is returned (otherwise the regain permission is kept). Thus, implicitly the permission of ServerThread changes from before to after. The use of the special regain permission ensures that the wr permission is not duplicated: once the wr permission has been returned by isDone, the regain permission is destroyed, and thus the isDone method cannot be called anymore for this task. A variant of this specification would be that the method would only return the wr permission, if the method was called with the regain permission (which is then destroyed). Alternatively, one could also imagine a specification where method isDone requires only a part of the regain permission and returns the corresponding part of the wr permission (see method run of the Matrix example in the next paragraph for such an example).

Fork/Join

Algorithms. An important class of concurrent algorithms are fork/join algorithmsthe concurrent variation of divide-and-conquer algorithms -i.e. each thread spawns off several other threads to do sub-computations. It waits for all these threads to finish, and then combines their results into a single result. We show how our annotation system can be used to show that the fork-join implementation in Figure 9 of the Fibonacci function [6, §4.4.1]4 does not contain data races.

Class Fib contains a permission p, that allows to write the field number 5 . The implementation of the method coInvoke from the class FJTask is such that it implicitly behaves as f1.fork(); f2.fork(); f2.join(); f1.join();. Since the fork method starts a new thread, that will execute its run method, the precondition of the run method is propagated to also be the precondition of the fork method. Further, each runnable object is supposed to contain an implicit join permission, that is returned to the creator of the class (and can be passed around). The join method has the following (implicit) specification: where Q run is the specified postcondition of the method run (notice that to get even more precision and flexibility, we could specify that only a fraction of the join permission is required, and that exactly this fraction of the permissions in run's postcondition are returned). Thus, when a thread creates objects f1 and f2, it obtains the permission to write their number fields, plus a special join permission. Within coInvoke, the permission on number is given to the forked threads, and after joining, they regain the permission to access number -however, they have lost the implicit join permission. As a more complex example, we sketch a program where we split and recombine permissions. Suppose we have a class Matrix with permission wr, that allows to write the elements in the matrix (but not to change the shape of the matrix), and with appropriate method annotations 6 . To be able to initialise the matrices, we need write permission on the elements. However, if we suppose that we only need read access to do the matrix computations, we can have the following programming pattern (where we have matrices a and b): This specifies that a fraction of the wr permission on a and b is needed for the thread to start, and that after joining the finished thread, exactly the same permission is given back. It is crucial here that exactly the same fraction of the permission is given back: this allows to conclude that after joining all threads, the main thread holds the complete wr permission again, and thus that the matrices can be re-initialised. 5 Notice that we do not require number to be volatile, in contrast with Lea's implementation. With our annotation system, and knowledge of the new Java Memory Model [START_REF] Manson | The Java memory model[END_REF], we can show that there will be no data races, thus there is no need for this variable to be volatile. 6 Where the grammar of permission declarations is extended to array declarations in the obvious way.

class Subject{ //@ permission p = {this}; //@ r = v : W; //@ lock {this} = v : W, obs : W (\split(q,1)); Observer obs; int v; /*@ requires (\exists int i,j. \has(this,i,j) & i >= 0 & j > 0) * @ obs.s == this; @ ensures (\exists int i,j. \has(this,i,j) & i >= 0 & j > 0) * @ obs.s == this; @*/ public void notifyObs(){ obs.update(); } /*@ requires (\exists int i,j. \has(this,i,j) & i >= 0 & j > 0) * @ obs.s == this; @ ensures (\exists int i,j. \has(this,i,j) & i >= 0 & j > 0) * @ obs.s == this; @*/ public void setState(int v){ this.v = v; notifyObs();} //@ requires \part(r); //@ ensures \part(r); int getState(){ return v; }} The last pattern we discuss is the subject observer pattern; a typical example of object-oriented programming [3, §5]. In this pattern, a single subject is observed by several observers: when the subject's state changes, the observers are notified, so that they can update their internal representation of the subject's state. Figures 10 and11 show a fragment of an implementation of this pattern. For clarity of presentation, we only have one observer per subject.

The Subject class has a permission p that allows to lock the considered subject. Once a subject is locked, write access to its state (v) and permission q on the observer are granted. Class Subject imposes client-side locking: to call a method on a subject, one must lock it beforehand. For example, notifyObs's precondition indicates that to call this method, a thread has to hold the lock on the subject: \has(this,i,j) with j greater than 0 is required.

The method described in paragraph 3.2 to handle reentrance works also well in this pattern. For example, when the synchronized(this) statement is executed in method update of class Observer, our system makes sure that the observer was not already locked (in this case, the permissions inside the observer is transferred to the executing thread), or that the thread considered already has the permissions inside the lock (i.e., write access to the field cache). Thus, after acquiring the observer a thread can execute cache = v;.

Note that a thread can lock the subject and the observer if it has permission p on the subject, while another thread can lock the observer if it has a part of permission q of this observer. Two different threads can simultaneously hold permission p and a part of permission q, because the subject only has half of the permission q on its observer (see annotation lock {this} = ..., obs : W class Observer{ //@ permission q = {this}, sub : R (); //@ lock {this} = cache : W; final Subject sub; int cache; //@ ensures q * s == sub; public Observer(Subject sub){ this.sub = sub; } //@ requires \part(q) * \part(sub.r); //@ ensures \part(q) * \part(sub.r); void update(){ int v = sub.getState(); synchronized(this){ cache = v; }} //@ requires \part(q); //@ ensures \part(q); public synchronized int getCache(){ return cache; }} (\split(q,1)); in class Subject). Thus, one can write a program where one thread updates the subject and the observer, while another thread simultaneously inspects the observer (to update a GUI for example). Finally, to notify the observer, one has to know that its underlying subject (field sub) is the subject where notifyObserver is called: this is specified (using standard JML) in notifyObs's precondition by obs.s == this. Because of this equality, the permission \split(v,0) (which is contained in \has(this,i,j) at notifyObs's entry) can be used to show that the precondition of update is satisfied.

CONCLUSION & FUTURE WORK

We have shown how our annotation system for permissions can be used to specify several common multithreaded programming patterns. The specifications are intuitive to understand, and at the same time highly expressive. They allow to prove absence of race conditions, and moreover the annotations can be used as auxiliary information for the further verification of the program: if we can deduce from the permission annotations that other threads cannot interfere at certain program points, then we do not have to consider the possible interleavings at these points (see [START_REF] Huisman | The stability problem for verification of concurrent object-oriented programs[END_REF] for details).

Currently, we have implemented a run-time checker for the annotation system, and we are working on the development of a static verification method. For this, we will extend an existing translation of JML-annotated programs into BoogiePL [START_REF] Barnett | Boogie: A modular reusable verifier for object-oriented programs[END_REF] with information about the permissions of the current threads. We will then use a verification condition generator for BoogiePL to generate appropriate proof obligations. order to establish new connections. This is one of the lessons that I now try to teach to my own PhD students. Happy birthday!

Figure 1 :

 1 Figure 1: Grammar for permission declarations

  PermPred ::= \split(Path.PermName,Nat) | \part(Path.PermName) | \has(Path, Nat, Nat) | \hasnot(Path) | Path.PermName | PermPred * PermPred Path ::= ArgAccess ? | ArgAccess ? .FieldAccess ArgAccess ::= this | VarName FieldAccess ::= FieldName | FieldAccess.FieldName

Figure 2 :

 2 Figure 2: Grammar for permission predicates

Figure 4 :

 4 Figure 4: Fragment of class Person, illustrating the lock splitting pattern

Figure 5 :

 5 Figure 5: Fragment of class Plotter, showing confinement across methods

Figure 7 :

 7 Figure 7: Fragment of class WorkerThread, part of the worker thread pattern

Figure 8 :

 8 Figure 8: Fragment of class Task, part of the worker thread pattern

Figure 9 :

 9 Figure 9: Fragment of fork-join implementation of Fibonacci function

  class Matrix{ final int[][] elems; //@ permission wr: w = elems : R[R[W]]; ...

  while(true){ initialise(a,b); // split matrix permissions for different threads coInvoke(...); // do fork-join matrix computations // recombine fractioned permissions into complete permission } provided that the run method for the matrix computation is specified as follows: requires \part(a.wr) * \part(b.wr); ensures \part(a.wr) * \part(b.wr) & (\forall int n.\old(\split(a.wr, n)) = \split(a.wr, n)) & (\forall int n.\old(\split(b.wr, n)) = \split(b.wr, n));

Figure 10 :

 10 Figure 10: Fragment of class Subject

Figure 11 :

 11 Figure 11: Fragment of class Observer

Notice that such a thread creation/join mechanism can simulate parallel composition, but parallel composition cannot simulate all programs using thread creation and joining.

We do not consider exceptional postconditions in this paper. Typically, exceptional permission postconditions would coincide with the method's normal permission postconditions.

To avoid possible race conditions, and therewith unexpected behaviours.

Even though, as pointed out by Lea, this is an unrealistic example, because there are much faster non-recursive solutions. However, because of its simplicity, it nicely illustrates the working of our annotation system.

Acknowledgements. It is our pleasure to dedicate this paper to Henk Barendregt. Henk has been one of the PhD supervisors of the first author, Marieke Huisman. During this period, Henk has taught me about the importance of looking at known facts with a different and fresh mindset, in

This work is funded in part by the IST programme of the EC, under the IST-FET-2005-015905 Mobius project and the French national research organisation (ANR), under the ANR-06-SETIN-010 ParSec project. Creative Commons


Michel Belliard

Numerical modeling of an in-vessel flow limiter using an Immersed Boundary Approach

Keywords: In-vessel flow limiter, Thermal-hydraulic system safety, CATHARE, CFD, GENEPI, Immersed Boundary Method

This work is in the context of the mitigation of the consequences of a large-break loss of coolant accident in a pressurized water reactor. To minimize the flow leaving the vessel and prevent or delay the uncovering of the core, CEA has devised a device, named in-vessel flow limiter, limiting the flow of fluid from the vessel to the break. The goal is to interfere as little as possible with the nominal operation flow and maximize the fluid retained in the event of this kind of accident.

In order to quickly perform a series of 3D-CFD simulations to optimize this device, it is imperative to have a simulation tool that provides sufficiently accurate results in a reasonable time. For this goal, an immersed boundary condition approach is retained. The solid obstacles constituted by the fins of the device are not extruded from the fluid domain, but included in the calculation domain itself. Their presence is considered by a local forcing term.

Through 3D/1D up-scaling of CFD global quantities, local pressure-drop coefficients, induced by the in-vessel flow limiter, can be provided to thermalhydraulic system safety codes. It allows safety studies of the thermal-hydraulic system taking into account the in-vessel flow limiter presence in a more realistic way.

Introduction -Context

The context of this work is set in the domain of Generation II and III nuclear power plants. Generation II reactors are the class of commercial reactors that was built by the end of the 1990s and that includes several types of design: PWR, BWR, CANDU, AGR and VVER. Generation III reactors are the innovative designs that are under construction or still in design phase: EPR, ATMEA1, AP1000, APR1400, ESBWR, . . . [START_REF] Belliard | Advanced Reactors Information System (ARIS)[END_REF]. More specifically, we focused on the light-water Pressurized Water Reactors (PWRs), which are the main type of reactors built and exploited in France. Nowadays passive safety systems are more and more included in the nuclear-reactor safety strategy to mitigate design basis accidents (for example AP600 and AP1000 [START_REF] Belliard | Advanced Reactors Information System (ARIS)[END_REF]; see also [START_REF]IAEA, Iaea-tecdoc-1624: Passive safety systems and natural circulation in water cooled nuclear power plants[END_REF]). A passive safety system is a system that activates itself without the need of mechanical or electrical actuation. The passive systems are divided into four main categories (A to D), depending on the particular phenomena/device that is not used for the activation of the structure [START_REF]IAEA-TECDOC-62: Safety related terms for advanced nuclear plants[END_REF]: For instance, the fuel cladding belongs to the category A (1, 2, 3 and 4) and the pressurizer surge line or the hydraulic diode -one-way flow reduction through vortex effect -to the category B (2, 3 and 4). The accumulators belong to the category C (3 and 4) and the SCRAM to the category D (4 only).

The interest of these particular systems is given by the possibilities that derive from their employment. Some of the main benefits are: the simplification of the pipe networks for the safety injection (SI) systems, the potential disappearance of some active elements such as some specific pumps and the economical saving (less active systems to be placed and operated). At CEA, some studies on passive safety systems have been done in the past years, notably for the in-vessel flow limiter (hydraulic diode) patented by the CEA [START_REF] Gautier | Dispositif limiteur de débit inverse de fluide[END_REF] designed to limit the amount of water lost during the short-term sequence of a Large-Break (LB) Loss Of Coolant Accident (LOCA), cf. Fig. 1.

French 900 MWe CP1 and 1,300 MWe P4 PWRs and low-pressure PWRs have been the reference reactors for these investigations [START_REF] Gautier | Development of the safety systems for a simplified, low pressure, medium sized PWR[END_REF]. An other example is the advanced accumulator with passive hydraulic-diode device considered in the Generation III projects (ATMEA1, AP1000, APR1400 , . . . ) [START_REF] Shiraishi | The Advanced Accumulator: A New Passive ECCS Component of the APWR[END_REF][START_REF] Chu | Development of Passive Flow Controlling Safety Injection tank for APR1400[END_REF]. The goal of this device is to set up a two-step injection regime. The first one is a high-rate injection of the amount of water needed to fill the vessel lower plenum and down-comer. Then, the second step is a low-rate injection limited to just maintain the water level. The expected goal is a better use of the water injected by the SI accumulators and a bigger delay for the on-set of the SI pumps.

Let us notice that it is important to assess the conjoint effect of these kind of passive devices. In fact, the effect of an elementary device can be increased or minimized in conjunction with other ones. For instance, In-vessel flow limiters and advanced accumulators contribute to strongly reduce the short-term primary-mass lost during LB-LOCA as demonstrated in [START_REF] Stratta | Thermal-hydraulic study of passive safety systems based on the hydraulic diode principle for the management of large-break loss of coolant accidents[END_REF] on a generic 3-loop middle-range electrical-power reactor of 1150 MWe, taking inspiration from the ATMEA1 reactor. In reference to the case without hydraulic diodes, the accumulator injection time is more than doubled (which means that the pumps can start with a bigger delay), the reflooding level is increased of almost +35% and the peak cladding temperatures are reduced of about -10% and -43% on the short and long term. These computations were done with the French reference thermal-hydraulic system safety code CATHARE [START_REF] Barré | The CATHARE code strategy and assessment[END_REF][START_REF] Geffraye | CATHARE 2v25 2: A single version for various applications[END_REF], originally devoted to the study of water-cooled reactor transients (standard operations or accidental transients from any kind of failures or size and location of breaks), that is based on 0D/1D and 3D modules using six-equation (mass, momentum and energy) two-fluid models. But, the relevance of these system-scale studies depends on the level of realism of the data introduced to model the hydraulic diodes.

In order to take into account the large-scale effect of hydraulic diodes in safety-system studies, we need information coming from small-scale experimental or numerical experiments. For instance, results from 3D Computational Fluid Dynamic (CFD) simulations can be used to up-scale relevant characteristics as the global pressure drop induced by the hydraulic diode. With this information, it is possible to model the effect of an in-vessel flow limiter through the definition of the pressure-drop coefficient in a 1D CATHARE model. This approach is already classic in this context, cf. document [START_REF]Assessment of CFD Codes for Nuclear Reactor Safety Problems -Revision 2[END_REF] for CFD codes in the nuclear reactor safety problems.

Moreover, CFD studies can be involved in the optimization process of the flowlimiter geometry to maximize the global pressure drop during a cold-leg LB-LOCA. The value of the global pressure-drop coefficient is linked to the geometry of the fins, cf. Fig. 1. The optimization process consists to design a fin geometry minimizing (respectively, maximizing) as much as possible the pressure drop during nominal operations (respectively, a cold-leg LB-LOCA). As many geometric scales are presented at the same time (typically several meters for the downcomer radial scale and one centimeter for the fin thickness), a CFD simulation of the two-phase flow inside the flow-limiter can be time consuming. It is incompatible with a big number of simulations needed to optimize the design. Instead of this conventional approach (i.e. body-fitted approach), we are motivated by fictitious domain approach [START_REF] Saul'ev | On the solution of some boundary value problems on high performance computers by fictitious domain method[END_REF][START_REF] Rukhovets | A remark on the method of fictive domains[END_REF] allowing a less precise but faster estimation of the pressure-drop coefficient. Following this way, we consider simulations over a full computational domain including the in-vessel flow-limiter fins and re-introduce their presence adding local external forces on the immersed boundaries. Moreover, a homogeneous relaxed equilibrium model of a liquid-vapour mixture [START_REF] Clerc | Numerical Simulation of the Homogeneous Equilibrium Model for Two-phase Flows[END_REF] can be considered. For same space discretization, this kind of three-balance-equation model generally run faster than a six-balance-equation model. Once a particular geometry exhibited as a good candidate, a reduced number of body-fitted CFD computations can be done to refine the design.

The aim of this paper is to present a methodology concerning the design of a fast-running two-phase CFD model of the in-vessel flow-limiter device and the 3D/1D up-scaling of these CFD results in order to provide some useful input data in 1D thermal-hydraulic system safety code. The numerical/experimental validation of this simulation tool is not the goal of this paper and only brief elements of validation are given here. In the future, a full validation will be issued based on an improved numerical scheme including, in particular, a second-order in space interpolation scheme across the Immersed Boundary (IB) [START_REF] Introïni | A Second Order Penalized Direct Forcing for Hybrid Cartesian/Immersed Boundary Flow Simulations[END_REF] and an IB-condition consistent fractional-step method [START_REF] Belliard | Penalized direct forcing and projection schemes for navier stokes[END_REF]. This paper is structured as follow. The two-phase fluid CFD model and the first-order in space immersed boundary models are first presented in Sections 2 and 3. Then, the conditions of the in-vessel flow limiter study are precised in Section 4 with the GENEPI code as CFD tool. The process concerning the 3D/1D up-scaling of the global pressure-drop coefficient is described in 5. Re-sults are given in Section 6 providing a range of values that may be used in input of safety-system codes and allowing us to compare with the pressure-drop coefficient values used in CATHARE studies [START_REF] Stratta | Thermal-hydraulic study of passive safety systems based on the hydraulic diode principle for the management of large-break loss of coolant accidents[END_REF].

Finally, we give some words concerning the perspectives of this work in Section 7 and draws some perspectives about future work concerning passive safety systems and expected gains allowed by their introduction.

The two-phase CFD model

The considered two-phase CFD model is the GENEPI one [START_REF] Grandotto | Calculs des écoulements diphasiques dans les échangeurs par une méthode aux éléments finis[END_REF][START_REF] Grandotto | Steam Generator Two-Phase-Flow Numerical Simulation with Liquid and Gas Momentum Equations[END_REF], designed for the steam-generator two-phase flow steady-state 3D computations through the resolution of three balance equations for a water liquid/steam mixture. It is based on a homogeneous relaxed equilibrium model with thermodynamic equilibrium of the two phases. But closure laws take into account the liquid/steam momentum disequilibrium. This code incorporates the possibility to model thin no-penetration obstacles using Immersed Boundary Conditions (IBCs) [START_REF] Belliard | Fictitious domain methods for two-phase flow energy balance computations in nuclear components[END_REF].

After averaging the local and instantaneous mass, momentum and energy equations for each phase, these are merged to obtain a mixture description of the two-phase flow. Provided that the following assumptions hold: (i) surface tension, viscous and turbulent dissipation are neglected and pressure terms are neglected in the enthalpy balance equation, (ii) same pressure for steam and liquid, (iii) eddy viscosity model, one obtains for two-phase flows the following mass and momentum balance equations:

∇.G = 0, ( 1
)
ρ∂ t V + G. ∇V + ∇.(x(1 -x)ρV R ⊗ V R ) = ρg -∇P -ΛV + ∇.µ T ( ∇V + ∇T V), ( 2
)
ρ∂ t H + G.∇H + ∇.(x(1 -x)ρLV R ) = ∇.(χ T ∇H), (3)
with µ T the turbulent dynamic viscosity, χ T the turbulent diffusion coefficient, Λ the singular-obstacle tensor and V R the relative velocity given by the driftflux Lellouche-Zolotar model [START_REF] Lellouche | Mechanistic Model For Predicting Two-Phase Void Fraction For Water in Vertical Tubes, Channels, and Rod Bundles[END_REF] and based on the Zuber-Findlay approach [START_REF] Zuber | Average Volumetric Concentration in Two-Phase Flow Systems[END_REF].

The density ρ, the static quality x and the latent heat L are determined through the equation of state of the water as a function of the pressure and mixture specific enthalpy. We solve in H, P and V variables. The time term presence in Eq. ( 2) allows to search the steady-state regime through a transient computation of a thermally dilatable fluid, ∇.G = 0, cf. Eq. [START_REF] Belliard | Advanced Reactors Information System (ARIS)[END_REF]. For that, the Chorin- Petrov-Galerkin (SUPG) method is applied to correct the convective term [START_REF] Brooks | Streamline upwind/petrov-galerkin formulations for convection dominated flows with particular emphasis on the incompressible navier-stokes equations[END_REF].

A conjugated gradient method, preconditioned by the diagonal, is used to solve the arising linear systems.

In GENEPI, the turbulence modeling is done thank to a very simple local scalar model, the Schlichting model [START_REF] Schlichting | Boundary Layer Theory[END_REF]:

µ T = a S |G|L T (4)
where L T is a characteristic length and a S a coefficient. This turbulence model is known to be quite diffusive. The turbulent diffusion coefficient χ T for the enthalpy balance equation is defined via the Prandtl number P r = µ T χ T .

130

According to the hyperbolic nature of the flow equations, Dirichlet boundary conditions are used at the inlets of the domain (mass flux and enthalpy) and

Neumann boundary conditions at the outlets (pressure). The other boundaries of the domain are impermeable walls. Generally, these are considered adiabatic and with no shear stress. 135

The Immersed Boundary Model

In the fictitious domain approach, introduced in the fifties by Hyman [START_REF] Hyman | Non-iterative numerical solution of boundary-value problems[END_REF] and the Russian's school [START_REF] Saul'ev | On the solution of some boundary value problems on high performance computers by fictitious domain method[END_REF][START_REF] Rukhovets | A remark on the method of fictive domains[END_REF], the original domain Ω is embedded in a fictitious domain Ω which is geometrically bigger and generally simpler-shaped.

Doing this, some immersed boundary Σ appears such that Ω

= Ω ∪ Σ ∪ Ω e ,
where Ω e is the complementary or 'exterior' domain (as the fins). The spatial discretization is now performed in Ω, independently of the shape of the original domain Ω. Then, the resolution of the new problem in Ω will be faster and simpler. The main issue is to enforce the original boundary conditions on the immersed interface Σ which is non-aligned with the mesh.

The ISI method

In this work, among the numerous fictitious domain methods (see [START_REF] Belliard | Fictitious domain methods for two-phase flow energy balance computations in nuclear components[END_REF] for a short introduction), we consider an element of the set of Immersed Boundary Methods: the fictitious domain method with Immersed Spread Interface (ISI) [START_REF] Ramière | A fictitious domain approach with spread interface for elliptic problems with general boundary conditions[END_REF][START_REF] Belliard | Fictitious domain methods for two-phase flow energy balance computations in nuclear components[END_REF]. The fictitious problem to be solved in Ω is built from the original problem in Ω, but an additional term takes into account the immersed boundary conditions. For velocity Dirichlet boundary conditions, the singular-obstacle tensor Λ of Eq. ( 2) will play this role. It allows us to take into account the nopenetration condition of the flow limiter fins in an implicit way during the first step of the Chorin-Gresho method.

Let us (u, v, w) be the local basis linked to a given obstacle (i.e. fin). The vectors u and v are tangential to the obstacle and the vector w is normal to the obstacle. For the element e, we define the singular-obstacle tensor by:

Λe = ρ e A e Ω e      Λ u 0 0 0 Λ v 0 0 0 Λ w      (5)
with A e the measure of this obstacle (area, m 2 ) intercepted by the element e

and Ω e the measure (volume; m 3 ) of this element. Λ u , Λ v and Λ w are the tensor coefficients in, respectively, the directions u, v and w. Here, we consider no-penetration obstacles in the normal direction, Λ w = 1/ with 0 < << 1, and slip conditions in the tangential directions, Λ u = Λ v = 0.

Immersed interfaces Σ, as the flow-limiter fins, are modeled by a collection of linear plane surfaces intercepting elements of Ω. In each intercepted element, the measure of the intercepted surface and the external normal vector are known.

As the singular-obstacle tensor is defined by element, all the nodes belonging to this element are concerned: i.e. we have a spread interface approximation of Σ.

Convergence order and elements of validation

As for the L 2 -penalty methods [START_REF] Saul'ev | On the solution of some boundary value problems on high performance computers by fictitious domain method[END_REF][START_REF] Angot | A penalization method to take into account obstacles in incompressible viscous flows[END_REF], that the ISI method generalizes, the parameter is known as the penalty parameter. Regardless to the Navier-Stokes solving method, the theoretical rate of convergence of the penalized solution toward the body-fitted one is comprised in the range [O( 1/4 );

O( 1 )] in L 2 (Ω)
norm [START_REF] Angot | A penalization method to take into account obstacles in incompressible viscous flows[END_REF]. Let us notice that for Dirichlet boundary conditions and elliptic problems, the theoretical rate of convergence in space of the Q 1 -finite element method with non-boundary-fitted meshes is O(h 1 ) in L 2 (Ω) norm, with h the space step [START_REF] Ramière | Convergence analysis of the Q1-finite element method for elliptic problems with non-boundary-fitted meshes[END_REF].

As a whole, contributions to the validation of this IB approach can be found in [START_REF] Belliard | Fictitious domain methods for two-phase flow energy balance computations in nuclear components[END_REF] and [START_REF] Introïni | A Second Order Penalized Direct Forcing for Hybrid Cartesian/Immersed Boundary Flow Simulations[END_REF]. On one side, in the context of dilatable two-phase flow elliptic problems, the work mentioned in [START_REF] Belliard | Fictitious domain methods for two-phase flow energy balance computations in nuclear components[END_REF] validates the ISI method with respect to body-fitted finite-element computations and to the JEBC method (an IB method using a finite-volume discretization). A first-order rate of convergence in space is numerically reached.

On the other side, in the context of incompressible one-phase flow Navier-Stokes equation, the work mentioned in [START_REF] Introïni | A Second Order Penalized Direct Forcing for Hybrid Cartesian/Immersed Boundary Flow Simulations[END_REF] gives elements of validation for a finitevolume first-order penalty method very similar to the ISI method. Again, a first-order rate of convergence in space is numerically reached on the test case of a laminar flow around a static cylinder of diameter D (Reynolds number = 20). These results are in very good agreement with those proposed in the literature, cf. Table 1. In a lesser degree, it is also true for our own GENEPI results using the ISI method (about 10% on the drag coefficient and 30% on the recirculation length) giving confidence in the ability to catch the magnitude of an obstacle's drag coefficient. Lw: recirculation length. [START_REF] Introïni | A Second Order Penalized Direct Forcing for Hybrid Cartesian/Immersed Boundary Flow Simulations[END_REF] (base) refers to the first-order penalty method of [START_REF] Introïni | A Second Order Penalized Direct Forcing for Hybrid Cartesian/Immersed Boundary Flow Simulations[END_REF]. ISI: 18 cells in the diameter D and = 10 -5 .

GENEPI

CFD studies of the in-vessel flow limiter

With the CATHARE code, even it allows coarse 3D meshing, it is not possible to run simulations including the effective geometry of the fins and the turbulence that derives from it. So it is really tough to evaluate the exact value to be set as pressure-drop coefficient and only parametric studies can be envisaged. This section is devoted to a first insight into the hydraulic of the in-vessel flow limiter through CFD simulations using the GENEPI code. In this preliminary work, we mainly restrict ourselves to liquid one-phase flows. But it is not a limitation and illustrations of simulations with two-phase flows are provided below.

Computational domain and meshing

The CFD computational domain Ω is a simplified planar geometry of 4. CATHARE computations in [START_REF] Stratta | Thermal-hydraulic study of passive safety systems based on the hydraulic diode principle for the management of large-break loss of coolant accidents[END_REF] highlight the benefice brought by a out-flow pressure-drop coefficient K out equals to 11.7, corresponding to a multiplying factor K = 15 of the initial pressure-drop coefficient. The motivation of this 3D/1D up-scaling approach is to answer to the question: is this value realistic ?

Methodology

Let us recall that the CFD computational domain is based on the CATHARE space discretization. To compare the CATHARE results to the CFD ones, we do not directly compare the value of the CATHARE local pressure-drop coefficients K out to a CFD estimation. But we consider the global pressure-drop coefficients K gl computed between the in-flow/out-flow section of the CFD computational domain marked in red in Fig. 2.

Taking as reference the normal mixture mass-flux at the beginning of the cold-

leg nozzle, G cl • n cl = ρ cl V cl
• n cl , we define K gl as:

K gl = 2 < P in > -< P out > < ρ cl >< |V cl • n cl | > 2 (6)
where the symbol < • > denotes an area average and n cl the surface normal.

A positive value of K gl denotes a pressure loss: < P out > ≤ < P in >. From a discrete-space point of view, considering the three meshes M i = M 1 . . . M 3 , we denote by K M i the computed approximations of K gl . Let's emphasize that the variables in Eq. ( 6) are defined for a liquid/steam mixture.

Similarly to the CFD global pressure-drop coefficient K gl , we can define a global pressure-drop coefficient for the CATHARE computation K gl,cath :

K gl,cath = 2 P in -P out ρ cl V 2 cl . (7)
Here the density ρ cl is computed using an upstream approximation and the pressures P in and P out using a downstream approximation. Obviously, the global 260 pressure-drop coefficient K gl is linked to the local pressure-drop coefficients K in and K out of the down-comer volume element.

Let us notice that the global pressure-drop coefficient can be coarsely estimated using macroscopic analytic expressions, as Borda-Carnot law [START_REF] Massey | Mechanics of Fluids[END_REF], or experimental correlations, as Idel'cik [START_REF] Idel'cik | Mémento des pertes de charge[END_REF]. The Borda-Carnot law defines the ir-265 reversible losses for incompressible flows through head loss coefficient 0 < ξ:

E 1 = E 2 + ξ/2 ρ(V 1 -V 2 )
2 with the head defined as 0 < E = P + 1/2 ρV 2 and location 2 downstream to location 1.

Results

As the GENEPI's turbulence model is quite rough, a parametric study is performed on the coefficient a S and the turbulence characteristic length L T , cf. Eq (4): a S ∈ {a=0.047, a/10} and L T ∈ {0.3, 1.0, 2.0}. The value of the penalty parameter is set to = 10 -5 .

We consider that the GENEPI steady state is reached when the relative L 2 -norm difference of the variables between two consecutive time iterations is less than 5.10 -3 δt. Here the variables are the pressure, the mass flux and the enthalpy (if computed) and δt is the time step. Computations were ran on a 2.67 GHz-Xeon workstation with 2 Go RAM. Whatever the configuration is, it takes typically about 5,000 time steps (δt ≈ 10 -3 s) and 5 to 6 hours of CPU time to reach the steady state on the finest mesh M 3 . These CPU times are compatible with a big number of simulations needed to optimize the flow limiter geometry. Although the mesh convergence is not fully reached, the trend of the evolution of the pressure-drop coefficient versus the space-step index is globally caught. Balancing Tensor Diffusivity method [START_REF] Gresho | A modified Finite Element Method for Solving the Time-Dependent Incompressible Navier Stokes Equations (Part 1: Theory)[END_REF] used instead of the SUPG method. as well as the velocity perturbations, cf. Fig. 8(b). As a whole, in the central in-flow region around Y = 0 m, the original velocity profiles and magnitudes are conserved when the flow limiter is added. Let's notice that, in the down-comer, the general direction of the incoming flow from the primary pipe is essentially axial (small velocity X-components) with or without the flow limiter.

Nominal-operation condition

Fig. 9 shows the mesh convergence concerning the pressure and velocity profiles. The profiles obtained with the M 3 (164,160 elements) and M 2 (48,640 elements) meshes are relatively closed in comparison with the profiles obtained with the coarse mesh M 1 (6,080 elements). Even if the mesh convergence is not completely reached, it gives confidence in the physical tendencies deduced from the fine-mesh result analysis. - Table 2. Obviously the geometry of the fins has to be optimized to enhance this effect while limiting the flow-limiter impact during nominal operations.

Moreover, considering Fig. 7, we conjuncture that the impact on the flow of the fins located at the bottom of the flow limiter is much greater during LB-LOCA condition than during nominal-operation one. (a S = a; L T = 1.0 m). Fig. 13 shows the mesh convergence concerning the pressure and velocity profiles. Here again, the computation using the coarse mesh M 1 is clearly underresolved in space. In contrast, the profiles obtained with the meshes M 2 and M 3 are rather close, which indicates that the mesh convergence is almost reached.

Two-phase flow

To illustrate the GENEPI capacities to deal with two-phase flows, a computation has been done for the case of a two-phase reverse flow (a S = a; L T = 1.0 m) without flow limiter and meshes M 1 and M 3 . For this computation, the enthalpy balance equation has been solved considering a turbulent Prandtl number of 0.5. CATHARE results provide the out-flow pressure value of 43.65 bar and the in-flow mixture specific enthalpy value of 1,291 kJ/kg under a void fraction equals to 0.42. The in-flow mass flow rate is unchanged. Figs 14 and15 illustrate the void fraction field and Fig. 16 shows a comparison of the velocity fields for liquid and two-phase flows. As expected, as the mixture density is lower than the liquid one, the velocity increases.

In this case, the global GENEPI-computed two-phase pressure-drop coefficients are about K M 1,2ϕ = 3.2 and K M 3,2ϕ = 3.0; a little less than (but not very different from) the one-phase ones: K M 1 = 3.9 and K M 3 = 3.7, cf Table 2. Fig. 17 shows a comparison of the pressure and velocity profiles along the line (0. m, [-2.0; 2.0] m, 0.1 m) for two-phase and one-phase flows. As previously underlined, in the case of two-phase flows, the velocity magnitude is strongly increased as well as the pressure drop (about 3 bar). • For the reverse flow direction and without in-vessel flow limiter, the values of the global pressure-drop coefficients provided by GENEPI (K M 3 ∈ [3.6; 3.7] and K M 3,2ϕ ≈ 3.0) may suggest that the CATHARE local pressuredrop coefficient for reverse flow direction are about K out ≈ 7. This can be a guideline for CATHARE studies as ones done in [START_REF] Stratta | Thermal-hydraulic study of passive safety systems based on the hydraulic diode principle for the management of large-break loss of coolant accidents[END_REF].

• For the reverse flow direction and with in-vessel flow limiter, the CFD global pressure-drop coefficient is multiplied by almost a factor two (K M 3 ∈ [3.6; 3.7] → [5.4; 6.4]). This leads to CATHARE local pressure-drop coefficient K out ≈ 25. This is higher than the lower values exhibited by parametric CATHARE computations [START_REF] Stratta | Thermal-hydraulic study of passive safety systems based on the hydraulic diode principle for the management of large-break loss of coolant accidents[END_REF] to get a beneficial effect of the invessel flow limiter on LB-LOCA transients: K gl,cath ≈ 3.9 for K out ≈ 11. take into account the in-vessel flow limiter in a realistic manner. Moreover, the range of these estimated values is higher than the considered one in previously CATHARE studies [START_REF] Stratta | Thermal-hydraulic study of passive safety systems based on the hydraulic diode principle for the management of large-break loss of coolant accidents[END_REF], that is conservative from a point of view of safety analysis. Nevertheless, considering the limitations related to the turbulence model, this conclusion needs to be consolidated by body-fitted CFD studies with more precise turbulence models. Some perspectives can be outlined about the improvement of this numerical model in order to proceed to the geometry optimization of the design of the flow-limiter fins. In particular, we can mentioned the space-interpolation scheme across the boundary interface to reach the second order [START_REF] Introïni | A Second Order Penalized Direct Forcing for Hybrid Cartesian/Immersed Boundary Flow Simulations[END_REF] and the definition of immersed-wall laws for RANS/large-eddy simulations.

1 .

 1 No moving working fluid, 2. No moving mechanical part, 3. No input signal of "intelligence", 4. No external power input or forces.

Figure 1 :

 1 Figure 1: Scheme of in-vessel flow limiters (hydraulic diode) located between the cold legs and the downcomer [4].

3 m x 4 m

 4 to which was added the broken cold-leg nozzle starting about 1 m before the down-comer. It extends up to 2 m below the cold-leg axis, cf. Fig.2. Considering a sectorized down-comer as in Fig.1, we only model one third of the down-comer, including one cold-leg entry in the vessel and one hot-leg pipe, cf.

Fig. 3 .

 3 Fig. 3. As a whole, the geometrical dimensions are taken from the CATHARE model, except for the thickness of the planar part. A compromise between the volumes and the radius leads us to set this data to about 0.2 m. The measures of the down-comer and the cold-leg nozzle surfaces are equal to 0.9 m 2 and

Figure 2 :

 2 Figure 2: Part of the PWR geometry taken into account by the CFD study of the in-vessel flow limiter.

Figure 3 :

 3 Figure 3: Example of the mesh used for the CFD study of the in-vessel flow limiter (mesh M 2 ; 48640 elements). The walls are colored in blue.

Figure 4 :

 4 Figure 4: Example of preliminary meshing of the flow limiter.

4. 2 .

 2 Boundary conditions Mass flux Q in is imposed on the in-flow boundary and ad-hoc pressure P out = 50 bar on the out-flow boundary. The in-flow and out-flow surfaces are the down-comer and the cold-leg nozzle surfaces depending on the considered main-flow direction (default, as in nominal-operation condition, or reverse, as in LB-LOCA condition). Slip-wall boundary conditions are considered on the walls. The inlet mass fluxes are given by the steady-state CATHARE computation for the default flow direction and by CATHARE transient results at time t = 10 s after the guillotine break opening for reverse flow direction. For the reverse flow direction case, this instant is chosen so that the equations of state of GENEPI and CATHARE are compatible (primary pressure in the range [40-60] bar). At this time, CATHARE provides the pressure value of 64.6 bar. The inlet mass-flux values Q in are • Nominal-operation condition (default flow direction): Q in = 4,690 kg/s, in ↔ cold leg and out ↔ downcomer, • LB-LOCA condition (reverse flow direction): Q in = 5,200 kg/s, in ↔ downcomer and out ↔ cold leg. 4.3. Computation hypothesis For simplicity reason, we do not give attention to the gravitational acceleration term in the GENEPI or CATHARE computations (around the flow limiter, gravitational effects are negligible in comparison with the inertial ones during the fast-depressurization phase). As well, the pressure range considered here is only [40-60] bar. It differs from the pressure addressed in the CATHARE simulations of the nominal operation or of the LB-LOCA. But it is the pressure drop that is important at first, mainly for one-phase flow computations. Concerning the turbulence model, cf. Eq (4), the GENEPI-code standard value for the Schlichting coefficient is a S = 0.047 and the turbulence characteristic length L T is related to the biggest eddy structures. As reference, we choose L T ≈ 1 m for the nominal-operation flow direction (azimuth scale in the downcomer) and L T ≈ 0.3 m for reverse flow direction toward the broken cold leg (≈ radial scale in the down-comer).5. 3D/1D up-scaling of the global pressure-drop coefficient 5.1. The CATHARE model In our CATHARE model, the cold-legs are modeled by three axial elements and the down-comer by a single axial element connected to a volume element, cf. Fig. 2 where only one cold-leg is shown. On this volume are defined the local pressure-drop coefficients for the junctions that connect this element to the cold-legs: K in for fluid flowing into the volume (default flow direction, nominaloperation condition) and K out for fluid flowing out (reverse flow).

Figs 5 to 11

 11 Figs 5 to 11 present some field distributions concerning the mixture velocity, the mixture pressure and the local external forces taking into account the flowlimiter fins ΛV, cf. Eq. (2) for the two studied flow configurations (nominal condition and reverse condition). For each configuration, we compare the cases with and without flow limiter for the reference turbulence-model coefficients.

.

  The default flow direction is defined as the nominal-operation flow direction (from the cold leg to the down-comer). The reverse flow direction is defined as the opposite direction. (*): unsteady computation. BTD:

For

  Fig. 7(a). Introducing the flow-limiter device, the regular standard flow path becomes much more irregular with the emergence of flow channels between the fins, cf. Fig. 6(b). Also, the in/out-flow pressure difference is increased in reference to the case without in-vessel flow limiter: K M 3 ∈ [-0.3; +0.1] →

( a )

 a Without flow limiter. (b) With flow limiter.

Figure 5 :

 5 Figure 5: Comparison of the pressure in nominal-operation condition; mesh M 3 (164,160 elements); (a S = a; L T = 0.3 m). Pressure iso-values ranges from 4.8 to 5.2 bar.

Figure 6 :

 6 Figure 6: Comparison of the velocity in nominal-operation condition; mesh M 3 (164,160 elements); (a S = a; L T = 0.3 m).

  (a) Nominal-operation flow; (a S = a; L T = 0.3 m). (b) Reverse flow; (a S = a; L T = 1.0 m).

Figure 7 :

 7 Figure 7: Local external-force distribution; mesh M 3 (164,160 elements).

  m/s) Vx with fl. lim. Vx without fl. lim. Vy with fl. lim. Vy without fl. lim. Nominal-operation conditions Mesh M3, L_T=0.3 m Fins region Flow-limiter outer region Flow-limiter outer region (b) Velocity X/Y -components.

Figure 8 :

 8 Figure 8: Comparison of the pressure and velocity profiles in nominal-operation condition; mesh M 3 (164,160 elements); (a S = a; L T = 0.3 m).

  outer region Flow-limiter outer region (b) Velocity X-component.

  outer region Flow-limiter outer region (c) Velocity Y -component.

Figure 9 :

 9 Figure 9: Mesh convergence of the pressure and velocity profiles in nominal-operation condition; (a S = a; L T = 0.3 m).

6. 2 .

 2 Fig. 7(b). Contrary to the nominal-operation case, the range of CF-computed values without flow limiter K M 3 ∈ [3.6; 3.7] sensibly differs from the Idel'cik (conical collector with front wall and screen) or Borda-Carnot (sudden flow-section reduction) estimation one K gl ∈ [1.1; 1.2]. But the hydraulic path is quite complex and does not reduce to simple configurations. As shown, the introduction of the flow-limiter device clearly increases the fluid vortex at the entry of the broken cold leg, cf. Fig. 11(b), and the in/out-flow pressure drop, cf. Fig. 10(b). The global pressure-drop coefficient is multiplied by almost a factor two in case of flow limiter: K M 3 ∈ [3.6; 3.7] → [5.4; 6.4], cf.

Fig. 12

 12 Fig. 12 shows a comparison of the pressure and velocity profiles with and without the flow limiter along the line (0. m, [-2.0; 2.0] m, 0.1 m), cf. Fig. 4. Here again, in the central out-flow region around Y = 0 m, the original profiles and magnitudes of the pressure and the velocity are conserved when the flow limiter is added. Whatever the situation is, there is a vortex of out-coming flow moving toward the primary pipe exit. But with the flow limiter, due to the presence of the fins, the velocity rotation is reinforced, cf. the plots of the velocity X-component in the fins region (Y ∈ [-1.5; -1.1] m) displayed in Fig.12(b). Hence, the flow-limiter induced pressure drop magnitude is relatively high: about 3.5 bar.

  (a) Without flow limiter. (b) With flow limiter.

Figure 10 :

 10 Figure 10: Comparison of the pressure in LB-LOCA condition; mesh M 3 (164,160 elements); (a S = a; L T = 1.0 m). Pressure iso-values ranges from 4.7 to 5.9 bar.

Figure 11 :

 11 Figure 11: Comparison of the velocity in LB-LOCA condition; mesh M 3 (164,160 elements); (a S = a; L T = 1.0 m).

  m/s) Vx with fl. lim. Vx without fl. lim. Vy with fl. lim. Vy without fl. lim. LB-LOCA conditions Mesh M3, L_T=1.0 m Fins region Flow-limiter outer region Flow-limiter outer region (b) Velocity X/Y -components.

Figure 12 :

 12 Figure 12: Comparison of the pressure and velocity profiles in LB-LOCA condition; mesh M 3 (164,160 elements); (a S = a; L T = 1.0 m).

  outer region Flow-limiter outer region (b) Velocity X-component.

  outer region Flow-limiter outer region (c) Velocity Y -component.

Figure 13 :

 13 Figure 13: Mesh convergence of the pressure and velocity profiles in LB-LOCA condition;

Figure 14 :

 14 Figure 14: Two-phase flow CFD study of the in-vessel reverse flow without flow limiter: void fraction; mesh M 1 (6,080 elements); (a S = a; L T = 1.0 m).

Figure 15 :

 15 Figure 15: Two-phase flow CFD study of the in-vessel reverse flow without flow limiter: void fraction profile along the line (0. m, [-2.0; 2.0] m, 0.1 m); mesh M 1 (6,080 elements); (a S = a; L T = 1.0 m).

  (a) Liquid flow. (b) Two-phase flow.

Figure 16 :

 16 Figure 16: Two-phase flow CFD study of the in-vessel reverse flow without flow limiter: comparison of the velocity; mesh M 1 (6,080 elements); (a S = a; L T = 1.0 m).

  m/s) Vx two-phase flow Vx one-phase flow Vy two-phase flow Vy one-phase flow LB-LOCA conditions Mesh M1, L_T=1.0 m (b) Velocity X/Y -components.

Figure 17 :

 17 Figure 17: Two-phase flow CFD study of the in-vessel reverse flow without flow limiter: comparison of the pressure and velocity profiles; mesh M 1 (6,080 elements); (a S = a; L T = 1.0 m).

6. 4 .

 4 CATHARE local pressure-drop coefficients Here, we try to evaluate the ad hoc CATHARE local pressure-drop coefficients K in and K out leading to global pressure-drop coefficients K gl,cath compatible with the above global CFD estimated values. Let us recall that the coefficients K in and K out are located at the junction between the down-comer volume element and the cold-leg axial elements of the CATHARE model, cf Section 5.1. This estimation of the CATHARE local pressure-drop coefficients is done for the cases with and without in-vessel flow limiter. Results are summarized in Table 3. • For the default flow direction and without in-vessel flow limiter, the global pressure-drop coefficients provided by GENEPI and CATHARE are similar (K M 3 ∈ [-0.3; +0.1] versus K gl,cath ∈ [-0.3; -0.1]) giving confidence in the CATHARE local pressure-drop coefficient for nominaloperation flow direction defined by K in ≈ 1 -1.8 = -0.8 [8].

7 .

 7 Concerning the reverse-flow case (LB-LOCA condition), two-phase-flow GENEPI simulations should refine the estimation of the global pressure-drop coefficients and, consequently, the CATHARE local pressure-drop coefficient K out . model of a PWR. Parametric studies on the turbulence model lead to determine the range of the global pressure-drop coefficients depending on the direction of the flow and the presence of the flow limiter. Estimations of the CATHARE local pressure-drop coefficients have been done opening the way to

Table 1 :

 1 Hydrodynamic coefficients associated with the problem of steady flow around a static cylinder of diameter D (Reynolds number 20). B.F.: Body-fitted. C d : drag coefficient.

	References

Table 2 :

 2 Summary of the global pressure-drop coefficients K gl

Conclusions and perspectives

In this paper, we have presented a CFD model of an in-vessel flow limiter to mitigate the consequences of a large-break loss of coolant accident in a pressurized-water reactor. The principle of this safety device is based on fins designed to create a strong flow vortex increasing the pressure drop toward the broken cold leg. A rough CFD model using a homogeneous relaxed equilibrium model of a liquid-vapour mixture and an immersed boundary approach has been set-up, allowing a less precise but fast estimation of the pressure drop following the geometry of the fins. Through 3D/1D up-scaling of a global pressure-drop coefficient, local pressure-drop coefficients can be provided to thermal-hydraulic system safety codes (here CATHARE), allowing the study of the in-vessel flow limiter effect on the thermal-hydraulic system.

One-phase and two-phase fluid CFD simulations have been run using the GENEPI code with a computation domain defined in coherence with the CATHARE

Nomenclature and acronyms

Latin symbols

• a S : Schlichting model dimensionless constant

• A e : obstacle area intercepted by the element e (m 2 )

• AGR: Advanced Gas-cooled Reactor • G: mixture mass flux (= ρV)

• H: mixture specific enthalpy (J kg -1 )

• H ls : saturated liquid specific enthalpy (J kg -1 )

• IB: Immersed Boundary

• IBC: Immersed Boundary Condition

• ISI: Immersed Spread Interface

• L: latent heat (J kg -1 )

• L T : typical vortex length (m)

• L w : recirculation length (m)

• LB LOCA: Large Break Loss Of Coolant Accident

• M i : domain mesh

• P : pressure (Pa)

• PWR: Pressurized Water Reactor

• P4: 1,300 MWe French PWR (Palier P4)

• SCRAM: Safety Control Rod Axe Man

• SI: Safety Injection

• SUPG: Streamline Upwind Petrov-Galerkin method

• t: time (s)

• V: mixture velocity (m s -1 )

• v R : relative velocity (gas minus liquid, m s -1 )

• VVER: Vodo-Vodiano Energuetitcheski Reaktor

Greek symbols

• : Penalty parameter (0 < << 1)

• χ T : turbulent diffusion coefficient for the enthalpy balance equation (kg m -1 s -1 )

• Λ: two-phase friction tensor (s -1 )

595

• µ T : two-phase turbulent dynamic viscosity (N s m -2 )

• Ω e : elementary volume of the element e (m 3 )

• ρ: mixture density (kg m -3 )