Thibaut Vuillemin

François Goichon

Cédric Lauradoux

Guillaume Salagnac

Entropy transfers in the Linux Random Number Generator

Keywords: random number generation, operating systems, linux, entropy génération de nombres aléatoires, systèmes d'exploitation, entropie

One of the services provided by the operating system to the applications is random number generation. For security reasons, the Linux Random Number Generator is built upon the combination of a deterministic algorithm known as the cryptographic post-processing and an unpredictable physical phenomenon called an Entropy Source. While the various cryptographic post-processing algorithms and their properties are well described in the literature, the entropy collection process itself is little studied. This report first presents the different approaches to random number generation, and then details the architecture of the Linux Random Number Generator. Then, we present the experiments we performed to monitor entropy transfers. Our results show that the main source of randomness in the system is the behavior of the hard drive, and that most random numbers produced by the generator are actually consumed by the kernel itself.

Introduction

On a computer system, access to hardware devices like the processors, main memory, or network interfaces has to be properly shared between all the users. Not only has the operating system (OS) to provide software interfaces between these hardware resources and the user programs, but also the OS should enforce fairness and other desirable properties amongst users. Otherwise, a malicious user program may try and monopolize a particular resource like the disk drive or the CPU, which damages system performance and may lead to denial of service. Therefore, understanding and controlling the behavior of all shared resources is of critical importance to evaluate resiliency of the system.

One little-studied shared OS resource is the random number generator. Many applications need random numbers. For example, all security-oriented programs use random numbers to generate reliable cryptographic keys. Also, random numbers are key to implementing Monte-Carlo simulations. Because algorithmic generators are too weak to satisfy the randomness needs of all those applications, a majority of operating systems provide a system-level Random Number Generator (RNG). The Linux kernel is not an exception and contains such a generator: the Linux Random Number Generator (LRNG). This component is available both for use by the kernel itself and by user-space applications. A non-privileged user of a GNU/Linux OS obtains random numbers by reading the /dev/random or /dev/urandom devices.

The Linux Random Number Generator is not a mere deterministic algorithm producing numbers. In order to offer "good quality" randomness, the LRNG implements a so-called True Random Number Generator architecture, which collects so-called entropy from various parts of the system. These terms will be defined in Section 2. However, a major drawback of this architecture is that the generator sometimes has to wait for some entropy to arrive, thus yielding very long response times.

To illustrate this situation, we conducted a little experiment, with a program that repeatedly requests one byte from the /dev/random device. The results are presented on Fig. 1 below 1 . While the shortest request completion time is under 1µs, the longest response time is about 35 seconds. That is a seven order-of-magnitude difference.

The litterature about random number generation techniques, and about the different generator architectures, is numerous. However, to the best of our knowledge, little has been done on the entropy collection process itself. In this work, we investigate this subject by studying how Linux generates, collects, and handles entropy to be used by the LRNG. The report is organized as follows. Section 2 defines several terms related to random number generation (such as entropy) and presents the different RNG architectures. Section 3 gives an overview of the Linux Random Number Generator. Section 4 describes the monitoring system we implemented to understand to the LRNG and presents the results we obtained. Finally, the Section 6 summarizes our conclusions and perspectives.

Context: random numbers 2.1 Random numbers characterization

There are various ways to define or describe what is called randomness. Two important properties are often expected for a sequence of numbers to be called "random", that is uniform distribution, and statistical independence. Let X 1 , X 2 , • • • X n be a sequence of random variables defined from the sample space Ω realizing their values in a measurable state space X . Let ω ∈ Ω be an event of the sample space. Let x i ∈ X be the i-th realization of X i such that X i (ω) = x i (i.e the definition of a random variable). Let P r(X i = x i ) be the probability of the random variable X i realizing the value x i . Equation 1 defines that X i follows a uniform distribution. Equation 2 defines statistical independence between two random variables X i and X i-1 .

P r(X

i = x i ) = 1 card(Ω) (1) ∀i ∈ {2, ..n}, P r(X i = x i |X i-1 = x i-1 ) = P r(X i = x i ) (2)
However, not all random variables satisfy these two properties. When it is not the case, then the numbers we get are not "purely random" but may still contain some amount of uncertainty. This amount is often named entropy and can be described by the following formula.

Let X be a random variable. Let X be the sample space of this random variable. Let x be a realization of this variable, with x ∈ X . For any i ∈ X , let p(i) = P r(x = i) be the probability of x being i. The Shannon entropy of variable X is defined as:

H(X) = - ∀i∈X p(i) log 2 p(i).
In order to use this formula, one has to know the probability distribution function p. If this function is not known, then the Shannon entropy can only be estimated from several realizations of the variable. In the literature, this problem of entropy estimation is an open topic [START_REF] Lauradoux | Online entropy estimation for non-binary sources and applications on iphone[END_REF].

Inria

Random Number Generators

A Random Number Generator is a computer program intended to behave like a random variable. Random number generation always attracted attention from both computer science and mathematics communities. The work of Knuth [START_REF] Knuth | The art of computer programming[END_REF] is seminal in this field. In the literature, different classes of random number generators have been studied:

• pseudo-random number generators (PRNG);

• cryptographically secure pseudo-random number generators (CSPRNG);

• true random number generators (TRNG).

We choose to give a simple description for each of them. Our descriptions are close to the ones given by Viega in his book Secure Programming [START_REF] Viega | Secure Programming Cookbook for C and C++[END_REF].

Pseudo-Random Number Generators

A Pseudo-Random Number Generator (PRNG) is a finite state machine. The initial value of the internal state is often called the seed. Upon each request, a transition function computes the next internal state, and an output function (sometimes the identity function) produces the actual number.

A PRNG will deterministically produce a sequence of values depending only on the initial seed. This sequence will always be periodic, because there are only so many distincts states the algorithm can visit before rolling over.

Remark If the internal state is stored as a n-bit vector, then the maximum period for the PRNG is 2 n . Indeed, a finite state machine has no other internal memory than its state, so after at most 2 n transitions it will necessarily end up in some already visited state, and from then on it will perform the same transitions again and again.

Example Linear congruential generators [START_REF] Park | Random number generators: Good ones are hard to find[END_REF] are a common class of PRNGs. Let x t be the internal state of the PRNG at time t. Let a, b and m be arbitrary integer constants. The PRNG transition function is:

x t+1 = a • x t + b mod m.

Cryptographically Secure Pseudo-Random Number Generators

A Cryptographically Secure Pseudo-Random Number Generator (CSPRNG) is a PRNG that satisfies cryptographic conditions [R 05] in addition to statistical conditions. Two classes of CSPRNGs can be distinguished depending on their security model: so-called stream ciphers [START_REF] Matthew | New Stream Cipher Designs -The eSTREAM Finalists[END_REF] offer empirical security by resisting to the state of the art of cryptanalysis attacks, while provably secure CSPRNGs have more theoretical security properties derived by reduction to a difficult problem [START_REF] Blum | A simple unpredictable pseudo-random number generator[END_REF].

As stated by Viega [START_REF] Viega | Secure Programming Cookbook for C and C++[END_REF]:

Cryptographic pseudo-random number generators are still predictable if you somehow know their internal state. The difference is that, assuming the generator was seeded with sufficient entropy and assuming the cryptographic algorithms have the security properties they are expected to have, cryptographic generators do not quickly reveal significant amounts of their internal state. Such generators are capable of producing a lot of output before you need to start worrying about attacks.

True Random Number Generators

Opposed to a PRNG, which is implemented entirely in software, a hardware RNG is based on some unpredictable physical phenomenon called an entropy source [START_REF] Kelsey | Entropy sources[END_REF]. Sampling this source is the electronic equivalent of rolling a die or flipping a coin. Examples of such generators include the Bull Mountain RNG implemented in recent Intel CPUs [START_REF]Bull mountain software implementation guide[END_REF] which uses thermal noise within the silicon as an entropy source. However, this architecture is limited by the speed of the physical phenomenon itself. Naive hardware RNGs thus typically offer insufficient throughput for most applications. This is why people generally build so-called True Random Number Generators [START_REF] Barak | A model and architecture for pseudo-random generation with applications to /dev/random[END_REF] by cascading a hardware RNG with a PRNG: the slow-but-very-unpredictable hardware RNG is used to periodically refresh the internal state of the fast-but-deterministic PRNG, in a process called reseeding. In this context, the hardware RNG is often referred to as a entropy harvester [START_REF] Barak | Extracting randomness using few independent sources[END_REF].

The resulting behavior is between that of a pure hardware RNG and a PRNG: when sufficient entropy is harvested from the source, then the PRNG is reseeded often and the output numbers are really non-deterministic. However, if the TRNG fails to harvest enough entropy from the source, then it does its best and acts as a pure PRNG.

Context: the Linux Random Number Generator

The Random Number Generator implemented by the Linux kernel [MT12, GPR06] is a True RNG: it contains both a PRNG and entropy harvesting mechanisms. In this section, we present the architecture of the LRNG as well as entropy transfer mechanisms.

Output interfaces

The Linux Random Number Generator provides three output interfaces:

• /dev/random • /dev/urandom • get_random_bytes()

The get_random_bytes() function is available only from within the kernel. It always return the requested amount of random data. On the other hand, the two /dev/(u)random devices are accessible from user-space, even to a non-privileged user. Whereas /dev/random can block the requesting process, /dev/urandom never blocks.

The Linux source code [START_REF] Mackall | random.c -a strong random number generator[END_REF] reads that /dev/urandom will return as many bytes as are requested. If too many many bytes are required, this will result in "random numbers that are merely cryptographically strong, but they are still acceptable for many applications". On the other hand, "/dev/random is suitable for use when very high quality randomness is desired (for example, for key generation or one-time pads), but it will only return a maximum of the number of bits of randomness contained in the entropy pool".

Entropy pools

The internal state is separated into three arrays of unsigned integers called the entropy pools: the input pool, the blocking pool, and the non-blocking pool.

The blocking pool, which size is 1024 bits, is related to /dev/random. When bits are requested to /dev/random, the LRNG reads numbers is this pool. The non-blocking pool, which size is also 1024 bits, is related to /dev/urandom. When bits are requested to /dev/urandom, the LRNG Inria reads numbers is this pool. It is also used by another function named get_random_bytes(). This function can be called from inside the kernel to request random numbers. The input pool, which size is 4096 bits, is connected to the blocking pool, the non-blocking pool and the entropy sources.

Every pool is associated with an integer called entropy counter. It is an estimation of the amount of entropy contained in the pool. The entropy counter of a pool varies over time, as explained below:

Random numbers request on the blocking pool When n random bits are requested to /dev/random, the blocking pool entropy counter is decreased by n. If it falls to zero, the read operation blocks until more entropy is acquired.

Random numbers request on the non-blocking pool

When n random bits are requested to /dev/urandom or get_random_bytes(), the non-blocking pool entropy counter is decreased by n. However, if it falls to zero, it will not block the request and the device will act as a PRNG.

Entropy transfer When the blocking or the non-blocking pool has to provide n entropy bits but has an entropy counter lower than n, then entropy is transferred from the input pool to the demanding pool. The input pool counter is decreased by n while the counter of the demanding pool is increased by n. If the input pool entropy counter is lower than n, then the transfer is rejected.

Entropy harvesting As a TRNG, the LRNG uses entropy sources to harvest entropy (see Section 2). When n entropy bits are gathered from an entropy source, the input pool is refreshed using those bits and the input pool entropy counter is increased by n.

Entropy sources

The LRNG depends on the rest of the kernel to provide it with entropy. It exports three functions, which are intended to be called by various device drivers:

• add_disk_randomness(),

• add_input_randomness(), and • add_interrupt_randomness().

The add_disk_randomness() function is called by the hard disk driver. It uses the the seek time of block layer request events as an entropy source. The add_input_randomness() function is called by the drivers of input peripherals such as the mouse, the keyboard, but also a joystick or a tablet. For every event, it uses values representing the event (e.g. "KEY" event, "A" key, key released) as an entropy source. The add_interrupt_randomness() function can a priori be called in any part of the kernel that generates interrupt. For every interruption, it uses the inter-interrupt timing as an entropy source. Every interruption is not satisfying: the inter-interrupt timing of a timer is predictable. We provide a big picture of the LRNG architecture in figure 2, including the entropy sources, the entropy pools and the output interfaces.

Input

Entropy estimation

In Subsection 3.2, we stated that when n entropy bits are gathered from a source, the input pool entropy counter is increased by n. We now explain how the LRNG computes n, that is, how it estimates how much uncertainty was conveyed by the event. The LRNG has to do such an estimation in order to detect regularities and reject non-uncertain sequences.

Let us consider an entropy source that produces events for the LRNG. Let t 0 , t 1 , t 2 , ... be the jiffies associated with those events. A jiffy is the current value of the kernel timer when an event happens.

First, the estimator calculates the jiffies differences on three levels:

δ i = t i -t i-1 , δ 2 i = δ i -δ i-1 , δ 3 i = δ 2 i -δ 2 i-1 .
Then, the estimator takes the minimum of the differences absolute values:

∆ i = min(|δ i |, |δ 2 i |, |δ 3 i |).
Finally, it applies the following function [START_REF] Lacharme | The Linux Pseudorandom Number Generator Revisited[END_REF]:

H i =    0 if ∆ i < 2 11 if ∆ i ≥ 2 12 log 2 (∆ i ) otherwise.
The choice of this particular estimator and of its parameters are not justified in the Linux source commentaries. It is a crude but cheap entropy estimator. It was chosen for its cost, not for its accuracy [START_REF] Lacharme | The Linux Pseudorandom Number Generator Revisited[END_REF]. Still, it is pretty good at detecting regularities. A retro-engineering [START_REF] Pousse | An interpretation of the linux entropy estimator[END_REF] study proposes an interpretation based on Newton polynomial interpolation.

Let us give an example:

∆ i = min(|11|, |6|, |2|) = 2 H 1041 = log 2 (2) = 1
According to the LRNG estimation, the event at time 1041 brings 1 bit of entropy.

RR n°8060

Alternatives to the LRNG

The LRNG is not the only random number generator that can be used on a GNU/Linux operating system. Many programming languages like C or Python include a random number generator in their standard library, but they are generally PRNG, not TRNG.

The Entropy Gathering Daemon [START_REF] Warner | EGD: The Entropy Gathering Daemon[END_REF] is a user-space TRNG, that is intended to be used on a Linux system running applications that need random numbers but for some reason are denied access to /dev/(u)random.

HAVEGE [START_REF] Seznec | HAVEGE: A user-level software heuristic for generating empirically strong random numbers[END_REF] is a TRNG that uses the state of volatile hardware components (caches, branch predictors, long pipelines, instruction level parallelism) as entropy sources.

The kernel module frandom is a kind of add-on to the LRNG. It provides the system with two new user-space output interfaces: /dev/frandom is 10-50 times faster than /dev/urandom but may consume some of the LRNG entropy, /dev/erandom is slower but does not consume any LRNG entropy. Both of them are not meant to be used for cryptographic purposes.

Contribution

In the previous part, we studied the architecture of the LRNG, including its components and algorithms. We know the LRNG principles. Now we would like to investigate the behavior of the LRNG while in use. Consequently, we implemented a system that monitors the inputs, the outputs and the internal state of the LRNG, and ran it in several test scenarios.

Experimental setup

We used virtual machines in order to develop our system, because it is both faster and safer to test a kernel on a virtual machine rather than on a real computer. We use the QEmu emulator to host the virtual machines, because it allows us to change the kernel image very easily at each boot.

However, we used real machines for the actual experiments, in order to avoid any influence from the emulator on the results.

We chose to implement the monitoring system directly inside the kernel. We could have used a kernel debugger, but it would have been very difficult to know how it modifies the LRNG behavior. By using our system, we have both a better knowledge and a better control on what happens.

The easiest way to log something from within the kernel is to use the printk function. It works like printf, and writes the output in the operating system log files. We chose not to use this function, because it could have generated disk events that would have been used as an entropy source by the LRNG. We want our monitoring system to have as less influence as possible. Instead of printk, we chose to send packets over network.

To our knowledge, there are two different kernel APIs to send UDP packets, available via include/linux/net.h and include/linux/netpoll.h. The first one is pretty similar to the standard userspace socket API. The second one is more low-level and not commonly used. Still, we chose to use it, because it works even in IRQ contexts, which is required to instrument input events.

Inria

Tested scenarios

We ran three experiments, corresponding to different scenarios:

1. a desktop workstation, 2. a file server, 3. a computing server.

In the desktop scenario, we installed a graphical environment and a web browser on the virtual machine. Once ready, we rebooted it and started to monitor its LRNG activity. For one hour, we did some typical desktop tasks, mainly reading articles on the internet and installed a few programs. Finally, we shut it down properly, in order to include machine booting and halting activities in the results.

In the file server scenario, we did something similar (one hour, from boot to shutdown). There were no graphical packages installed on the system, and we ran wget to download a large data file from a Debian image mirror. Approximately 1.6 GB of data was transferred during the experiment.

In the computing server scenario, we launched a mathematical C program, which looks for integers n ∈ N such that 2 √ n -1 ∈ N and 3 √ n + 1 ∈ N. For instance, the number 26 is such an integer (26 -1 = 25 = 5 2 and 26 + 1 = 27 = 3 3 ).

We chose those scenarios because they were likely to stress different entropy sources of the LRNG. The desktop workstation is likely to generate mouse, keyboard and other input events; the file server mainly uses the disk and the network; and the computing server has a heavy CPU load.

Entropy inputs

In this experiment, we study the origin of the entropy harvested by the LRNG. For that we measure the proportions of entropy that comes from the three functions:

• add_disk_randomness(),

• add_input_randomness(), and

• add_interrupt_randomness()

We summarize the results in the table below, and then discuss the lessons learned from our observations. A graphical representation of those results can be found in Fig. 3. The term "generic input" represents an input event that cannot be associated with a particular device. Secondary result: In our three scenario, we did not observe a single "IRQ event". We assume than our experiments did not include any activity likely to generate interrupts that can be used by the add_interrupt_randomness() function, such as an interaction with a USB-key.

Desktop workstation scenario

Secondary result:

In the server scenarios, the only entropy provider is the disk. In those scenarios, all the interactions have been performed remotely, without using the machine mouse or keyboard, so there was no reason to observe any entropy generated but an input device.

Secondary result:

The desktop workstation obviously produces far more entropy than the server scenarios. The amount of collected entropy in this scenario is almost 5 times higher compared the file server scenario, and more than 15 times higher compared the computing scenario. As expected, the more user activity there is, the more entropy the system collects.

Entropy outputs

In this experiment, we measure the quantity of entropy that flows through the three output interfaces:

• /dev/random,

• /dev/urandom, and

• get_random_bytes()

We summarize the results in the table below, and a graphical representation of those results can be found in Fig. 4.

Desktop workstation scenario

Primary result: In our three scenarios, most entropy consumption is caused by requests on get_random_bytes(). It means the kernel is the first entropy consumer.

Secondary result: The total amount of consumed entropy almost equals the total amount of produced entropy, in every scenario. This result was to be expected, because 1) obviously the system cannot consume more entropy than it has collected, and 2) the LRNG does not produce useless entropy : once the entropy counter for the input pool has reached its maximum value (4096 bits), the LRNG starts ignoring any incoming events in order to save CPU resources. Thus, it collects no more entropy.

Entropy consuming applications

In this experiment, we study the userspace clients requesting data on the /dev/(u)random devices, in order to find out which applications need random numbers and how much entropy they consume. Due to major differences between the results in each scenario, we choose to present them one after another. We summarize the results in the tables below, and graphical representations are available in figures 5, 6 and 7. The main entropy consumer is the ELF binary loader, which is part of the kernel (implemented in fs/binfmt_elf.c). The ELF file format is the default format for executable binaries on x86 Linux systems. The loader needs random numbers to implement Address Space Layout Randomization2 . Other important consumers are :

• SVN, because we were using an HTTPS repository ;

• the Chrome browser, mainly for use by various subsystems like TLS, SSL, (the networkrelated part of) ffmpeg, and python (for seeding its own "random" module)

We also note that a small but non-negligible part of the entropy is consumed by a process named "php5". This process is not part of a web server, but a mere php interpreter invoked during the start-up phase to run some system script.

The remaining 5% of entropy consumers include both userspace processes related to the desktop environment and kernel functions related to network. Once again, the main entropy consumer is the ELF binary loader. The php5 consumer is still there ; it is the second consumer, but the amount of consumed bytes is equivalent.

46% 26% 21% 2% 5% [K] load_elf_binary [U] svn [U] chromium-browse [U] php5
The third consumer is a kernel function called inet_frag_secret_rebuild() and located in the net/ipv4/inet_fragment.c file, which is part of the IPv4 stack. The implementation uses a hash table to store IP datagram fragments. If the hash function was the always the same, then it would be practical for an attacker to set up a denial of service attack just by sending well-known specific datagrams which caused hash collisions happen. That's why a random seed is used during the initialization in order to make the hash function unpredictable from the outside.

The fourth consumer is the Apache web server, even if there was no HTTPS request during the experiment. The server needs random numbers when it starts in order to initialize a seed used by the SSL protocol. In this scenario we find only two consumers, and both are kernel functions. The ELF binary loader is still the main consumer, eating an overwhelming majority of the entropy. Only 5% are required by the second consumer, a function of the IPv4 stack which has been presented before.

Inria 95% 5% load_elf_binary+0xad5 inet_frag_secret_rebuild+0x2c
Figure 7: Randomness consumers -Computing server

Subsection conclusion

In this subsection, we studied how entropy is consumed. For each request on /dev/(u)random, we identified the corresponding process. For each call to get_random_bytes(), we analyzed the corresponding kernel function. In this conclusion, we remind the results that seem really significant.

Primary result: In our three scenarios, the vast majority of entropy is consumed by one kernel function : the ELF binary loader.

Secondary result:

The other consumers vary greatly according to the scenario.

Entropy counter of the input pool

In this experiment, we measure the evolution over time of the input pool entropy counter. We present the results on Fig. 8.

The vertical axis represents the input pool level. It goes from 0 (empty pool) to around 4096, which is the size of the input pool bit field. The level cannot increase over this size : if all the bits of the bit field are unpredictable, no additional entropy can be stored.

Moreover, we denote that the entropy counter does not fall to zero exactly, and does not cap to 4096 exactly. The behavior is due to thresholds. When the entropy counter goes over 3584 bits, the LRNG starts dropping most samples to avoid wasting CPU time and reduce lock contention. On the other hand, if the entropy counter falls under 64 bits, entropy transfers are blocked until the input pool level reaches 64 bits, so that the transfer is significant enough.

The horizontal axis represents the time, expressed in CPU cycles. Even if the scenario duration is always the same (one hour), the number of CPU cycles is different according to the scenario. It is about 10 10 cycles in the server scenario, ten times higher in the computing scenario and again 10 times higher in the desktop workstation scenario. These results are logical regarding the to the recent CPUs sleep mode. On the one hand, the CPU load is extremely heavy in the computing scenario. One the other hand, it is very light in the file server scenario, because the computer is mainly waiting for the network. The desktop workstation scenario stands for an average between the two precedent cases.

Primary result:

The entropy counter varies greatly in the desktop scenario. It can go from its maximum to its minimum very quickly. The unstable behavior is due the constant user interaction. A long entropy shortage seems to be almost impossible in such a scenario.

Primary result:

In the file server scenario, the entropy counter increases continuously and slowly until the computer halts. We assume that the slow increase is due to the reduced number of CPU cycles, because the number of jiffies is also reduced and the LRNG estimator is based on the jiffies differences. If we would have represented the desktop workstation scenario and the file server scenario on the same X axis, the slopes would have been similar. We assume that a moderate entropy consumption would not cause any shortage.

Primary result:

In the computing scenario, the entropy counter remains extremely low during all the experiment, because both entropy production and consumption are really basic. The situation is almost a constant entropy shortage.

Inria

LRNG customers

Thanks to our experiments, we acquired information about the LRNG "customers", who use the LRNG random numbers and thus consume entropy. Mainly, we have got valuable information about the first of the customers: the kernel itself. However, even if we have some information about other customers, it is highly specific to our scenarios.

In this last part, we try to cross-check our results with the source code and the documentations of widely used applications:

• google_authenticator,

• openssl,

• openssh.

The Google authenticator is a enforced authenticator with a two-step verification, available for GMail or GApps users who feel concerned with security. It uses /dev/urandom to get random numbers. When someones suggests to switch from /dev/urandom to /dev/random be cause the latter is supposed to be more secure3 , the developers oppose a categorical refusal. One the one hand, the don't want to rely on a device that can block indefinitely. One the other hand, they argue that no practical attack has been shown on /dev/urandom despite its theoretical weaknesses.

The OpenSSH library only uses /dev/urandom, even when generating keys with ssh-keygen.4 . Once again, the developers don't want to rely on a blocking device. According to them, it is so much easier to compromise a key by social engineering rather than by an attack on /dev/urandom, it would be a useless loss of time to use /dev/random instead.

The OpenSSL library is able to use both devices, but it will always try to use the non-blocking one. Although no explanation is given in the documentation, we assume that the reasons are similar to the previous ones.

Primary result: Even if entropy shortage seems to be a common phenomenon according to our results, the end users never notice them because the LRNG "customers" only use the non-blocking interfaces: /dev/urandom and get_random_bytes().

Conclusions and perspectives

In this work, we studied the randomness resource in the GNU/Linux operating system with a system-oriented approach. Thanks to a built-in monitoring system, we did several measures on the random number generator of the Linux kernel in order to understand the underlying mechanism that transform unpredictable interactions like disk seeks of mouse events into reliable random numbers, meant to be used by security programs.

Even when entropy is rare, it does not prevent the kernel from running normally. This is due to the fact that the /dev/random device was never used in our experiments, whereas /dev/urandom and get_random_bytes() are commonly used and non-blocking. Between the two methods, i.e /dev/urandom and get_random_bytes(), the latter one is the most often called meaning that the kernel is consuming a lot of random numbers. Inside the kernel, the most entropy is consumed by the function used to load the binary files in memory. It needs random numbers in order to randomize the address space of the binaries for security reasons. In all our experiments, the disk is a significant entropy provider even in the desktop scenario which includes interactions with a graphical environment.

If we compare our experiments, we observe that the total amount of entropy gathered is far higher in the desktop scenario than in the other two. From the perspective of the entropy gathered, a file server produces only a bit more entropy than a computing server. However, the input pool level for the file server is far higher than in the case of the computing server. Therefore, a quantitative comparison is not fair. The desktop workstation scenario and the computing scenario are both likely to have entropy shortage. They are very brief for the desktop workstation and should not cause any problem. However, the computing scenario is almost always in situation of entropy shortage.

To our knowledge, no attack on /dev/urandom has been shown in the state of the art. However, the systematic use of /dev/urandom instead of /dev/random relies on the postulate that "there will be enough entropy soon". The supposition is quite true in the desktop workstation scenario; but according to our experiments, it is absolutely false in the situation of a computing server scenario, and not certain in a file server scenario. Besides, the overwhelming majority of Linux systems are servers, not desktop workstations. And on servers, the only entropy provider is the disk, whereas this device plays a less and less important role in modern systems. In the future, it might be not acceptable anymore to consider that there will always be enough entropy soon.

Figure 1 :

 1 Figure 1: Response time of /dev/random for 1000 successive one-byte requests. The average response time is 264.518 ms with a standard deviation of 1675.535 ms.

Figure 2 :

 2 Figure 2: Global architecture of the the LRNG. Lozenge boxes are entropy sources ; rectangle boxes are the three entropy pools, together making up the internal state of the generator. Rounded boxes are the three output interfaces. An arrow between two boxes represents the possibility of an entropy transfer.

Figure 3 :

 3 Figure 3: Randomness inputs

Figure 4 :

 4 Figure 4: Randomness outputs

Inria

  4

  Figure 6: Randomness consumers -File server

Figure 8 :

 8 Figure 8: Evolution of the input pool entropy counters over time. The top graph corresponds to the desktop workstation scenario. The middle graph corresponds to the file server scenario. The bottom graph corresponds to the computing server scenario.

3

  Random numbers characterization . . . . . . . . . . . . . . . . . . . . . . . . . . 3 2.2 Random Number Generators . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 5 Context: the Linux Random Number Generator 6 3.1 Output interfaces . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 6 3.2 Entropy pools . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 6 3.3 Entropy sources . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 7 3.4 Entropy estimation . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 9 3.5 Alternatives to the LRNG . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 10 4 Contribution 10 4.1 Experimental setup . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 10 4.2 Tested scenarios . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 11 4.3 Entropy inputs . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 11 4.4 Entropy outputs . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 13 4.5 Entropy consuming applications . . . . . . . . . . . . . . . . . . . . . . . . . . . . 14 4.5.1 Desktop Workstation scenario . . . . . . . . . . . . . . . . . . . . . . . . . 15 4.5.2 File server scenario . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 16 4.5.3 Computing server scenario . . . . . . . . . . . . . . . . . . . . . . . . . . . 16 4.5.4 Subsection conclusion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 17 4.6 Entropy counter of the input pool . . . . . . . . . . . . . . . . . . . . . . . . . . 17 4.7 Proportion of events generating entropy . . . . . . . . . . . . . . . . . . . . . . . 19 4.8 Differences between VMs and real machines . . . . . . . . . . . . . . . . . . . . . 19

.5.1 Desktop Workstation scenario

	Desktop Workstation scenario
	Input	Entropy (bits) Entropy (%)
	[K] load_elf_binary	-11896	45.59
	[U] svn	-6728	25.78
	[U] chromium-browse		21.12
	[U] php5	-640	2.45
	[K] inet_frag_secret_rebuild	-320	1.23
	[U] gnome-terminal	-256	0.98
	[U] gconfd-2	-128	0.49
	[U] gconftool-2	-128	0.49
	[U] gdm3	-128	0.49
	[K] nl_pid_hash_rehash	-128	0.49
	[U] python	-104	0.40
	[K] rt_cache_invalidate	-64	0.25
	[K] rt_genid_init	-64	0.25

This test is performed on a Ubuntu 12.04 machine, running the default desktop environment. The processors is an Intel Core i3 (2.27 GHz x 4), and the machine has 3.7 GB of RAM. The version of the Linux kernel is 3.2.0. We measure the completion time for each of the 1000 requests.RR n°8060

RR n°8060

Address Space Layout Randomization (ASLR) is a computer security method which involves randomly arranging the positions of key data areas, usually including the base of the executable and position of libraries, heap, and stack, in a process's address space. cf http://en.wikipedia.org/wiki/Address_space_layout_randomization RR n°8060

http://code.google.com/p/google-authenticator/issues/detail?id=107

https://groups.google.com/forum/?fromgroups#!topic/comp.security.ssh/UG8KDZBAM7g%5B1-25%5D Inria

Proportion of events generating entropy

In this experiment, for each event type, we compare the number of events that actually generate entropy with the total number of event, in order to find the proportion of "useful" events. We summarize the results in the table below.

Desktop workstation scenario

Input

All events (bits)

Secondary result:

The less disk events there are, the more entropy-efficient the events are.

Differences between VMs and real machines

Where as the final experiments have been run on physical computers, a huge part of the development has been led on virtual machines. Consequently, we have noticed the differences between the results on virtual and physical machines and we present them below.

Secondary result:

In every one-hour scenario we played on virtual machines, the number of cycles elapsed was almost the same. Besides, the number of cycles seems to be a maximum, like there is no sleep mode on the CPU. This specificity has an unexpected consequence: the slopes on the input pool level graph tend to be higher, especially for the file server scenario. Surprisingly, it means that the virtualization would have a rather positive effect on the entropy harvesting process.

Secondary result:

The repartition of entropy inputs and outputs are more or less the same.

Only the disk looks a bit more efficient in the case of a virtual machine.

LRNG customers 20

6 Conclusions and perspectives 20


Annette Bieniusa

Marek Zawirski

Nuno Preguiça

Marc Shapiro

HASLab Carlos Baquero

Valter Balegas

Sérgio Duarte

An optimized conflict-free replicated set

Keywords: Data replication, optimistic replication, commutative operations Réplication des données, réplication optimiste, opérations commutatives

à la diffusion de documents scientifiques de niveau recherche, publiés ou non, émanant des établissements d'enseignement et de recherche français ou étrangers, des laboratoires publics ou privés.

Introduction

Eventual consistency of replicated data supports concurrent updates, reduces latency and improves fault tolerance, but forgoes strong consistency (e.g., linearisability). Accordingly, several cloud computing platforms implement eventually-consistent replicated sets [START_REF] Decandia | Dynamo: Amazon's highly available key-value store[END_REF][START_REF] Sovran | Transactional storage for georeplicated systems[END_REF]. Eventual Consistency, allows concurrent updates at different replicas, under the expectation that replicas will eventually converge [START_REF] Vogels | Eventually consistent[END_REF]. However, solutions for addressing concurrent updates tend to be either limited or very complex and error-prone [START_REF] Saito | Optimistic replication[END_REF].

We follow a different approach: Strong Eventual Consistency (SEC) [START_REF] Shapiro | Conflict-free replicated data types[END_REF] requires a deterministic outcome for any pair of concurrent updates. Thus, different replicas can be updated in parallel, and concurrent updates are resolved locally, without requiring consensus. Some simple conditions (e.g., that concurrent updates commute with one another) are sufficient to ensure SEC. Data types that satisfy these conditions are called Conflict-Free Replicated Data Types (CRDTs). Replicas of a CRDT object can be updated without synchronization and are guaranteed to converge. This approach has been adopted in several works [START_REF] Gene | Efficient solutions to the replicated log and dictionary problems[END_REF][START_REF] Thomas | A majority consensus approach to concurrency control for multiple copy databases[END_REF][START_REF] Preguiça | A commutative replicated data type for cooperative editing[END_REF][START_REF] Weiss | Logoot-undo: Distributed collaborative editing system on P2P networks[END_REF][START_REF] Sovran | Transactional storage for georeplicated systems[END_REF].

The set is a pervasive data type, used either directly or as a component of more complex data types, such as maps or graphs. This paper highlights the semantics of sets under eventual consistency, and introduces an optimized set implementation, Optimized Observed Remove Set.

Principle of Permutation Equivalence

The sequential semantics of a set are well known, and are defined by individual updates, e.g., {true}add(e){e ∈ S} (in "{pre-condition} computation {post-condition}" notation), where S denotes its abstract state. However, the semantics of concurrent modifications is left underspecified or implementation-driven.

We propose the following Principle of Permutation Equivalence [START_REF] Bieniusa | Brief annoucement: Semantics of eventually consistent replicated sets[END_REF] to express that concurrent behaviour conforms to the sequential specification: "If all sequential permutations of updates lead to equivalent states, then it should also hold that concurrent executions of the updates lead to equivalent states." It implies the following behavior, for some updates u and u :

{P }u; u {Q} ∧ {P }u ; u{Q } ∧ Q ⇔ Q ⇒ {P }u u {Q}
Specifically for replicated sets, the Principle of Permutation Equivalence requires that {e = f }add(e) remove(f ){e ∈ S∧f / ∈ S}, and similarly for operations on different elements or idempotent operations. Only the pair add(e) remove(e) is unspecified by the principle, since add(e); remove(e) differs from remove(e); add(e). Any of the following post-conditions ensures a deterministic result: where < CLK compares unique clocks associated with the operations. Note that not all concurrency semantics can be explained as a sequential permutation; for instance no sequential execution ever results in an error mark.

!!"# "##$!%# !!$%"# !"#$%&'#()*+ !%"# "##$%%# !!$%"# &!'$!%# "##$(%# !!$%$("# !%"# !"#$%&'#()*+ !!$%$("# !!$%$("# &'()*+,#-# &'()*+,#.# (a) Dynamo shopping cart !""!#"# !"#$%&'#()*+ #$%&%'%(%)% $#%# !""!#"# &'%% *#+!#"# ()*+,-.#/# ()*+,-.#&# #$%&%'%(%)% $#%# #$%&%,%(%)% $#%# #$%&%,%(%)% $#%# #$%&%,%(%,% $%# !"#$%&'#()*+ !""!#"# #$%&%-%(%,% $#%# &'## *#+.#/% !""!#"# #$%&%,%(%,% $%# #$%&%-%(%,% $#%# #$%&%0%(%0% $%# #$%&%0%(%0% $%# (b) C-Set !""!#"!# !""$!#%&"# !"#$%&'#()*+ $#%# !""!#"#!# !""$!#%&'"# (#)!#"#!# (#)$!$!#%&"&!#%&'"%"# '()*+,-#.# '()*+,-#/# *0$+#%&,%# 10$%# $#%# $#%# $#%# $%# !"#$%&'#()*+ !""!#"#!" !""$!#%&'''"# $#%# (#)!#"-!# (#)$!$!#%&"&!#%&'"%"# !""!#"#!- !""$!#%&''"# $%# $#%# $#%# $#%# *0$!#%&'"%# 10$%# *.$!#%&"&!#%&'"%# 10$%# *.$!#%&"&!#%&'"%####### 10$%# *./0-------------------- 10$!#%&"&!#%&'"%# *./0------------------- 10$!#%&"&!#%&'"%# *0$!#%&''"%######### 10$!#%&"&!#%&2"%# *./!#%&''',0------- 10$!#%&"&!#%&2"%# *./!#%&''"&!#%&'''"%- 10$!#%&"&!#%&2"%# *./!#%&''"&!#%&'''"%- 10$!#%&"&!#%&2"%# (c) OR-Set

A Review of Existing Replicated Set Designs

In the past, several designs have been proposed for maintaining a replicated set. Most of them violate the Principle of Permutation Equivalence (Fig. 1). For instance, the Amazon Dynamo shopping cart [START_REF] Decandia | Dynamo: Amazon's highly available key-value store[END_REF] is implemented using a register supporting read and write (assignment) operations, offering the standard sequential semantics. When two writes occur concurrently, the next read returns their union. As noted by the authors themselves, in case of concurrent updates even on unrelated elements, a remove may be undone (Fig. 1(a)). Sovran et al. and Asian et al. [START_REF] Sovran | Transactional storage for georeplicated systems[END_REF][START_REF] Aslan | C-Set: a commutative replicated data type for semantic stores[END_REF] propose a set variant, C-Set, where for each element the associated add and remove updates are counted. The element is in the abstraction if their difference is positive. C-Set violates the Principle of Permutation Equivalence (Fig. 1(b)). When delivering the updates to both replicas as sketched, the add and remove counts are equal, i.e., e is not in the abstraction, even though the last update at each replica is add(e).

Add-wins Replicated Sets

In Section 2 we have shown that when considering concurrent add and remove operations over the same element, one among several post-conditions can be chosen. Considering the case of an add wins semantics we now recall [START_REF] Shapiro | Conflict-free replicated data types[END_REF] the CRDT design of an Observed Remove Set, or OR-Set, and then introduce an optimized design that preserves the OR-Set behaviour and greatly improves its space complexity. These CRDT specifications follow a new notation with mixed state-and operation-based update propagation. Although the formalization of this mixed model, and the associated proof obligations that check compliance to CRDT requisites, is out of the scope of this report the notation is easy to infer from the standard CRDT model [START_REF] Shapiro | Conflict-free replicated data types[END_REF][START_REF] Shapiro | A comprehensive study of Convergent and Commutative Replicated Data Types[END_REF][START_REF] Shapiro | Convergent and commutative replicated data types[END_REF].

E := E \ R T := T ∪ R compare (A, B) : boolean b let b = ((A.E ∪ A.T ) ⊆ (B.E ∪ B.T )) ∧ (A.T ⊆ B.T ) merge (B) E := (E \ B.T ) ∪ (B.E \ T ) T := T ∪ B.T
System model synopsis: We consider a single object, replicated at a given set of processes/replicas. A client of the object may invoke an operation at some replica of its choice, which is called the source of the operation. A query executes entirely at the source. An update applies its side effects first to the source replica, then (eventually) at all replicas, in the downstream for that update. To this effect, an update is modeled as an update pair (p, u) that includes two operations such that p is a side-effect free prepare(-update) operation and u is an effect(-update) operation; the source executes the prepare and effect atomically; downstream replicas execute only the effect u. In the mixed state-and operation-based modelling, replica state can both be changed by applying an effect operation or by merging state from another replica of the same object. The monotonic evolution of replica states is described by a compare operation, supplied with each CRDT specification.

Observed Remove Set

Figure 2 shows our specification for an add-wins replicated set CRDT. Its concurrent specification {P }u 0 . . . u n-1 {Q} is for each element e defined as follows:

• ∀i, u i = remove(e) ⇒ Q = (e / ∈ S) • ∃i : u i = add(e) ⇒ Q = (e ∈ S).
To implement add-wins, the idea is to distinguish different invocations of add(e) by adding a hidden unique token n, and effectively store (e, n) pair. A pair (e, n) is removed by adding it to a tombstone set. An element can be always added again, because the new pair (e, n ) uses always a fresh token, different from the old one, n = n. If the same element e is both added and removed concurrently, the update-prepare of remove concerns only observed pairs (e, n 1 ), (e, n 2 ), . . . and not the concurrently-added unique pair (e, n ). Therefore the add wins by adding a new pair. We call this object an Observed Remove Set, or OR-Set. As illustrated in Figure 1(c), OR-Set is immune from the anomaly that plagues C-Set.

Space complexity: The payload size of OR-Set is at any moment bounded by the number of all applied add (effect-update) operations.

Optimized Observed Remove Set

The OR-Set design uses extensively unique identifiers and tombstones, as other CRDTs [START_REF] Preguiça | A commutative replicated data type for cooperative editing[END_REF][START_REF] Weiss | Logoot-undo: Distributed collaborative editing system on P2P networks[END_REF][START_REF] Shapiro | A comprehensive study of Convergent and Commutative Replicated Data Types[END_REF]. We now show how to make CRDT practical by minimizing the required metadata.

Immediately discarding tombstones: When comparing two payloads P and P , respectively containing some element e and the other not, it is important to know if e has been recently added to P , or if it was recently removed from P . The presented add-wins set uses tombstones to unambiguously answer this question, even when updates are delivered out of order or multiple times.

Tombstones accumulate (as a consequence of the monotonic semilattice requirement); if they cannot be discarded, memory requirements grow with the number of operations. To address this issue, Wuu's 2P-Set [START_REF] Gene | Efficient solutions to the replicated log and dictionary problems[END_REF] garbage-collects tombstones that have been delivered everywhere, basically by waiting for acknowledgements from each process to every other process. This adds communication and processing overhead, and requires all processes to be correct. We devise a novel technique to eliminate tombstones without these limitations and offer conflict-free semantics at an affordable cost. We present our solution using add-wins as the example.

To recapitulate, in OR-Set, adding an element e creates a new unique (e, n) pair to the E part of the payload. Removing an element moves all pairs containing e observed at the source from E to T . 1 Note that adding some pair (e, n) always happens-before removing the same pair (e, n). If updates are delivered only in causal order,once, the add always executes before any related removes, and the tombstone set T is not necessary when executing operations. However, we also need to support state-based merge, which joins two replicas possibly unrelated by happens-before. When merging two replicas in which only 1 A practical implementation will just set a mark bit on the representation of the removed pair and will deallocate any other associated storage. Consider for instance the extension of OR-Set to a map: a key will have some associated value, e.g., E would contain triples (e, n, value). When the key is removed, value can be discarded, but the corresponding (e, n) pair(s) must remain in T . prepare (e) let r = myID()

INRIA

--r = source replica let c = v[r] + 1 effect (e, c, r) pre causal delivery if c > v[r] then let O = {(e, c , r) ∈ E|c < c} v[r] := c E := E ∪ {(e, c, r)} \ O update remove (element e)
prepare (e) --Collect all unique triples containing e let R = {(e, c, i) ∈ E} effect (R) --Remove triples observed at source pre causal delivery one replica has some pair (e, n), we need to know if the pair has been added to the replica that contains it or if it was removed in the other replica. We leverage these observations to propose a novel remove algorithm that discards a removed pair immediately and works safely with merge. It compactly records happensbefore information to summarizes removed elements. Figure 3 presents Optimized OR-Set (Opt-OR-Set) based on this approach.

E := E \ R compare (A, B) : boolean b let R = {(c, i)|0 < c ≤ A.v[i] ∧ ∃e : (e, c, i) ∈ A.E} let R = {(c, i)|0 < c ≤ B.v[i] ∧ ∃e : (e, c, i) ∈ B.E} let b = A.v ≤ B.v ∧ R ⊆ R merge (B) let M = (E ∩ B.E) let M = {(e, c, i) ∈ E \ B.E|c > B.v[i]} let M = {(e, c, i) ∈ B.E \ E|c > v[i]} let U = M ∪ M ∪ M let O = {(e, c, i) ∈ U |∃(e, c , i) ∈ U : c < c } E := U \ O v := [max(v[0], B.v[0]), . . . , max(v[n], B.v[n])]
Each replica i maintains a vector v [START_REF] Parker | Detection of mutual inconsistency in distributed systems[END_REF] to summarize the unique identifiers it has already observed. Entry v[j] = n at replica i indicates that this replica has observed n successive identifiers generated at j: (1, j), (2, j), . . . , (n, j). Replica i maintains its local counter as the i-th entry in the vector v[i], initially 0. A replica generates new unique identifiers (c, i) by incrementing its local counter. Note that to summarize successive identifiers in a vector, OptORSet requires causal delivery of updates. 2When add is invoked, the source associates it with a unique identifier made of the next local counter value and source replica identifier. When the add is delivered to a downstream replica, it should have an effect only if it has not been previously delivered; for this, it checks if the unique identifier is incorporated in the downstream replica's vector. When mergeing payloads, an element should be in the merged state only if: either it is in both payloads (set M in Figure 3), or it is in the local payload and not recently removed from the remote one (set M ) or vice-versa (M ) -an element has been removed if it is not in the payload but its identifier is reflected in the replica's vector.

This approach can be generalized to any CRDT where elements are added and removed, e.g., a sequence [START_REF] Preguiça | A commutative replicated data type for cooperative editing[END_REF][START_REF] Weiss | Logoot-undo: Distributed collaborative editing system on P2P networks[END_REF] or a graph [START_REF] Shapiro | Convergent and commutative replicated data types[END_REF].

Coalescing repeated adds: Another source of memory growth in the original OR-Set is due to the elements added several times. Similarly to tombstones, they pollute the state with unique identifiers for every add. We observe that for every combination of element and source replica, it is enough to keep the identifier of the latest add, which subsumes previously added elements. The OptORSet specification leverages this observation in add and merge definitions, by discarding unnecessary identifiers (set O).

Space complexity: The payload size of OptORSet set is bounded by O(|elements|n+n) at any moment, where n is the number of processes in the systems and |elements| is the number of elements present in the set. The first component corresponds to the maximum number of timestamps in set E and the second captures the size of the vector v. In the common case, where the number of processes repeatedly invoking adds can be considered a constant, the payload size is O(|elements| + n).

Conclusions

Conflict-Free Replicated Data Types (CRDTs) allow a system to maintain multiple replicas of data that are updated without requiring synchronization while guaranteeing Strong Eventual Consistency. This allows, for example, a cloud infrastructure to maintain replicas of data in data centers spread over large geographical distance and still provide low access latency by choosing the closest, to client, data center.

In this paper we reviewed existing replicated set designs and contrasted then with the CRDT OR-Set design, under the principle of permutation equivalence. Having in mind that the base OR-Set favored simplicity at the expense of scalability, we introduced a new optimized design, Optimized OR-Set, that greatly improves its scalability and should favor efficient implementations of sets and other CRDTs that share the OR-Set design techniques.

Figure 1 :

 1 Figure 1: Examples of anomalies and a correct design.

  INRIA payload set E, set T --E: elements; T : tombstones --sets of pairs { (element e, unique-tag n), . . . } initial ∅, ∅ query contains (element e) : boolean b let b = (∃n : (e, n) ∈ E) query elements () : set S let S = {e|∃n : (e, n) ∈ E} update add (element e) prepare (e) let n = unique() -unique() returns a unique tag effect (e, n) E := E ∪ {(e, n)} \ T -e + unique tag update remove (element e) prepare (e) --Collect pairs containing e let R = {(e, n)|∃n : (e, n) ∈ E} effect (R) --Remove pairs observed at source

Figure 2 :

 2 Figure 2: OR-Set: Add-wins replicated set

  payload set E, vect v --E: elements, set of triples (element e, timestamp c, replica i) -v: summary (vector) of received triples initial ∅, [0, . . . , 0] query contains (element e) : boolean b let b = (∃c, i : (e, c, i) ∈ E) query elements () : set S let S = {e|∃c, i : (e, c, i) ∈ E} update add (element e)

Figure 3 :

 3 Figure 3: Optimized OR-Set (Opt-OR-Set).

RR n°8083

It is easy to extend this solution for updates delivered out of happens-before order by using instead a version vector with exceptions[START_REF] Malkhi | Concise version vectors in WinFS[END_REF].INRIA

* This research is supported in part by ANR project ConcoRDanT (ANR-10-BLAN 0208), by ERDF, COMPETE Programme, by Google European Doctoral Fellowship in Distributed Computing received by Marek Zawirski, and FCT projects #PTDC/EIA-EIA/104022/2008 and #PTDC/EIA-EIA/108963/2008. type de données ensemble répliqué sans conflit


E Birmelé

P Crescenzi

R Ferreira

R Grossi

V Lacroix

A Marino

N Pisanti

G Sacomoto

M.-F Sagot

Efficient bubble enumeration in directed graphs

Polymorphisms in DNA-or RNA-seq data lead to recognisable patterns in a de Bruijn graph representation of the reads obtained by sequencing. Such patterns have been called mouths, or bubbles in the literature. They correspond to two vertex-disjoint directed paths between a source s and a target t. Due to the high number of such bubbles that may be present in real data, their enumeration is a major issue concerning the efficiency of dedicated algorithms. We propose in this paper the first linear delay algorithm to enumerate all bubbles with a given source.

Introduction

In recent papers [START_REF] Sacomoto | KisSplice: de-novo calling alternative splicing events from rna-seq data[END_REF][START_REF] Peterlongo | Identifying SNPs without a reference genome by comparing raw reads[END_REF], algorithms for identifying two types of polymorphim, respectively SNPs (Single Nucleotide Polymorphisms) in DNA, and alternative splicing in RNA-seq data were introduced. Both correspond to recognisable patterns in a de Bruijn graph (DBG) built from the reads provided by a sequencing project. In both cases, the pattern corresponds to two vertex-disjoint paths between a pair of source and target vertices s and t. Properties on the lengths or sequence similarity of the paths then enable to differentiate between different types of polymorphism.

Such patterns have been studied before in the context of genome assembly where they have been called bulges [START_REF] Pevzner | De novo repeat classification and fragment assembly[END_REF] or bubbles [START_REF] Robertson | De novo assembly and analysis of RNA-seq data[END_REF][START_REF] Simpson | ABySS: A parallel assembler for short read sequence data[END_REF][START_REF] Zerbino | Velvet: Algorithms for de novo short read assembly using de bruijn graphs[END_REF]. However, the purpose in these works was not to enumerate all these patterns, but "only" to remove them from the graph, in order to provide longer contigs for the genome assembly. More recently, ad-hoc enumeration methods have been proposed but are restricted to non-branching bubbles [START_REF] Iqbal | De novo assembly and genotyping of variants using colored de bruijn graphs[END_REF], i.e., each vertex from the bubble has in-degree and out-degree 1, except for s and t. Furthermore, in all these applications [START_REF] Robertson | De novo assembly and analysis of RNA-seq data[END_REF][START_REF] Simpson | ABySS: A parallel assembler for short read sequence data[END_REF][START_REF] Iqbal | De novo assembly and genotyping of variants using colored de bruijn graphs[END_REF][START_REF] Pevzner | De novo repeat classification and fragment assembly[END_REF][START_REF] Zerbino | Velvet: Algorithms for de novo short read assembly using de bruijn graphs[END_REF] since the patterns correspond to SNPs or sequencing errors, the authors only considered paths of length smaller than a constant. On the other hand, bubbles of arbitrary length have been considered in the context of splicing graphs [START_REF] Sammeth | Complete alternative splicing events are bubbles in splicing graphs[END_REF]. However, in this context, a notable difference is that the graph is a DAG. Additionally, vertices are coloured and only unicolour paths are then considered for forming bubbles. Finally, the concept of bubble also applies to the area of phylogenetic networks [START_REF] Gusfield | Optimal, efficient reconstruction of phylogenetic networks with constrained recombination[END_REF], where it corresponds to the notion of a recombination cycle. Again for this application, the graph is a DAG.

In this paper, we adopt the term bubble, which is being most used in the community, and this will denote two vertex-disjoint paths between a pair of source and target vertices with no condition on the path length or the degrees of the internal nodes. We then consider the more general problem of enumerating all bubbles in a arbitrary directed graph. That is, our solution is not restricted to acyclic or de Bruijn graphs. This problem is quite general but it was still an open question whether a polynomial-delay algorithm could be proposed for solving it. The algorithm presented in [START_REF] Sacomoto | KisSplice: de-novo calling alternative splicing events from rna-seq data[END_REF] was an adaptation of Tiernan's algorithm for cycle enumeration [START_REF] Tiernan | An efficient search algorithm to find the elementary circuits of a graph[END_REF] which does not have a polynomial delay, in the worst case the time elapsed between the output of two solutions is proportional to the number of paths in the graph, i.e. exponential in the size of the graph. It was not clear at the time if more efficient cycle enumeration methods in directed graphs such as Tarjan's [START_REF] Tarjan | Enumeration of the elementary circuits of a directed graph[END_REF] or Johnson's [START_REF] Johnson | Finding all the elementary circuits of a directed graph[END_REF] could be adapted to efficiently enumerate bubbles in directed graphs.

The aim of this paper is to show a non trivial adaptation of Johnson's cycle (what he called elementary circuit) enumeration algorithm to identify all bubbles in a directed graph in the same theoretical complexity. Notably, the method we propose enumerates all bubbles with a given source with O(|V | + |E|) delay, where V , resp. E, is the set of vertices, resp. arcs, of the graph. The algorithm requires an initial transformation, described in Section 3, of the graph for each source s that takes O(|V |+|E|) time and space. Moreover, we briefly describe, in Section 6, a slightly more complex version of the algorithm (but with the same overall complexity) that is more space and time efficient in practice.

De Bruijn graphs and bubbles

A de Bruijn graph (DBG) is a directed graph G = (V, E) whose set of vertices V are labelled by k-mers, i.e. words of length k. An arc in E links a vertex u to a vertex v if the suffix of length k -1 of u is a prefix of v. Given two vertices s and t in G, an (s, t)-path is a path from s to t. By an (s, t)-bubble, we mean two vertex-disjoint (s, t)-paths that only shares s and t.

In the case of next generation sequencing (NGS) data, the k-mers correspond to all words of length k present in the reads (strings) of the input dataset, and only those. In relation to the classical de Bruijn graph for all possible words of size k, the DBG for NGS data may then not be complete. Vertices may also be labelled by the number of times each k-mer is present in the reads. In general a vertex will be labelled by both a k-mer and its reverse complement, and the DBG used in practice will thus be a bi-directed multigraph. Figure 1 gives an example of a portion of a DBG that corresponds to a bubble generated by a SNP or a sequencing error. In this paper, we ignore all details related to the treatment of NGS data using De Bruijn graphs that are not essential for the algorithm described, and consider instead the more general case of finding all (s, t)-bubbles in an arbitrary directed graph.

Turning bubbles into cycles

Let G = (V, E) be a directed graph, and let s ∈ V . We want to find all (s, t)bubbles for all possible target vertices t. We transform G into a new graph

G s = (V s , E s ) where |V s | = 2|V | and |E s | = O(|V | + |E|). Namely, V s = {v, v | v ∈ V } E s = {(u, v), (v, u) | (u, v) ∈ E and v = s} ∪ {(v, v) | v ∈ V and v = s} ∪ {(s, s)}
Let us denote by V the set of vertices of G s that were not already in G, that is V = V s \ V . The two vertices x ∈ V and x ∈ V are said to be twin vertices.

Observe that the graph G s is thus built by adding to G a reversed copy of itself, where the copy of each vertex is referred to as its twin. The arcs incoming to s (and outgoing from s) are not included so that the only cycles in G s that contain s also contain s. New arcs are also created between each pair of twins: the new arcs are the ones leading from a vertex u to its twin ū for all u except for s where the arc goes from s to s. An example of a transformation is given in Figure 2.

We define a cycle of G s as being bipolar if it contains vertices of both V and V . As the only arc from V to V is (s, s), then every bipolar cycle C contains also only one arc from V to V . This arc, which is the arc (t, t) for some t ∈ V , is called the swap arc of C. Moreover, since (s, s) is the only incoming arc of s, all the cycles containing s are bipolar. We say that C is twin-free if it contains no pair of twins except for (s, s) and (t, t).

Definition 1 (Bubble-cycle).

A bubble-cycle in G s is a twin-free cycle of size greater than four1 . Proposition 1. Given a vertex s in G, there is a one-to-two correspondence between the set of (s, t)-bubbles in G for all t ∈ V , and the set of bubble-cycles of G s . Fig. 2: Graph G and its transformation G s . We have that s, e, e, b, a, s, s is a bubble-cycle with swap arc (e, e) that has a correspondence to the (s, e)-bubble composed by the two vertex-disjoint paths s, e and s, a, b, e .

Proof. Let us consider an (s, t)-bubble in G formed by two vertex-disjoint (s, t)paths P and Q. Consider the cycle of G s obtained by concatenating P (resp. Q), the arc (t, t), the inverted copy of Q (resp. P ), and the arc (s, s). Both cycles are bipolar, twin-free, and have (t, t) as swap arc. Therefore both are bubble-cycles.

Conversely, consider any bubble-cycle C and let (t, t) be its swap arc. C is composed by a first subpath P from s to t that traverses vertices of V and a second subpath Q from t to s composed of vertices of V only. By definition of G s , the arcs of the subpath P form a path from s to t in the original graph G; given that the vertices in the subpath Q from t to s are in V and use arcs that are those of E inverted, then Q corresponds to another path from s to t of the original graph G. As no internal vertex of Q is a twin of a vertex in P , these two paths from s to t are vertex-disjoint, and hence they form an (s, t)-bubble.

Notice that there is a cycle s, v, v, s for each v in the out-neighborhood of s. Such cycles do not correspond to any bubble in G, and the condition on the size of C allows us to rule them out.

The algorithm

Johnson [START_REF] Johnson | Finding all the elementary circuits of a directed graph[END_REF] introduced a polynomial delay algorithm for the cycle enumeration problem in directed graphs. We propose to adapt the principle of this algorithm, the prunned backtracking, to enumerate bubble-cycles in G s . Indeed, we use a similiar prunning strategy, modified to take into account the twin nodes. Proposition 1 then ensures that running our algorithm on G s for every s ∈ V is equivalent to the enumeration of (twice) all the bubbles of G. To do so, we explore G s by recursively traversing it while maintaining the following three variables. We denote by N + (v) the set of out-neighbors and N -(v) as the set of in-neighbors of v.

1.

A variable stack which contains the vertices of a path (with no repeated vertices) from s to the current vertex. Each time it is possible to reach s from the current vertex by satisfying all the conditions to have a bubblecycle, this stack is completed into a bubble-cycle and its content output. 2. A variable status(v) for each vertex v which can take three possible values:

f ree: v should be explored during the traversal of G s ; blocked: v should not be explored because it is already in the stack or because it is not possible to complete the current stack into a cycle by going through v -notice that the key idea of the algorithm is that a vertex may be blocked without being on the stack, avoiding thus useless explorations; twinned: v ∈ V and its twin is already in the stack, so that v should not be explored. 3. A set B(v) of in-neighbors of v where vertex v is blocked and for each vertex w ∈ B(v) there exists an arc

(w, v) in G s (that is, w ∈ N -(v)).
If a modification in the stack causes that v is unblocked and it is possible to go from v to s using free vertices, then w should be unblocked if it is currently blocked.

Algorithm 1 enumerates all the bubble-cycles in G by fixing the source s of the (s, t)-bubble, computing the transformed graph G s and then listing all bubblecycles with source s in G s . This procedure is repeated for each vertex s ∈ V . To list the bubble-cycles with source s, procedure CYCLE(s) is called. As a general approach, Algorithm 3 uses classical backtracking with a pruned search tree. The root of the recursion corresponds to the enumeration of all bubble-cycles in G s with starting point s. The algorithm then proceeds recursively: for each free out-neighbor w of v the algorithm enumerates all bubble-cycles that have the vertices in the current stack plus w as a prefix. If v ∈ V and v is twinned, the recursion is also applied to the current stack plus v, (v, v) becoming the current swap arc. A base case of the recursion happens when s is reached and the call to CYCLE(s) completed. In this case, the path in stack is a twin-free cycle and, if this cycle has more than 4 vertices, it is a bubble-cycle to output.

The key idea that enables to make this pruned backtracking efficient is the block-unblock strategy. Observe that when CYCLE(v) is called, v is pushed in the stack and to ensure twin-free extensions, v is blocked and v is twinned if v ∈ V . Later, when backtracking, v is popped from the stack but it is not necessarily marked as free. If there were no twin-free cycles with the vertices in the current stack as a prefix, the vertex v would remain blocked and its status would be set to free only at a later stage. The intuition is that either v is a deadend or there remain vertices in the stack that block all twin-free paths from v to s. In order to manage the status of the vertices, the sets B(w) are used. When a vertex v remains blocked while backtracking, it implies that every out-neighbor w of v has been previously blocked or twinned. To indicate that each out-neighbor w ∈ N + (v) (also, v ∈ N -(w) is an in-neighbor of w) blocks vertex v, we add v to each B(w). When, at a later point in the recursion, a vertex w ∈ N + (v) becomes unblocked, v must also be unblocked as possibly there are now bubble-cycles that include v. Algorithm 2 implements this recursive unblocking strategy.

An important difference between the algorithm introduced here and Johnson's is that we now have three possible states for any vertex, i.e. free, blocked Algorithm 1: Main algorithm and twinned, instead of only the first two. The twinned state is necessary to ensure that the two paths of the bubble share no internal vertex. Whenever v is twinned, it can only be explored from v. On the other hand, a blocked vertex should never be explored. A twin vertex v can be already blocked when the algorithm is exploring v, since it could have been unsuccessfully explored by some other call. In this case, it is necessary to verify the status of v, as it is shown in the graph of Figure 3a. Indeed, consider the algorithm starting from s with (s, a) and (a, b) being the first two arcs visited in the lower part. Later, when the calls CYCLE(c) and CYCLE( b) are made, since ā is twinned, both b and c remain blocked. When the algorithm backtracks to a and explores (a, c), the call CYCLE(c) is made and c is already blocked.

1 for s ∈ V do
Another important difference with respect to Johnson's algorithm is that there is a specific order in which the out-neighborhood of a vertex should be explored. In particular, notice that the order in which Algorithm 3 explores the neighbors of a vertex v is: first the vertices in N + (v) \ {v} and then v. A variant of the algorithm where this order would be reversed, visiting first v and then the vertices in N + (v) \ {v}, would fail to enumerate all the bubbles. Indeed, intuitively a vertex can be blocked because the only way to reach s is through a twinned vertex and when that vertex is untwinned the first one is not unblocked. Indeed, consider the graph in Figure 3b One way to address the problem above would be to modify the algorithm so that every time a vertex v is untwinned, a call to UNBLOCK(v) is made. All the bubble-cycles would be correctly enumerated. However, in this case, it is not hard to find an example where the delay would then no longer be linear. Intuitively, visiting first N + (v) \ {v} and, then v, works because every vertex u that was blocked (during the exploration of N + (v) \ {v}) should remain blocked when the algorithm explores v. Indeed, a bubble would be missed only if there existed a path starting from v, going to s through u and avoiding the twinned vertices. This is not possible if no path from N + (v)\{v} to u could be completed into a bubble-cycle by avoiding the twinned vertices, as we will show later on.

Proof of correctness and complexity analysis

Proof of correctness: Algorithm 3 enumerates all bubbles with source s

Lemma 1. Let v be a vertex of G s such that status(v) = blocked, S the set of vertices currently in the stack, and T the set of vertices whose status is equal to twinned. Then S ∪ T is a (v, s) separator, that is, each path, if any exists, from v to s contains at least one vertex in S ∪ T .

Proof. The result is obvious for the vertices in S ∪T . Let v be a vertex of G s such that status(v) = blocked and v / ∈ S ∪ T . This means that when v was popped for the last time, CYCLE(v) was equal to false since v remained blocked.

Let us prove by induction on k that each path to s of length k from a blocked vertex not in S ∪ T contains at least one vertex in S ∪ T .

We first consider the base case k = 1. Suppose that v is a counter-example for k = 1. This means that there is an arc from v to s (s is an out-neighbor of v). However, in that case the output of CYCLE(v) is true, a contradiction because v would then be unblocked.

Suppose that the result is true for k -1 and, by contradiction, that there exists a blocked vertex v / ∈ S ∪ T and a path (v, w, . . . , s) of length k avoiding S ∪ T . Since (w, . . . , s) is a path of length k -1, we can then assume that w is free. Otherwise, if w were blocked, by induction, the path (w, . . . , s) would contain at least one vertex in S ∪ T , and so would the path (v, w, . . . , s).

Since the call to CYCLE(v) returned false (v remained blocked), either w was already blocked or twinned, or the call to CYCLE(w) made inside CYCLE(v) gave an output equal to false. In any case, after the call to CYCLE(v), w was blocked or twinned and v put in B(w).

The conditional at line 16 of the CYCLE procedure ensures that when untwinned, a vertex immediately becomes blocked. Thus, since w is now free, a call to UNBLOCK(w) was made in any case, yielding a call to UNBLOCK(v). This contradicts the fact that v is blocked.

Theorem 1. The algorithm returns only bubble-cycles. Moreover, each of those cycles is returned exactly once.

Proof. Let us first prove that only bubble-cycles are output. As any call to UN-BLOCK (either inside the procedure CYCLE or inside the procedure UNBLOCK itself) is immediately followed by the popping of the considered vertex, no vertex can appear twice in the stack. Thus, the algorithm returns only cycles. They are trivially bipolar as they have to contain s and s to be output.

Consider now a cycle C output by the algorithm with swap arc (t, t). Let (v, w) in C with v = s and v = t. If v is free when v is put on the stack, then v is twinned before w is put on the stack and cannot be explored until w is popped. If v is blocked when v is put on the stack, then by Lemma 1 it remains blocked at least until v is popped. Thus, v cannot be in C, and consequently the output cycles are twin-free.

So far we have proven that the output produces bubble-cycles. Let us now show that all cycles C = {v 0 = s, v 1 , . . . , v l-1 , v l = s, v 0 } satisfying those conditions are output by the algorithm, and each is output exactly once.

The fact that C is not returned twice is a direct consequence of the fact that the stack is different in all the leaves of a backtracking procedure. To show that C is output, let us prove by induction that the stack is equal to {v 0 , . . . , v i } at some point of the algorithm, for every 0 ≤ i ≤ l -1. Indeed, it is true for i = 0. Moreover, suppose that at some point, the stack is {v 0 , . . . , v i-1 }.

Suppose that v i-1 is different from t. As the cycle contains no pair of twins except for those composing the arcs (s, s) and (t, t), the path {v i , v i+1 , . . . , v l } contains no twin of {v 0 , . . . , v i-1 } and therefore no twinned vertex. Thus, it is a path from v i to s avoiding S ∪ T . Lemma 1 then ensures that at this point v i is not blocked. As it is also not twinned, its status is free. Therefore, it will be explored by the backtracking procedure and the stack at some point will be {v 0 , . . . , v i }. If v i-1 = t, v i = t is not blocked using the same arguments. Thus it was twinned by the call to CYCLE(t) and is therefore explored at Line 17 of this procedure. Again, the stack at some point will be {v 0 , . . . , v i }.

Analysis of complexity: Algorithm 3 has linear delay

As in [START_REF] Johnson | Finding all the elementary circuits of a directed graph[END_REF], we show that Algorithm 3 has delay O(|V | + |E|) by proving that a cycle has to be output between two successive unblockings of the same vertex and that with linear delay some vertex has to be unblocked again. To do so, let us first prove the following lemmas.

Lemma 2. Let v be a vertex such that CYCLE(v) returns true. Then a cycle is output after that call and before any call to UNBLOCK.

Proof. Let y be the first vertex such that UNBLOCK(y) is called inside CYCLE(v). Since CYCLE(v) returns true, there is a call to UNBLOCK(v) before it returns, so that y exists. Certainly, UNBLOCK(y) was called before UNBLOCK(v) if y = v. Moreover, the call UNBLOCK(y) was done inside CYCLE(y), from line 37, otherwise it would contradict the choice of y. So, the call to CYCLE(y) was done within the recursive calls inside the call to CYCLE(v). CYCLE(y) must then return true as y was unblocked from it.

All the recursive calls CYCLE(z) made inside CYCLE(y) must return false, otherwise there would be a call to UNBLOCK(z) before UNBLOCK(y), contradicting the choice of y. Since CYCLE(y) must return true and the calls to all the neighbors returned false, the only possibility is that s ∈ N + (y). Therefore, a cycle is output before UNBLOCK(y). Lemma 3. Let v be a vertex such that there is a (v, s)-path P avoiding S ∪ T at the moment a call to CYCLE(v) is made. Then the return value of CYCLE(v) is true.

Proof. First notice that if there is such a path P , then v belongs to a cycle in G s . This cycle may however not be a bubble-cycle in the sense that it may not be twin-free, that is, it may contain more than two pairs of twin vertices. Indeed, since the only constraint that we have on P is that it avoids all vertices that are in S and T when v is reached, then if v ∈ V , it could be that the path P from v to s contains, besides s and s, at least two more pairs of twin vertices. An example is given in Figure 2b. It is however always possible, by construction of G s from G, to find a vertex y ∈ V such that y is the first vertex in P with y also in P . Let P be the path that is a concatenation of the subpath s y of P , the arc (y, y), and the subpath y s in P . This path is twin-free, and a call to CYCLE(v) will, by correctness of the algorithm, return true. Theorem 2. Algorithm 3 has linear delay.

Proof. Let us first prove that between two successive unblockings of any vertex v, a cycle is output. Let w be the vertex such that a call to UNBLOCK(w) at line 37 of Algorithm 3 unblocks v for the first time. Let S and T be, respectively, the current sets of stack and twinned vertices after popping w. The recursive structure of the unblocking procedure then ensures that there exists a (v, w)-path avoiding S ∪ T . Moreover, as the call to UNBLOCK(w) was made at line 37, the answer to CYCLE(w) is true so there exists also a (w, s)-path avoiding S ∪ T . The concatenation of both paths is a again a (v, s)-path avoiding S ∪ T . Let x be the first vertex of this path to be visited again. Note that, if no vertex in this path is visited again there is nothing to prove, since v is free, CYCLE(v) needs to be called before any UNBLOCK(v) call. When CYCLE(x) is called, there is a (x, s)-path avoiding the current S ∪ T . vertices. Thus, applying Lemma 3 and then Lemma 2, we know that a cycle is output before any call to UNBLOCK. As no call to UNBLOCK(v) can be made before the call to CYCLE(x), a cycle is output before the second call to UNBLOCK(v).

Let us now consider the delay of the algorithm. In both its exploration and unblocking phases, the algorithm follows the arcs of the graph and transforms the status or the B lists of their endpoints, which overall require constant time. Thus, the delay only depends on the number of arcs which are considered during two successive outputs. An arc (u, v) is considered once by the algorithm in the three following situations: the exploration part of a call to CYCLE(u); an insertion of u in B(v); a call to UNBLOCK(v). As shown before, UNBLOCK(v) is called only once between two successive outputs. CYCLE(u) cannot be called more than twice. Thus the arc (u, v) is considered at most 5 times between two outputs. This ensures that the delay of the algorithm is O(|V | + |E|).

Practical Speedup

Speeding up preproccessing. In Section 3, the bubble enumeration problem was reduced to the enumeration of some particular cycles in the transformed graph G s for each s. It is worth observing that this does not imply building from scratch G s for each s. Indeed, notice that for any two vertices s 1 and s 2 , we can transform G s1 into G s2 by: (a) removing from G s1 the arcs (s 1 , s 1 ), (s 2 , s 2 ), (v, s 2 ), and (s 2 , v) for each v ∈ N -(s 2 ) in G; (b) adding to G s1 the arcs (s 1 , s 1 ), (s 2 , s 2 ), (v, s 1 ), and (s 1 , v) for each v ∈ N -(s 1 ) in G.

Avoiding duplicate bubbles. The one-to-two correspondence between cycles in G s and bubbles starting from s in G, claimed by Proposition 1, can be reduced to a one-to-one correspondence in the following way. Consider an arbitrary order on the vertices of V , and assign to each vertex of V the order of its twin. Let C be a cycle of G s that passes through s and contains exactly two pairs of twin vertices. Denote again by t the vertex such that (t, t) is the arc through which C swaps from V to V . Denote by swap predecessor the vertex before t in C and by swap successor the vertex after t in C. Proposition 2. There is a one-to-one correspondence between the set of (s, t)bubbles in G for all t ∈ V , and the set of cycles of G s that pass through s, contain exactly two pairs of twin vertices and such that the swap predecessor is greater than the swap successor.

Proof. The proof follows the one of Proposition 1. The only difference is that, if we consider a bubble composed of the paths P 1 and P 2 , one of these two paths, say P 1 , has a next to last vertex greater than the next to last vertex of P 2 . Then the cycle of G s made of P 1 and P 2 is still considered by the algorithm whereas the cycle made of P 2 and P 1 is not. Moreover, the cycles of length four which are of the type {s, t, t, s} are ruled out as s is of the same order as s.

Conclusion

We showed in this paper that it is possible (Algorithm 3) to enumerate all bubbles with a given source in a directed graph with linear delay. Moreover, it is possible to enumerate all bubbles, for all possible sources (Algorithm 1), in O((|E| + |V |)(|C| + |V |)) total time, where |C| is the number of bubbles. This required a non trivial adaptation of Johnson's algorithm [START_REF] Johnson | Finding all the elementary circuits of a directed graph[END_REF].

Fig. 1 :

 1 Fig.1: Bubble due to a substitution (gray letter).

2 stack:=∅; 3 for v ∈ G s do 4 statusAlgorithm 2 : 6 UNBLOCK

 23426 Procedure U N BLOCK(v) 1 /* recursive unblocking of vertices for which popping v creates a path to s */ 2 status(v):=free; 3 for w ∈ B(v) do 4 delete w from B(v); 5 if status(w)==blocked then

Algorithm 3 :

 3 and the twin-first variant starting in s with (s, a) and (a, b) being the first two arcs explored in the lower part of the graph. When the algorithm starts exploring b the stack contains s, a, b . After, Procedure CYCLE(v)

Fig. 3 :

 3 Fig. 3: (a) Example where the twin v is already blocked when the algorithm starts exploring v. By starting in s and visiting first (s, a) and (a, b), the vertex c is already blocked when the algorithm starts exploring c. (b) Counterexample for the variant of the algorithm visiting first the twin and then the regular neighbors. By starting in s and visiting first (s, a) and (a, b), the algorithm misses the bubble-cycle s, a, c, c, b, s .

The only twin-free cycles in of size four in G s are generated by the outgoing edges of s. There are O(|V |) of such cycles.

This work was supported by the french ANR MIRI BLAN08-1335497 Project and the ERC Advanced Grant Sisyphe held by Marie-France Sagot. Partially supported by Italian project PRIN AlgoDEEP (2008TFBWL4) of MIUR. The second author received additional support from the Italian PRIN project 'DISCO'.