Mickaël Hoerdt
email: hoerdt@clarinet.u-strasbg.fr

Frédéric Beck
email: beck@clarinet.u-strasbg.fr

Damien Magoni

Jean-Jacques Pansiot
email: pansiot@crc.u-strasbg.fr

A Source Discovery Protocol for ASM Applications in SSM Networks

The Single Source Multicast (SSM) model has received considerable attention from the research community as it could finally bring a scalable solution for some multicast applications. In the SSM model it is the responsibility of the applications to tell the network what are their sources of interest. As a consequence, the SSM model is well suited for one to many multicast applications (such as TV and radio streaming) because there is only one source to discover and it is supposed to be well known. However there still does not exist an SSM session layer protocol and middleware in order to help Any Source Multicast (ASM) applications to discover the presence or departure of the SSM sources. In this paper, we propose a detailed architecture for such a source discovery session layer protocol. We have implemented our protocol as a middleware called Libemu and freely available as an open source library.

I. INTRODUCTION

Classical multicast multi-sources applications use the Deering's [START_REF] Deering | Multicast routing in a datagram internetwork[END_REF] ASM group communication model where anyone, including non-members, can send packets with one local operation to a group identifier G with best effort guarantees that its packets will be transmitted to all group members: the network duplicates packets and manages the dynamics of sources and receivers. With the new Holbrook's [START_REF] Holbrook | A channel model for multicast[END_REF] SSM model it is the responsibility of applications to tell the network who are their sources of interest. In this paper we propose a detailed architecture and implementation of a source discovery session layer protocol for making ASM applications able to use SSM sources. Among the most important sections of this paper, section IV presents an overview of the problem and our solution to it. In section V we describe the ASM over SSM emulation part of the protocol. Section VI presents an optional source management mechanism to perform some control operations on an applicative group. Finally section VII shows the various states in which a source can be.

II. PREVIOUS WORK

In his thesis [START_REF] Holbrook | A channel model for multicast[END_REF], Holbrook considered the case of using multi-source applications in SSM network only environment. He defined two schemes, an hybrid approach of them and evaluated them: he proposed a sender advertisement scheme using a control channel and a session relaying scheme using a central node, and argues that both of them have advantages and drawbacks. He measured the startup delay of the sessions and simulated both propositions but did not define an architecture for it. That's why we propose a detailed architecture based on his idea in this paper.

In paper [START_REF] Chesterfield | An extensible rtcp control framework for large multimedia distributions[END_REF] Chesterfield and Schooler propose to modify RTP/RTCP to have it working under SSM environment. It is aimed at large scale multimedia distribution and control and dot not consider the case of medium sized multiparty conferences. Their solution could work under specific applications [START_REF] Mccanne | vic:Afl e xible framework for packet video[END_REF], [START_REF] Sasse | Successful multiparty audio communication over the internet[END_REF] but reduces the information provided by RTP/RTCP protocol for conference applications and could only work with multicast applications using RTP/RTCP. We try here to define an architecture aimed at multi-source applications, as independent as possible from the application and only providing independent services for an applicative group.

An architecture based on session relaying has been proposed in [START_REF] Zappala | Using ssm proxies to provide efficient multiple-source multicast delivery[END_REF]. They propose to put multiple proxies on the network to improve the deployment of SSM. They do not consider the problem of multi-source in SSM networks.

A more recent work proposes to modify the SSM model to allow multiple source to emit toward a channel [START_REF] Sarac | Ssm extensions: Network layer support for multiple senders in ssm[END_REF]. They argue that this require little modifications in the forwarding state entries, improve the multicast state scalability, and eliminate the single point of failure problem present in the solutions proposed in [START_REF] Holbrook | A channel model for multicast[END_REF]. We believe that their solution will require a very high deployment effort as it require to change the IGMP and PIM-SSM specification, in the hosts and in the routers. This paper describes a protocol which could be used in the current SSM architecture to achieve the source discovery in the existing SSM Network. We have implemented the protocol as a freely available library for the research community [START_REF] Beck | Libemu 1.0[END_REF].

III. MOTIVATION AND GOALS

As a first point, we would like to have a quickly usable architecture, according to the latest network and application research efforts. We would like to be able to use multisource multicast at the inter-domain Internet scale too. This means using SSM service, because it is simple to maintain for operators and because it is believed to be a viable solution for large scale Internet multicast.

Second, this proposition must be scalable regarding the number of announced sources and fault tolerant regarding the source dynamics. This means using soft-state retransmission and UDP for control messages.

For application performances and ease of portability, the architecture must be as transparent as possible with respect to the classical ASM model.

On the implementation viewpoint, this architecture must be portable, independent of the IP version stack used and easy to use for an application developer.

IV. PROBLEM OVERVIEW

In this section, we will describe an architecture that solves the problems cited, and permits the use of channels in multisender applications. This architecture is inspired by Hugh W. Holbrook's thesis [START_REF] Holbrook | A channel model for multicast[END_REF].

A. Problems to solve

We aim to use channels in multi-sender applications, but to do so, we have several problems to solve.

Firstly, in the ASM model described by Deering, a multicast group is only identified by its IP address. But here, we will work with SSM multicast channels, which are described by a tuple (S,G) composed of the unicast IP address of a source S and an IP multicast address G where only S is allowed to put data in. But then, knowing that at the moment, in most of the multi-source multicast applications a multicast group is identified by both its IP multicast address and the port used to send data to, the problem of identifying a multicast group appears.

Secondly, on the one hand, in ASM when a receiver joins a group, he discovers automatically the sources when they send data to the group. In the same way, when a new source appears, the receivers learn it simply by receiving the packets sent. On the other hand, with the SSM model, no such mechanism exists, because we work with channels (S,G) where S is the single source being able to emit on this channel.

Consequently, the receivers must have a way of discovering the active sources and being informed of the arrival of new ones. To achieve these goals, we have two choices as exposed in [START_REF] Holbrook | A channel model for multicast[END_REF]:

• Session Relaying.

• Sender Advertising. We will now explain a solution, based on sender advertising.

B. Solution proposed

This solution uses at least a channel per source, and a control channel, used for relaying various administration information and source advertising for the group. This control channel identifies the multicast group, and runs on top of UDP for sending data. Therefore, to announce a group, it is only necessary to advertise this channel and the UDP port used, receivers having only to join the multicast channel and bind it to obtain all information necessary to take part in the session. The source of this control channel will be the controller of the group, and will perform various operations for this group, like advertising new sources, giving the list of the existing ones on demand and optionally making some network group management.

As it is possible for the same controller to have authority for several groups, it has been decided to use a single controller that takes care of multiple groups. This controller uses a single UDP socket listening on an announced port to communicate with the sources or the receivers which want to learn the existing sources. The source advertising messages and other administration messages are sent to all the receivers via the control channel. One control channel corresponds to one ASM group at the network level. This controller takes care on a controller S of all the control channels identified by a channel (S,...). These channels do not necessarily use the same UDP port to sent data to.

This controller is run at the application level of the OSI model, because it takes care of the interaction between the network and the applications used in multicast sessions. The advantage of having this controller running on the application layer, is that it is possible to choose its location so that the performances are the best, and it introduces the possibility of doing group management at the application layer. Moreover to manipulate UDP socket is very easy and allows quick implementation on the hosts. This paper describes a protocol which could be used between the sources or the receivers and the controller in SSM networks to manage sources.

V. ASM OVER SSM EMULATION

In this section, we describe the ASM over SSM emulation part of the protocol, which means the way the receivers learn the active sources for a group at a given time, and the new ones which could appear during the session.

A. Protocol description for source-controller interaction

Now we will describe the interaction between a source for a group, and the controller (that means the controller responsible of source announcing ) for this group. We'll only take a look at the source advertising, the group management control is optional and is described in Section VI.

1) With a new source: When a new source wants to announce itself, it must begin by informing the controller. The controller then stores the information about the source in its cache. To achieve robustness and scalability, this announcement is sent periodically with a message composed of the protocol header with an optional SDP payload [START_REF] Handley | Session announcement protocol[END_REF], [START_REF] Handley | Sdp: Session description protocol[END_REF]. We call this message ON. It is sent to the UDP socket of the controller's process. This message must contain all the necessary information for a receiver to join the channel used by the source to send its data to the group. It is composed of the new source channel information and possibly transport or protocol over IP information.

It must be noted that the transport information is only destined to applications and not the protocol described here. In the previous example, the application behavior is the same than the ASM one. It means that UDP sources have to send to the same UDP announced port and that each receiver have to bind on the same UDP port for data reception.

Once this message ON has been treated, the controller announces the source to the receivers and optionally send an ON RESP to the new source to inform it about various session parameters, depending on the group management policy. To do so, it forwards this message to the whole group via the control channel, if the source is allowed to send data to the group. Therefore, all the receivers learn that this source exists, and are able to join the channel it uses to send data if desired. Figure 1 illustrates this mechanism. 2) With an existing source: Periodically, a source must send a message ON to the controller, in order to inform it that it is still active, and the controller forwards this message on its control channel to the receivers. For each source, the controller keeps a timer ON validity which gives the time the source entry in his cache is valid. If the controller does not receive such a message before the timer associated with the source ends, it considers that the source is not active anymore, removes the corresponding entry in his cache and sends a message OF F to inform the receivers that this source no longer exists.

The sources and the receivers also have to maintain timers in order to take care of the validity of the various information and messages. For example, a source has to run a timer for the validity of a message ON in order to take care of the periodical sent of this message.

A source can also explicitly announce that it stops its activity by sending a message OF F to the controller, which removes the corresponding entry in its cache and forwards this message on the corresponding control channel in order to inform the receivers, which then leave this channel. The mechanism is illustrated on figure 2. Note that some upper layer protocols already have a function to advertise the end in a session participation(BYE packets from RTCP protocol). this mechanism may be redundant for some applications but one important thing of our approach is that the control plane is over UDP and the data plane is over IP. This means that any multicast applications working over IP could use this source discovery protocol.

The message OF F is almost identical to the message ON. It is sent to the announced port of the control channel on an UDP socket and is a message composed of this protocol's header with payload containing the identifier of the group and the identifier of the channel used to send data. These two parameters are sufficient to leave the channel corresponding to the source. If the message OF F is lost, the source timer will expire in the controller cache because it is no more updated with the periodic ON messages. As a consequence a OF F message will be sent in the control channel and the receivers will unsubscribe for this source.

B. Protocol description for receiver-controller interaction

Here we describe the interaction between a receiver and the group's controller. We distinguish the case when a receiver demands information about the sources which are active for the group, and the case when a source has announced itself.

1) Information Request: The controller acts here as a directory of sources. A receiver which wants to learn the list of the existing sources for a group on which the controller has authority, sends him a unicast message on a known port from the receiver. This message is an UDP packet with a payload containing the protocol header and the group identifier, which is the tuple (S,G). We call this message a message INFO REQ .

When the controller receives such a message, it takes a look at its cache, and responds to the receiver in unicast with a message INFO RESP giving this list of sources and enough information for the receiver to join the corresponding channels if wanted.

If the controller wants to be able to give these information, it has to maintain a cache containing all these. That means that the cache must contain the group identifier (S,G), if there are more than one controller, the control channels used by these controllers, the dates of beginning and ending of validity of the session, and the source list. That gives the table I.

The source list is composed of zero or several entries. Each entry describes a source and gives the information that a The flags KICKED and MUTED are used if source management operations are performed, and are set to 0 by default, and the coding information are optional.

This cache is used by the controller to store all the information for a group, but the receivers have to maintain the same cache, but limited to the joined applicative group. This cache is updated thanks to the messages sent on the control channel.

The message INFO RESP is an UDP packet with as payload this protocol's header with a SDP description describing the session and the sources. The SDP part contains the group identifier, the dates of beginning and end of validity of the session, and a list of zero or more channels and ports used by the sources to send data to the group, and possibly more information about these sources which the receivers could need.

When a receiver receives a packet INFO RESP ,i fi th a s not done so yet, it joins the control channel, and does the same with the channels used by the sources which it has interest in. See Figure 3. 2) Source advertising: When a new source announces its activity to the controller. The controller has to inform all the receivers. To do so, it forwards the message ON received from the source on the control channel.

In the same way, when the controller detects the explicit or implicit departure of a source, it forwards or sends, according to the case, a message OF F for this source on the control channel.

VI. SOURCE MANAGEMENT CAPABILITIES

This protocol allows the possibility to perform some control operations on an applicative group. These source management operations are optional. In this section, we give two examples of such operations, the KICK and MUTE mechanisms, inspired by the control operations available in the IRC protocol [START_REF] Oikarinen | Internet relay chat protocol[END_REF].

A. KICK

The KICK operation consists in forbidding a source to send data to the group. As we are using channels, the controller cannot block the traffic from the source, but it has the possibility to advice the receivers to not receive it. Indeed, if a source is not wanted anymore for a group, the flag KICKED is set in the entry corresponding from the controller's cache, and a message KICK is sent on the control channel for this source. Then, all the receivers learn that the source is not desired anymore and have the possibility to leave the corresponding channel, after having set the flag KICKED in their cache.

As the source is still active, it sends periodically ON messages to the controller, but as long as the flag KICKED is set for this source, the controller does not forward these messages, but instead sends a KICK message.

If the source is authorized to emit again, it is enough to withdraw the KICKED flag and to propagate a message ON for this source on the control channel, the receivers unsetting the KICKED flag in their cache will have the possibility to join again the corresponding channel if it has been left. This mechanism is transparent for the sources: this operation is done by the receivers which are supposed to leave the channel used by the source.

B. MUTE

This mechanism is less violent than the KICK one. Here the source is not undesirable, but in order to save resources or to maintain session integrity, it is wished that only specified sources transmit data at a given time.

To do so, the controller sends a message MUTE containing the group identifier to the sources that should not send data. When they receive such a message, sent in unicast to their announced port UDP (the address is known from the channel used to send data to), the sources stop sending data on their channels, but keep doing the other things like sending periodical ON messages to the controller.

To authorize the source to emit again, a message UNMUTE containing the group identifier is send to the source.

If a source has been asked to mute, the flag MUTED is set in the corresponding entry in the controller's cache, and this flag is removed if the source is authorized to send data to the group. As long as the flag MUTE is set, the controller periodically sends MUTE message to the source. This operation is totally transparent for the receivers, unlike the KICK one where they have to perform JOIN and LEAVE actions. But this operation needs the cooperation of the source which has to stop sending data to the group. That implies that sources also run a process which listens on this port and performs the needed operations if MUTE or UNMUTE messages are received.

C. Timers

To perform all these operations, both the source and the controller have to keep running timers, in order to take care of the validity of KICK and MUTE messages.

If the timer associated to the MUTE message expires, the source decides that it can send data again, removes the MUTE flag and sends again on its data channel.

A receiver has to run a timer giving the validity of a KICK message. If this timer runs out, the KICKED flag is removed and the channel data corresponding to the source is joined (if it was left).

In the same way, the controller has to maintain two timers giving the duration before the next sent of a KICK or MUTE message. The timers used for KICK and MUTE periodical messages have the same value than the timer used for the ON message.

VII. STATES OF THE SOURCE

In this section, we describe the different states in which a source can be, depending on the action realized thanks to this protocol. We will distinguish the ASM over SSM emulation part and the Control operations.

A. ASM over SSM emulation states

The different in which a source can be in this case are:

• ACTIVE: the source is considered to be active, notifying it periodically to the controller with ON messages.

• INACTIVE: the source has never been active or has sent a OF F message. The following state diagram illustrated in figure 4 describes the different states and the transition between them for the ASM over SSM emulation part, which means the source advertising.

B. Source management states

If the optional source management operations are performed, some new states are introduced. The different states in which a source can be in this case are:

• ACTIVE: the source is normally active.

• MUTED: The source has been asked by the controller to stop sending data on its channel.

• KICKED: The source has been KICKED by the controller (but it is still active).

• INACTIVE: The source has announced that it is not active anymore. The MUTE state is transparent for the receivers, the muted source's channel (multicast tree) is kept. This state is only available for sources. In the same way, the KICKED state has only sense for the receivers, the channel of a kicked source is removed (actually it is virtually removed by the controller's KICK message, the source keeps its data channel and keeps on sending data on it). We can imagine that we send a message KICK to the source to advertise it that it has been kicked from the session, but it is not an obligation as the receivers will not get the data anymore. The INACTIVE state is reached when the receiver has received the OF F message forwarded on the control channel.

The following state diagram illustrated in figure 5 describes the different states and the transition between them for the source management operations.

VIII. BREAKDOWNS AND CONSEQUENCES

The basic idea to prevent breakdowns is setup a backup controller. The master controller S sends on its control channel periodical ALIVE messages to show that it is always active. This message is a UDP packet with this protocol's header and a SDP description containing the group identifier for which it is the controller.

If the backup controller S' does not receive any ALIVE message from S on (S,G) after the configured time, it considers that S is down and begins sending ALIVE packets and source advertising. A source which was kicked by the other controller will be announced as active the next time it sends a message ON unless the new one kicks it. In the same way, a source that has been asked to mute, after the timer of the validity of the MUTE message has expired, will begin to send data again. To avoid this situation, it is possible to distribute this by adding a KICK flag in each receiver's cache. Then, if the master controller kicks a source, it announces it with a message KICK, which causes this KICK to be set in all receivers, including the second master. This flag would be valid until a ON message for this source is sent. Such a mechanism gives us at the application level the possibility for the application based on this protocol to ask the user if he still wants to receive data from the source even if the controller recommends to kick it, or if he wants to follow the controller's order in case of a MUTE message.

Supposing that the master controller S breaks down, S' detects it and begins assuming this role as already described earlier. A receiver that has already joined the group is a member of both control channels, and then receives the control data on (S',G'), and except the time until S' takes the control, nothing special happens. In the same way, a receiver that wants to join the session joins both control channels, and has the information he needs from S'.

IX. TIMERS AND THEIR DEFAULT VALUES

The messages ON, KICK and MUTE are sent periodically in order to add robustness to the protocol. Therefore, we have to use several timers to define the validity of the different states and information.

In this section, we describe the different timers used by this protocol and their default value. Some of them are configurable, the other ones depend on each others.

• ON Validity : This timer describe the validity of a source entry in the controller's cache. It is the reference used to define most of other timers. If we allow its value to being configurable, we could distribute it to the sources and the receivers thanks to the ON message. Default: 15 seconds.

• ON F orwarding : This timer describes the period the controller uses to send ON messages on the control channel. Default: 1/3 * ON Validity (5 seconds).

• ON Cycle : The ON Cycle is the period between two ON messages sent by a source. In order to be robust to one message ON lost, we must have 2 * ON Cycle > ON Validity Default: 1/3 ON Validity (5 seconds).

• ON S rc: The ON Src is the time of validity of a ON message for a receiver. Default: ON Validity (15 seconds).

• KICK: It is the time of validity of a KICK message sent by the controller, and is used at the receiver's level. Default: ON Validity (15 seconds).

• MUTE: It is the time of validity of a MUTE message sent by the controller, and is used at the source's level. Default: ON Validity (15 seconds).

• ALIV E: It is the time of validity of an ALIVE message sent by the controller, and is used at the receiver's level.

Actually, it only makes sense for the other controllers if we have redundancy in order to add robustness. Default: 2*ON Validity (30 seconds).

• ALIV E Cycle : It is the time between two ALIVE messages sent by the controller. Default: 3/4 * ON Validity (12 seconds).

X. SCALABILITY CONSIDERATIONS

An important question is the scalability of this architecture, especially at the controller point. The controller is acting as a sort of rendez-vous point ala PIM, and a relay node for control messages between sources and receivers in the emulated ASM group. It is then a central point of failure, sensible to the bandwidth generated by the different announcing sources.

We consider here the problem of bandwidth scalability at the controller in the worst case with the default timers values. What is the amount of required download bandwidth at the controller point in function of the number of synchronized announcing sources?

If we take the messages defined in [START_REF] Beck | Source discovery protocol in ssm network[END_REF], an ON message which contains one channel will take 158 bytes on the wire, considering that the source is transmitting on Ethernet and IPv6 and that the announced channel is IPv6. With IPv4 it will take 90 bytes considering that the announced channel is IPv4.

In a critical case, if there are N synchronized sources announcing themselves and located at 100 different hosts and that all the messages arrive together at the controller, the mean bandwidth required is a simple linear function, N*158 byte/sfor IPv6 and N*90 byte/s for IPv4 as illustrated in figure 6. The generated bandwidth is not too high as we are considering this case and that applications with more that 100 applicative sources at the same time are very specific and will be certainly based on the ASM network service model if it is finally deployed. Moreover thanks to the acknowledgement mechanism, it is possible to rate limit the sources at the controller. Note that once the messages arrive at the controller, the Controller has complete control over the bandwidth usage of its control channel.

XI. SECURITY CONSIDERATIONS

At this moment, this protocol does not provide security mechanism at all. The only security warrants are based on IP. Because this protocol is based on top of a connectionless protocol (UDP), it's very easy for someone to abuse the security mechanisms offered by IP with techniques like IP spoofing. This is why it is very important on the receiver part to bind on the two part of the (S,G) control channel. According to the SSM service model, and how the SSM multicast tree construction is done at the IP level, it is more difficult to spoof a SSM channel than an IP.

More security mechanisms could be imagined, like using a single session identifier per group which is only distributed between the controllers.

But the security problem is not limited to the controller. We can imagine that a source usurps another one's identity and carries out false advertisements for this one. By doing so, a source can be announced as no more active whereas it still sends data. To avoid this, we could use an appropriate authentication mechanism, which can be done thanks to the acknowledgment messages for example.

XII. IMPLEMENTATION

Our architecture is further defined in an Internet draft available to the research community [START_REF] Beck | Source discovery protocol in ssm network[END_REF]. In particular, this draft contains the full description of the protocol messages, the format of the packets as well as the definition of all headers and flags. We have implemented our middleware as a library called Libemu that can be used by applications running over Windows, Linux or *BSD operating systems depending on the IP stack version. The library is implemented in C++ and consists in more than 8000 lines of code. It is freely available for download at [START_REF] Beck | Libemu 1.0[END_REF]. It's worth noticing that our library makes use of the GNU Common C++ library (freely available) and needs an host implementation of the multicast protocols IGMPv3 and MLDv2. These protocol implementations are provided by the Kame's project under *BSD and can be found in the latest GNU/Linux Kernel version 2.4.22 or in the new 2.6.x versions.

XIII. CONCLUSION

Our work aims at providing a multicast library for multisource applications wanting to connect to Source Specific Multicast only networks. Source Specific Multicast is a new multicast over IP distribution paradigm: by pushing the source selection in the hosts, it is a good candidate for an Internet wide multicast service offer. It is well fitted for TV and radio streaming applications but it does not provide automatic source discovery as the original Any Source Multicast paradigm from Steeve Deering does. Considering multi-source applications with source dynamics and old applications based on the ASM model, a multicast library (usable as a middleware) implementing a source discovery protocol was missing. In this paper we have filled this gap by designing an architecture able to provide such a discovery service. We have defined a new protocol and implemented its corresponding middleware in order to enable ASM applications linked against our library to discover SSM sources. Our next task consists in assessing the performance of our proposal through extensive simulations. We also plan to evaluate our implementation in the environment of the 6NET project and to define a proxy architecture based on this protocol in order to use unmodified ASM applications on SSM networks.

Fig. 1 .

 1 Fig. 1. Interaction Source-Controller for a message ON

Fig. 2 .

 2 Fig. 2. Interaction Source-Controller for a message OF F

Fig. 3 .

 3 Fig. 3. Interaction Receiver-Controller

Fig. 4 .

 4 Fig. 4. State diagram for ASM over SSM emulation

Fig. 5 .

 5 Fig. 5. State diagram for source management

Fig. 6 .

 6 Fig. 6. Required download bandwidth vs number of announcing sources.

  to join the channel used for sending data to the group. Such an entry is illustrated on table II.

			TABLE I
		CONTROLLER'S CACHE INFORMATION
	control channel	second control channel	other control channels	beginning	ending	sources
	receiver may need

XIV. ACKNOWLEDGMENTS

This work was partly sponsored by the 6net European project under contract number IST-2001-32603.

Active


The temperature spectrum generated by frictional heating in isotropic turbulence

In every turbulent flow with non-zero viscosity, heat is generated by viscous friction. This heat is then mixed by the velocity field. We consider how heat fluctuations generated this way are injected and distributed over length scales in isotropic turbulence. A triadic closure is derived and numerically integrated. It is shown how the heat fluctuation spectrum depends on the Reynolds and Prandtl numbers.

Introduction

Turbulence allows for an efficient conversion of kinetic energy into heat. In general, in the presence of external heat-sources, the influence of this heat-production is small, and few studies have considered how these heat fluctuations are exactly redistributed over scales by the turbulent flow itself. In a recent work, De Marinis et al., 2013, addressed the question how this heat generation influences the wavenumber spectrum of heat fluctuations, using the Eddy-Damped Quasi-Normal Markovian (EDQNM) model, combined with an ad-hoc, spectrally local model. Unfortunately, in the important case of a spatially uniform temperature distribution, their model predicts zero frictional heat production, which is at odds with the physics of turbulent flows, where heat fluctuations are produced as soon as the instantaneous viscous dissipation is not uniform in space. In the present manuscript we will derive an EDQNM model of the heat fluctuation production which is more generally applicable, since it is independent of the initial distribution of temperature fluctuations.

Analysis of the problem

The problem we will consider is conceptually the simplest possible setting in which a non-trivial effect due to the heat generation by turbulent dissipation can be observed: statistically isotropic and stationary incompressible turbulence, governed by the Navier-Stokes equations. The viscous dissipation of kinetic energy is a heat source, and the advection of this scalar, the temperature, will be considered passive, ignoring the influence of thermal expansion of the fluid elements. It will be shown later that the temperature fluctuations are in most situations small, so that this assumption seems very reasonable. The heat diffusion coefficient α, dynamic viscosity µ and density ρ are assumed to be uniform and constant. The temperature is then governed by the following equation,

∂Θ ∂t + u i ∂Θ ∂x i = α ∂ 2 Θ ∂x 2 i + 1 ρc p τ ij D ij , (2.1)
with c p the specific heat, which we will assume to be constant, and where the viscous stress tensor τ ij is given, for a Newtonian fluid, by

τ ij = 2µD ij , (2.2)
and

D ij = 1 2 ∂u i ∂x j + ∂u j ∂x i .
(2.3)

Combining these expressions we obtain [START_REF] Landau | Fluid Mechanics, Second Edition[END_REF],

∂Θ ∂t + u i ∂Θ ∂x i = α ∂ 2 Θ ∂x 2 i + ν c p ∂u i ∂x j ∂u i ∂x j + ∂u i ∂x j ∂u j ∂x i , (2.4)
with ν = µ/ρ the kinematic viscosity. The velocity u i and temperature Θ are functions of space x and time t, but this dependence will be omitted from the notation unless more than one different space locations or time-instants appear in the same expression.

Introducing the Reynolds decomposition, Θ = Θ + θ, with Θ the average temperature and θ the fluctuation, we find the following equation for the mean temperature.

∂ Θ ∂t = ν c p ∂u i ∂x j ∂u i ∂x j = ǫ c p , (2.5)
with ǫ the mean dissipation rate and where the overlined quantities denote averages.

Clearly, the average temperature increases since no heat-sinks exist in our system. In a statistically stationary turbulent flow, where ǫ is constant, the mean temperature will thus increase linearly in time.

The equation for the fluctuations is

∂θ ∂t + u i ∂θ ∂x i = α ∂ 2 θ ∂x 2 i - ν c p ∂u i ∂x j ∂u i ∂x j + ν c p ∂u i ∂x j ∂u i ∂x j + ∂u i ∂x j ∂u j ∂x i . (2.6)
The variance of the temperature fluctuations is then given by

∂θ 2 ∂t = -2α ∂θ ∂x i ∂θ ∂x i + 2 ν c p θ ∂u i ∂x j ∂u i ∂x j + θ ∂u i ∂x j ∂u j ∂x i .
(2.7)

Note that the second term in expression (2.11) of De Marinis et al., 2013, is zero due to homogeneity and the last term is not complete since their basic equation (2.1) is not correct. This does not change the rest of their analysis, which did not explicitly use this information, but it is important if one wants to derive a model analytically from the governing equations, as we will do in the following section.

The temperature distribution over different lenghtscales is given by the temperature spectrum, defined such that

E θ (k)dk = θ 2 , (2.8)
which implies, due to homogeneity and isotropy,

E θ (k)δ(k + k ′ ) = 4πk 2 θ(k)θ(k ′ ), (2.9)
with θ(k) the Fourier transform of the temperature fluctuation. The equation for E θ (k) is straightforwardly derived from (2.6),

∂ ∂t + 2αk 2 E θ (k) = -4πk 2 δ(k -p -q)ik i u i (p)θ(q)θ(-k)dpdq -4πk 2 ν c p δ(k -p -q)p i q i u j (p)u j (q)θ(-k)dpdq -4πk 2 ν c p δ(k -p -q)p i q j u j (p)u i (q)θ(-k)dpdq, (2.10)
with u i (k) the Fourier transform of the velocity fluctuation. This equation is unclosed and three triple correlations need to be determined before we can integrate (2.10). In De Marinis et al., 2013, the first correlation, u i (p)θ(q)θ(-k) was modeled by classic EDQNM theory [START_REF] Herring | A comparative assessment of spectral closures as applied to passive scalar diffusion[END_REF][START_REF] Vignon | Thermal spectral calculation using eddy-damped quasinormal Markovian theory[END_REF]. The second correlation was modeled by the ad-hoc model,

-4πk 2 ν c p δ(k -p -q)p i q i u j (p)u j (q)θ(-k)dpdq ∼ ν c p k 5/2 E(k) E θ (k). (2.11)
As stated in the introduction, model (2.11) cannot be used in the case of a uniform temperature distribution. In that case the temperature fluctuations will stay zero forever, because the above model for the production is proportional to the spectrum E θ (k). One could argue that the model can, perhaps, give correct results for the case in which the temperature spectrum is initially nonzero. This would however be fortuitous, since the temperature is supposed to be passive within our assumptions. It should therefore not have any influence on the injection rate of heat fluctuations. We will in the following section derive an EDQNM model for the viscous production term of heat fluctuations. Subsequently we will carry out simulations to show the temperature spectrum that results from this production term. We will compare these results to those obtained by De Marinis et al., 2013.

Derivation of an EDQNM model for the production of viscous heat fluctuations

We use the procedure outlined in [START_REF] Bos | Reduction of mean-square advection in turbulent passive scalar mixing[END_REF] to derive the EDQNM model for the triple correlations in equation (2.10): first we derive a Direct Interaction Approximation (DIA, [START_REF] Kraichnan | The structure of isotropic turbulence at very high Reynolds numbers[END_REF] of the triple correlation. Then we modify the time-dependence to obtain a closure expression of the EDQNM type. Note that Orszag already pointed out this relation between the Eulerian DIA equations and the EDQNM equations for isotropic turbulence in the original paper [START_REF] Orszag | Analytical theories of Turbulence[END_REF]. The idea behind the DIA is the following: in the triple correlations that we are interested in, the turbulent fluctuations are expanded around a state u (0) , θ (0) in which all but one particular triad are retained. The influence of this particular triad only is denoted u (1) , θ (1) and is treated as a perturbation around the state u (0) , θ (0) . The expansion is then truncated at the lowest non-vanishing order. The expressions for u (1) and θ (1) are in the present case,

u (1) i (k, t) = -i t 0 G(k, t, s) k a δ ib + k b δ ia -2 k i k a k b k 2 u a (p, s)u b (q, s)ds (3.1)
and

θ (1) (k, t) = -i t 0 G θ (k, t, s) (k a u a (p, s)θ(q, s) + k a u a (q, s)θ(p, s)) ds -2 t 0 G θ (k, t, s) ν c p p i q i u j (p, s)u j (q, s)ds -2 t 0 G θ (k, t, s) ν c p p i q j u j (p, s)u i (q, s)ds (3.2)
in which G(k, t, s) and G θ (k, t, s) are the nonlinear response functions of the velocity and temperature fluctuations, respectively. Expression (3.1) and the first line of (3.2) are the standard DIA contributions for the velocity and the scalar [START_REF] Kraichnan | The structure of isotropic turbulence at very high Reynolds numbers[END_REF][START_REF] Roberts | Analytical theory of turbulent diffusion[END_REF]. The second and third line of (3.2) are the DIA contribution to the scalar stemming from the frictional heat source nonlinearity. Substituting these perturbed quantities into the expression for the triple correlations we obtain, after invoking the weak-dependence hypothesis [START_REF] Kraichnan | The structure of isotropic turbulence at very high Reynolds numbers[END_REF],

∂E θ (k) ∂t = T θ (k) + P θ (k) -D θ (k), (3.3) with T θ (k) = 8πk 2 δ(k -p -q)k i k j t 0 Φ ij (p, t, s) × [G θ (k, t, s)Φ θ (q, t, s) -G θ (q, t, s)Φ θ (k, t, s)] dsdpdq, P θ (k) = 16πk 2 ν c p 2 δ(k -p -q) t 0 G θ (k, t, s) × (p m q m ) 2 Φ ij (p, t, s)Φ ij (q, t, s)
+ 2p m q m p i q j Φ aj (p, t, s)Φ ia (q, t, s) + p i p j q m q n Φ mn (p, t, s)Φ ij (q, t, s) ] dsdpdq, (3.4)

D θ (k) = 2αk 2 E θ (k),
where the letters T, P, D denote Transfer, Production and Diffusion, respectively. In this expression we have introduced the quantities

Φ ij (k, t, s)δ(k + k ′ ) = u i (k, t)u j (k ′ , s), Φ θ (k, t, s)δ(k + k ′ ) = θ(k, t)θ(k ′ , s).(3.5)
We consider the isotropic, mirror-symmetric case so that, introducing the time correlation functions R(k, t, s) and R θ (k, t, s),

Φ ij (k, t, s) = Φ ij (k, t, t)R(k, t, s), Φ θ (k, t, s) = Φ θ (k, t, t)R θ (k, t, s), (3.6)
we can write

Φ ij (k, t, s) = E(k) 4πk 2 δ ij - k i k j k 2 R(k, t, s), Φ θ (k, t, s) = E θ (k) 4πk 2 R θ (k, t, s). (3.7)
The expressions for T θ (k) and P θ (k) can be rewritten, after some algebraic manipulations, as

T θ (k) = ∆ Θ I (k, p, q)(1 -z 2 )k 3 q 2 E(p) E θ (q) q 2 - E θ (k) k 2 dpdq pq , P θ (k) = 2 ν c p 2 ∆ Θ II (k, p, q)(4x 4 -3x 2 + 1)kp 2 q 2 E(p)E(q) dpdq pq (3.8) in which Θ I (k, p, q) = t 0 G θ (k, t, s)R θ (q, t, s)R(p, t, s)ds Θ II (k, p, q) t 0 G θ (k, t, s)R(q, t, s)R(p, t, s)ds. (3.9)
and in which ∆ indicates the domain in the pq-plane in which k, p, q can form a triangle (in other words |p -q| k |p + q|) and x, y, z are given by x = -p i q i /(pq)

y = k i q i /(kq) z = k i p i /(kp). (3.10)
Assuming a fluctuation-dissipation relation for the two-time quantities and exponentially decaying time-correlations,

R(k, t, s) = G(k, t, s) = exp(-η(k)(t -s)) for t s, (3.11) R θ (k, t, s) = G θ (k, t, s) = exp(-η θ (k)(t -s)) for t s, (3.12)
these expressions simplify, in the long-time limit, to

Θ I (k, p, q) = 1 η θ (k) + η θ (q) + η(p) Θ II (k, p, q) = 1 η θ (k) + η(p) + η(q) . (3.13)
It is at this point that the choice of the form of the typical frequencies η and η θ will determine the inertial range scaling of the spectra. We will choose the standard forms [START_REF] Herring | A comparative assessment of spectral closures as applied to passive scalar diffusion[END_REF],

η(k) = 1.3 k 0 s 2 E(s)ds + νk 2 and η θ (k) = αk 2 . (3.14)
We note here that the expression for T θ (k) is the standard EDQNM closure for scalar fluctuations [START_REF] Herring | A comparative assessment of spectral closures as applied to passive scalar diffusion[END_REF][START_REF] Vignon | Thermal spectral calculation using eddy-damped quasinormal Markovian theory[END_REF]. The expression for the viscous heat production term P θ (k) is new. In principle, the presence of the frictional heating term could change the behaviour of the scalar frequency η θ (k). This possibility is not explored in the present investigation.

The equation for the energy spectrum is equivalently closed by the EDQNM model [START_REF] Orszag | Analytical theories of Turbulence[END_REF],

∂ ∂t + 2νk 2 E(k) = ∆ Θ(k, p, q)(xy + z 3 )k 3 p 2 E(q) E(p) p 2 - E(k) k 2 dpdq pq + F (k), (3.15)
with F (k) a constant forcing term acting in the small wavenumber range and

Θ(k, p, q) = 1 η ′ (k) + η ′ (p) + η ′ (q) . (3.16) with η ′ (k) = 0.36 k 0 s 2 E(s)ds + νk 2 .
(3.17)

It was shown in the past that the variance of the dissipation-rate fluctuations, involving fourth order correlations of the velocity field, was not correctly captured by DIA related approaches [START_REF] Chen | Non-Gaussian statistics in isotropic turbulence[END_REF]. The present model involves the correlation between the dissipation-rate fluctuation and the temperature fluctuation. This is a mixed velocity-temperature third order correlation. In general, it seems that DIA more reliably represents third than fourth order correlations, but a validation of the model (3.8) by direct numerical simulation would be valuable, in particular since the modeled correlation involves the fluctuation of dissipation rate.

Heat production at unit Prandtl number

Preliminary analysis of the closure expression

The goal is here to first derive expressions for the heat-fluctuation spectrum generated by frictional heating, by considering the evolution equation (3.3) and the EDQNM closure for the frictional heating term P θ (k) in this equation. Then we will verify these analytical results by numerically integrating the closure equations. For our analytical considerations we will assume the energy spectrum to be constant in time and to obey Kolmogorov scaling, E(k) ∼ ǫ 2/3 k -5/3 , from the wavenumber corresponding to the integral scale down to the Kolmogorov wavenumber k η ∼ ǫ 1/4 ν -3/4 . To evaluate the convolution integral in (3.4), there are three distinct cases to consider. Firstly, the convolution integral is dominated by local interactions. Secondly, the integral is dominated by infrared nonlocal interactions p ≪ k ≈ q. Thirdly the integral is dominated by ultraviolet nonlinear interactions k ≪ p ≈ q. We have verified the convolution integral analytically for all three types of contributions and we have found that in the case of a very long inertial range, the nonlocal interactions contribute insignificantly compared to the local interactions. We will therefore only consider the approximation of (3.3) for the case of local interactions. In that case (3.4) can be approximated by its dimensional estimate,

P θ (k) ∼ ν c p 2 k 5 E(k) 2 τ (k), (4.1)
where τ (k) is the local estimate of the triad interaction time Θ II (k, p, q). If we suppose

τ (k) ∼ ǫ -1/3 k -2/3 , (4.2)
and we assume Kolmogorov-scaling we find

P θ (k) ∼ ν c p 2 ǫk. (4.3)
In the case of unit and large Prandtl number, the diffusion will be small in the inertialconvective range, and we expect a balance to establish between the production term and the transfer term,

P θ (k) ≈ T θ (k). (4.4)
This balance will indeed be observed in the numerical evaluation of the closure. In this range, it was shown [START_REF] Rubinstein | Reassessment of the classical closures for scalar turbulence[END_REF]) that the scalar transfer can be qualitatively predicted by a Kovaznay-type of gradient model, we find from (4.3), (4.4) and (4.5),

T θ (k) ∼ ∂ ∂k E θ (k)E(k) 1/2 k 5/2 , ( 4
E θ (k) ∼ ν c p 2 ǫ 2/3 k 1/3 . (4.6)
Integrating this spectrum up to the kolmogorov scale we have,

θ 2 ∼ ǫν c 2 p ∼ R -2 λ . (4.7)
In a steady state, the total production will be balanced by the scalar dissipation, so that,

P θ (k)dk = ǫ θ ≡ D θ (k)dk. (4.8)
Assuming a scaling of the production spectrum of the form (4.3) up to the Kolmogorov scale k η , we have

ǫ θ ∼ ǫ 3/2 ν 1/2 c 2 p , ( 4
.9) which shows that in the limit R λ → ∞,

ǫ θ ∼ R -1 λ . (4.10)
This means that in the limit of vanishingly small viscosity, both the mean-square heat fluctuations and their dissipation tend to zero.

Numerical integration of the closure for unit Prandtl number

Equations (3.3) and (3.15) are integrated numerically on a logarithmically spaced grid with 22 wavenumbers per decade for Reynolds number R λ = 100, 300, 1000. The specific heat c p is kept constant, since its variation only induces a trivial change on the results. Indeed, since the advection equation is linear and depends linearly on c -1 p , the temperature spectrum will simply be proportional to c -2 p . The initial temperature spectrum is zero. All wavenumber spectra in the remainder of this results section are normalized using Kolmogorov-Corrsin-Obukhov variables, i.e., using the variables ν, ǫ and ǫ θ . The normalized spectra are indicated by a tilde. For instance, the normalized energy spectrum is given by Ẽ and k is normalized by the Kolmogorov lengthscale η = ν 3/4 ǫ -1/4 . Since the velocity field contains a large scale forcing, the energy spectrum will develop towards a steady state. This long-time asymptotic energy spectrum is shown in Figure 1, left. The energy spectra, normalized using Kolmogorov variables, collapse perfectly at high wavenumbers. The temperature spectrum which develops asymptotically due to the viscous heat production is shown in Figure 1, right. It is an increasing function of the wavenumber, roughly proportional to k 1/3 and the spectra rapidly fall off in the far dissipation range. For different Reynolds numbers the temperature spectra collapse, using Kolmogorov-Corrsin-Obukhov variables,

(k) = ν -5/4 ǫ -1/4 E(k), (
Ẽθ (k) = (ǫ θ ν 5/4 ǫ 3/4 ) -1 E θ (k).
(4.12)

In Figure 2 we show the three different contribution terms in the evolution equation of the temperature spectrum: transfer T θ (k), production by frictional heating P θ (k) and dissipation of temperature fluctuations D θ (k). We observe that the main part of the production is locally (in wavenumber) dissipated by the diffusion. The transfer term allows for a small transfer of temperature variance towards smaller scales. But since the production takes mainly place in the diffusive range (at unit Prandtl number), the fluctuations do not survive long enough to be strained significantly by the velocity field, and the transfer is therefore confined to a small wavenumber span in the dissipation range. This picture will change at large Prandtl number. In Figure 2, right, we show the same con-tributions, but in a logarithmic representation. This allows to more closely focus on the small and intermediate wavenumber ranges. In the inertial range, the production term is proportional to k. The production is largely dominant over the diffusion D θ (k) in this range and it is therefore the transfer term which balances the production, as assumed in expression (4.4). The wavenumber dependence of P θ (k) and E θ (k) are as predicted in expressions (4.3) and (4.6), which were obtained by assuming that the convolution integral in (3.4) is dominated by local interactions. A further consequence of this assumption is that the mean-square temperature fluctuations should be proportional to R -2 λ . The dissipation of temperature fluctuations should be inversely proportional to the Reynolds number. This is indeed observed in Figure 3.

5. Temperature spectrum in the inertial-diffusive and viscous-convective range 5.1. Phenomenology for small and large Prandtl numbers The energy spectrum rapidly decays for k > k η . Since the convolution integral in (3.8) contains interactions with E(p) and E(q) under the constraint that |p -q| k |p + q|, negligible production will be observed beyond k = 2k η . For very high Prandtl number we therefore expect a Batchelor-type phenomenology in which the scalar fluctuations are strained by the smallest available velocity scale, around k η [START_REF] Batchelor | Small-scale variation of convected quantities like temperature in turbulent fluid. Part 1. conductivity[END_REF]. This scaling is given by

E θ (k) ∼ ǫ θ ν ǫ k -1 . (5.1)
Substituting ǫ θ by expression (4.9), we have

E θ (k) ∼ ǫν c 2 p k -1 . (5.2)
For the opposite case, in which the Prandtl number is very small, the production term will be balanced by the diffusion term. The production term, was approximated by

P θ (k) ∼ ν c p 2 k 5 E(k) 2 τ (k), (5.3)
assuming dominance of the local interactions in the generation of the heat fluctuations. However, the triad-correlation timescale is in this case not given by (4.2), since at these scales the dynamics are now dominated by the diffusion of heat, so that

τ (k) ∼ 1/(αk 2 ), (5.4)
and the production spectrum in the inertial-diffusive range will then be given by

P θ (k) ∼ ν c p 2 ǫ 4/3 α k -1/3 . (5.5)
Since the diffusion spectrum is

D θ (k) = 2αk 2 E θ (k), (5.6)
we find from

P θ (k) ≈ D θ (k), (5.7) that E θ (k) ∼ ν c p 2 ǫ 4/3 α 2 k -7/3 .
(5.8)

Numerical integration of the closure for small and large Prandtl number

The results for a large Prandtl number simulation (P r = 100, R λ = 100) are shown in Figure 4. The k -1 spectrum starts to appear but is not yet fully developed at this Prandtl number. Since the peak of the production mechanism remains fixed around k η , independently of the Prandtl number, there is no reason why it should not become more pronounced at higher Prandtl number. The balance of diffusion and transfer is observed in Figure 4, right, for wavenumbers larger than k η where the energy spectrum drops rapidly.

Figure 5, left, shows that, at low Prandtl number, the temperature spectrum is indeed approximately proportional to k -7/3 , in agreement with expression (5.8). The balance between the diffusion term and the heat fluctuation source term holds, as is observed in Figure 5, right. The production spectrum is not exactly given by a powerlaw proportional to k -1/3 , which is consistent with the fact that the temperature spectrum is not described by a pure powerlaw either.

Comparison with previous model predictions

It was shown in the previous sections that according to the closure we derived, the generation of heat fluctuations is dominated by local interactions in wavenumber space. This allowed us to show that a simple, phenomenological model of the production is

P θ (k) ∼ ν c p 2 k 5 E(k) 2 τ (k), (6.1) with τ (k) ∼ ǫ -1/3 k -2/3 for P r 1, kη ≪ 1 1/(αk 2 ) for P r ≪ 1, kη ≪ 1 . ( 6
E θ (k) ∼ ν c p 2 ǫ 2/3 k 1/3 . (6.3)
If we consider the low Prandtl number case, balance (5.7) should hold in the inertial range and the temperature spectrum will be proportional to (6.4) as also proposed by De Marinis et al., 2013. The predictions of their model for these particular cases are therefore in agreement with the results of our EDQNM model. This might seem surprising, given the fact that the model does not give physical results for the basic case in which no initial fluctuations are present. We note however that all the models of the form

E θ (k) ∼ ν c p 2 α -2 ǫ 4/3 k -7/3 ,
P θ (k) ∼ E θ (k) γ ν c p 2(1-γ) E(k) 3/2-γ k 7/2-2γ (6.5)
will give the same steady-state inertial-convective range behavior, independent of the value of γ. The only choice predicting a non-zero production in the absence of initial fluctuations is γ = 0. The value that will give correct scaling in the inertial-diffusive range, at low Prandtl number is γ = 1/2 corresponding to the model of De Marinis et al., 2013. It is therefore not possible to find a model of the form (6.5) which gives behaviour in agreement with both (4.6) and (5.8) and which predicts nonzero production in the absence of initial fluctuations. The simplest phenomenological model that satisfies these criteria is (6.1), involving a timescale which can mimic both diffusive and inertial behaviour.

We note that as soon as initial heat fluctuations are present in the considered flow, as is a prerequisite for (2.11) to generate heat, the distribution of these fluctuations will determine the instantaneous heat generation and all transients will be affected by this. The influence of these unphysical transients can be very long if the initial heat fluctuations are present in the largest flow scales, which decay slowly.

Conclusion

An EDQNM model for the heat-production term due to viscous friction is derived. Temperature wavenumber spectra are deduced from this closure. In particular, it is shown that the heat-production is mainly generated through local triad interactions. The heat production spectrum is then proportional to k. Balancing this production spectrum against nonlinear transfer, we find the heat-spectrum for unit and large Prandtl number in the inertial-convective subrange. In the very diffusive case (very small Prandtl number), the production is in equilibrium with the diffusivity spectrum, leading to a different wavenumber dependence. The predictions are Let us consider typical values for c p and ν. For a gas such as air, c p ≈ 10 3 Jkg -1 K -1 , ν ≈ 10 -5 m 2 s -1 . In order to observe heat-fluctuations of the order of 10 -3 K, a dissipation rate of the order of 10 5 W kg -1 is needed. Combining this with equation (2.5), this implies that the mean temperature of the gas increases at a rate of 10 2 Ks -1 . For a typical liquid such as water, this increase is even larger. This example illustrates the smallness of the temperature variance in practice. It does not seem to us that in a practical situation in which the temperature of a fluid increases by 100K every second, fluctuations of the order of a millikelvin are of major importance. Note that the scalings (7.1) are observed approximately from the numerically integrated closure expressions. We stress that these are the results which are given within the framework of the EDQNM assumptions and the model cures for the defect of a previous model proposition which excluded heat-fluctuation generation in the absence of initial fluctuations [START_REF] De Marinis | Temperature dynamics in decaying isotropic turbulence with Joule heat production[END_REF]. A validation of the EDQNM model for the heat production and a verification of the proposed scalings (7.1) using high resolution direct numerical simulations of isotropic turbulence would certainly be valuable, in particular since DIA and EDQNM might have difficulties in representing correlations involving fluctuations of the dissipation rate.

E θ (k) ∼         

Figure 1 .

 1 Figure 1. Left: Energy spectrum, normalized by Kolmogorov variables, for three different Reynolds numbers. Right: corresponding temperature fluctuation spectrum, generated by frictional heating (P r = 1).

Figure 2 .Figure 3 .

 23 Figure 2. Left: different terms in the evolution equation of the temperature spectrum (P r = 1): transfer T θ (k), production by frictional dissipation P θ (k) and dissipation of temperature fluctuations D θ (k). Right: logarithmic representation.

Figure 4 .

 4 Figure 4. Left: E θ (k) for R λ = 100, P r = 100. Right, transfer T θ (k), production P θ (k) and dissipation spectrum D θ (k).

Figure 5 .

 5 Figure 5. Left: E θ (k) for R λ = 1000, P r = 10 -4 . Right, transfer T θ (k), production P θ (k) and dissipation spectrum D θ (k).

  .2)Comparing the form of the above model to the model of De Marinis et al., 2013, the main difference is the absence of the heat fluctuation spectrum in (6.1). In the case of the model of De Marinis et al., 2013, the production term cannot be determined independently from the temperature fluctuation distribution. All cases we considered in the present work would thus yield zero temperature fluctuations using that model. However, in the particular case of an initial heat-fluctuation spectrum, in the absence of further heat-sources and for very long times, it is possible to compare the predictions of the model of De Marinis et al., 2013 with the present case. For instance in the inertialconvective range we expect a balance between nonlinear-transfer of scalar fluctuations and production. Combining expression (4.4), (4.5) and (2.11), we find in this particular case for the temperature fluctuation spectrum in the inertial-convective range (assuming Kolmogorov scaling),

  have obtained these expressions for the heat-fluctuation spectra, it is possible to estimate how large the fluctuations are in reality. Expression (4.7) shows that, assuming Kolmogorov scaling, we have

Ackowledgments

The author thanks Robert Rubinstein for his invaluable comments and acknowledges discussions with Marcello Meldi, Sergio Chibbaro and Pierre Sagaut.


Muhammad Muzzamil Luqman
email: muhammad_muzzamil.luqman@univ-lr.fr

Petra Gomez-Krämer

Jean-Marc Ogier
email: jean-marc.ogier@univ-lr.fr

Mobile phone camera-based video scanning of paper documents

Mobile phone camera-based document video scanning is an interesting research problem which has entered into a new era with the emergence of widely used, processing capable and motion sensors equipped smartphones. We present our ongoing research on mobile phone camera-based document image mosaic reconstruction method for video scanning of paper documents. In this work, we have optimized the classic keypoint feature descriptor based image registration method, by employing the accelerometer and gyroscope sensor data. Experimental results are evaluated using optical character recognition (OCR) on the reconstructed mosaic from mobile phone camera-based video scanning of paper documents.

I. INTRODUCTION

In recent years, the availability of camera equipped, processing capable, inertial sensors fitted and moderate priced mobile phones (a.k.a. smartphones), has attracted the attention of the research community to employ them for complementing the classical document scanning devices. The use of these devices for document scanning provides interesting advantages over the traditional document scanning devices. They can be used to scan thick books, historical documents that are too fragile to touch, text in scenes (walls, whiteboards, etc.), and large sized documents [START_REF] Liang | Camera-based analysis of text and documents: a survey[END_REF]. However, the use of smartphones introduces new challenges to document scanning which are not faced by classical document scanning devices. These challenges include uneven lighting, perspective distortion, nonplaner surfaces, motion blur and low resolution of the cameras [START_REF] Liang | Camera-based analysis of text and documents: a survey[END_REF] [START_REF] Hannuksela | Document image mosaicing with mobile phones[END_REF].

In this paper we present our ongoing work on mobile phone camera-based video scanning of paper documents. The video scanning of a paper document is achieved by swiping the mobile phone camera over the paper document and recording the accelerometer sensor data along with capturing the video frames. During the video scanning the orientation of the mobile phone camera is obtained from the gyroscope data and the user is provided with visual feedback on the orientation of the phone to avoid perspective distortion. A complete mosaic image of the paper document is reconstructed from the captured video frames by employing an optimized keypoint feature descriptor-based image registration techniques. The optimization is achieved by employing the recorded accelerometer sensor data. The resulting reconstructed mosaic has a higher resolution than a simple photo of the document taken by same camera.

In literature the camera-based scanning of paper documents has been approached by various works which are mainly motivated by panorama reconstruction and image mosaicing techniques from the computer vision research community. In [START_REF] Hannuksela | Document image mosaicing with mobile phones[END_REF] first image feature-based technique is used to estimate the camera motion and to assist the user to capture images of patches of document. The estimated camera motion is used with a keypoint feature descriptor-based technique for registration of captured image patches and reconstruction of a mosaic of the document. In [START_REF] Sawhney | Robust video mosaicing through topology inference and local to global alignment[END_REF] an algorithm for 2D scanning of a planar scene is proposed. The topology of the video frames are inferred on a 2D manifold by alignment of successive video frames and overlapping video frames. The aligned frames are merged by using a multi-resolution method for constructing a seamless mosaic. In [START_REF] Nakai | Camera-based document image mosaicing using LLAH[END_REF] local likely arrangement hashing (LLAH), which is originally an image retrieval technique, is used for keypoint detection and feature description in frames. Images are aligned by matching LLAH feature descriptors and the feature correspondences are used for combining input frames for reconstruction of mosaic. In [START_REF] Liang | Mosaicing of cameracaptured document images[END_REF] first captured frames are rectified by removing perspective distortion using texture flow information. A hough transform-based voting scheme is used for finding translation and scaling between video frames. The reconstruction of the mosaic is achieved by a sharpness-based seamless composition of overlapping images. In [START_REF] Yang | Inertial sensors aided image alignment and stitching for panorama on mobile phones[END_REF] inertial sensors in mobile phones have been employed for constructing panoramas on mobile phones. In a first step the position and relative displacement of video frames are computed by inertial sensor data. Using the alignment estimation from inertial sensor data, a more precise alignment of the video frames is computed by using a keypoint feature descriptorbased technique and the mosaic image is constructed by using the feature correspondences.

The perspective distortion is very important to be handled in case of camera-based document scanning. A document image mosaicing technique should directly or indirectly rectify perspective distortion of the captured frames before reconstructing the mosaic image. A summary of methods for content based correction of the perspective distortion in camera captured document images is presented in [START_REF] Jagannathan | Perspective correction methods for camera based document analysis[END_REF].

In this paper we present our ongoing research on document image mosaicing. We are inspired by the work in [START_REF] Yang | Inertial sensors aided image alignment and stitching for panorama on mobile phones[END_REF] for employing the inertial sensors for document mosaic image reconstruction. However, we are working on elaborating a lightweight algorithm that could be implemented on smartphones. The two novel contributions of our work are the following: 1) We use the gyroscope sensor to give visual user feedback to avoid perspective distortion during the video scan of the paper document whereas perspective distorsion is not considered in [START_REF] Yang | Inertial sensors aided image alignment and stitching for panorama on mobile phones[END_REF]. 2) We compute the direction of swipe to optimize the keypoint feature descriptor-based image registration method by using only accelerometer sensor data whereas the authors of [START_REF] Yang | Inertial sensors aided image alignment and stitching for panorama on mobile phones[END_REF] optimize a keypoint feature descriptor-based image registration method by computing the displacement of the mobile phone from accelerometer and gyroscope sensor data.

The remainder of this paper is organized as follows. We present a detailed description of our method of video scanning of paper documents in Section II. In Section III we discuss the experimental evaluation and the results. In Section IV we conclude our work and present future directions of research.

II. MOBILE PHONE CAMERA-BASED VIDEO SCANNING

In this section we present a detailed description of our method for mobile phone camera-based video scanning of paper documents. We first describe the capturing of video frames and accelerometer sensor data along with the gyroscope based visual feedback for avoiding perspective distortion. This is followed by description of our methodology for finding the direction of swipe (of video scanning), from the accelerometer data recorded with the captured frames in the video sequence. Finally we describe the image registration of the frames of the captured video and the reconstruction of the complete mosaic image of the paper document.

A block diagram of our method for mobile phone camera based video scanning of paper documents is presented in Fig. 1.

A. Video scanning of paper documents

The video scanning of a paper document is achieved by a one-dimensional swipe of the mobile phone camera on the paper document. The swipe could either be from the top to the bottom of the document or from the bottom to the top of the document. During the video scanning we record the accelerometer sensor data along with capturing of the video frames.

Processing of accelerometer sensor data: An accelerometer in a smartphone measures the acceleration α = (α x , α y , α z ) of the phone in each direction of the X, Y and Zaxis. The accelerometers in smartphones are usually not of very high quality (because of cost constraints) and thus the obtained acceleration data is very noisy. Hence, the raw accelerometer readings are full of random noise and are unusable in their original form. The rise in temperature of the mobile phone (resulting from camera and screen heat) increases the random noise in the accelerometer data. In order to make sure that the accelerometer reading is as close to the real value as possible, we compute the calibration offset of the accelerometer sensor by placing the phone on a flat surface and averaging x readings along each of the three axes separately and independently. We then subtract the calibration offset from the future readings of the accelerometer; hence obtaining calibrated readings. For removing random noise from the accelerometer readings we smooth the accelerometer readings by using a Kalman filterbased running average. The smoothed accelerometer readings are recorded with the video frames. Fig. 2 presents respectively the raw (noisy) accelerometer readings for the whole video scanning of a paper document, the calibrated and smoothed accelerometer readings for the whole video scanning of a paper document, the raw accelerometer data associated to the captured frames and the calibrated and smoothed accelerometer data associated to the captured frames, for the X-axis of the phone accelerometer sensor.

Use of gyroscope sensor data for visual user feedback:

In order to avoid perspective distortion, we employ the gyroscope sensor data for obtaining the orientation of the mobile phone camera. Modern smartphones are fitted with 3D gyroscopes, which measure the angular velocity ω = (ω x , ω y , ω z ) along the X, Y and Z-axis of the phone. The angular velocity along the X-axis is termed pitch, along the Y-axis yaw and along the Z-axis roll. We obtain the angular velocity from the gyroscope as an angle of rotation around each of the three axes. We then employ it for providing a visual feedback to the user in order to keep the mobile phone camera parallel to the document plane and so to avoid the perspective distortion. We show three triangles on the screen for the visual feedback to adjust the orientation of the device. The three triangles are mapped to the gyroscope data along the three axis, respectively. A change in orientation of the device along an axis, rotates its respective triangle on the screen. By aligning any two of the three triangles, the user can keep the mobile phone parallel to the XY, YZ or ZX plane. To illustrate it pictorially, the Fig. 3 presents some screenshots of the interface during the video scanning of a paper document.

B. Finding the direction of swipe by using accelerometer data

The accelerometer sensor data measures the acceleration of the phone along the X, Y and Z-axis. We use the calibrated and smoothed accelerometer data for inferring the swipe direction of the mobile phone during the video scanning of the paper document. The accelerometer readings along each of the three axes is between -1g and +1g (where g = 9.8m/s 2 ). The document page is placed on a planar surface (e.g. a table) and we use the phone in landscape mode (as shown in the screenshots of Fig. 3) for video scanning of documents. This setup makes the mobile phone's X-axis as the primary axis of swipe. The sensor coordinate system, indicating the X, Y and Z-axis of the mobile phone that we use for our research, is shown in Fig. 4.

To compute the swipe direction of the mobile phone (from the top to the bottom or from the bottom to the top of the document page) we employ a very simple methodology. We count the number of positive and negative readings in the recorded accelerometer data of the frames of a video scan. If there are more negative values than positive ones, this means that the phone is swiped in negative direction of the X-axis. And if there are more positive values than the negative ones, this means that the phone is swiped in positive direction of the X-axis. A swipe in negative X-axis direction corresponds to a top to bottom video scan of a paper document whereas a swipe in positive X-axis direction corresponds to a bottom to top video scan. This simple methodology is robust, efficient and very useful for detecting the direction of swipe during the video scanning of paper documents.

C. Image registration of the video frames

The captured frames from video scanning of the paper document are registered by employing a keypoint feature descriptor-based technique. Successive frames in the video sequence are aligned and the projective transform or homography is computed between them. Afterwards, the computed homographies between successive frames are employed for calculating the homography of each frame to a reference frame in the video sequence.

Optimized alignment of successive video frames: For the alignment of successive frames in the video sequence, we use a keypoint feature descriptor-based alignment technique. As discussed by [START_REF] Szeliski | Image Alignment and Stitching[END_REF], the feature based image alignment techniques have some very interesting advantages over pixel-based image alignment techniques. Namely for mobile phone camera-based video scanning of paper documents, the feature based image alignment methods are more efficient and robust than the pixel based methods; specially in case of uneven lighting and scene motions. There are many keypoint detectors and feature descriptors in state-of-the-art. They include the famous SIFT [START_REF] Lowe | Distinctive image features from scale-invariant keypoints[END_REF], SURF [START_REF] Bay | Surf: Speeded up robust features[END_REF], FAST [START_REF] Rosten | Machine learning for high-speed corner detection[END_REF], ORB [START_REF] Rublee | ORB: An efficient alternative to SIFT or SURF[END_REF] and FREAK [START_REF] Alahi | FREAK: Fast Retina Keypoint[END_REF].

To align two successive frames, we first employ the SIFT keypoint detector to obtain a set of keypoints in two frames. Second, we extract the SIFT feature descriptors on each of the detected keypoints in the two frames. Third, we perform FLANN based feature matching [START_REF] Muja | Fast approximate nearest neighbors with automatic algorithm configuration[END_REF] between the two frames and employ RANSAC [START_REF] Fischler | Random sample consensus: a paradigm for model fitting with applications to image analysis and automated cartography[END_REF] for refining the initial correspondences obtained by FLANN.

In order to optimize keypoint detection, feature descriptor computation and feature matching, we use the direction of mobile phone swipe during video scanning. During the keypoint detection and feature descriptor computation phases we use the direction of swipe to avoid processing the complete image and to ignore the top and bottom parts of successive frames (or vice versa depending upon the direction of swipe). As a result, it reduces the search space to be exploited during feature matching. The size of top and bottom parts of successive frames respectively to be ignored considering the direction of swipe is controlled by a parameter which is computed automatically for a video sequence and it takes into account the resolution of video frames (in pixels) and the speed of the swipe in the sequence (assumed to be constant during the video scanning). The speed of swipe is estimated from the total number of video frames in the video scan of a document page. The parameter for ignoring the top and bottom parts of successive frames is computed as:

G = h n (1)
where G denotes the number of pixel rows to be ignored on the top and bottom of successive frames respectively, h is the height of captured video frames in pixels and n is the total number of frames in video sequence.

Homography computation for two successive frames: We compute the planar homography or projective transform between two successive frames by minimizing the backpropagation error and further refine the computed homography by using the Levenberg-Marquardt method to minimize the backpropagation error [START_REF] Szeliski | Image Alignment and Stitching[END_REF]. The homography or projective transform between two successive frames is represented by a homography matrix.

For a captured video sequence of n frames given by:

V = {f 1 , f 2 , f 3 , . . . , f n-1 , f n } (2)
the set of homographies between successive frames is:

h = {h (1,2) , h (2,3) , h (3,4) , . . . , h (n-1,n) } (3)
where h (i,j) is the homography between frames f i and f j .

Homography between non-successives frames: Employing the well established properties of matrices, the homographies computed for successive frames are employed in a cascade matrix multiplication, for computing the homographies between non-successive frames. For example for computing the homography h [START_REF] Liang | Camera-based analysis of text and documents: a survey[END_REF][START_REF] Liang | Mosaicing of cameracaptured document images[END_REF] between the frames f 1 and f 5 , the homographies h (1,2) , h (2,3) , h [START_REF] Sawhney | Robust video mosaicing through topology inference and local to global alignment[END_REF][START_REF] Nakai | Camera-based document image mosaicing using LLAH[END_REF] and h (4,5) are matrix multiplied. This permits us to define a homography between any pair of frames in the video sequence (whether successive or non-successive).

D. Reconstruction of the complete mosaic

The complete mosaic image of the video scanned paper document is constructed by a projection of the pixels in captured frames onto a reference frame. Thus, we first select a reference frame in the captured video sequence. For simplicity, we suppose here that the first frame f 1 is selected as the reference frame. Then, we compute the size of the complete mosaic image by projecting the four corners of each frame in the video sequence using the corresponding homography matrix of the frame. To avoid holes or missing pixels in the complete mosaic image, the construction of the complete mosaic image is achieved by inverse projection of each of its pixels onto the sequence of frames. If a pixel of the mosaic is projected onto a subpixel in a frame, we use bilinear interpolation for computing the subpixel intensity from the intensities of the neighboring pixels. If a pixel of the mosaic is projected onto more that one frames of the video sequence we use the median value of the intensities of corresponding pixel of those frames.

III. EXPERIMENTATION

In this section, we evaluate our method for mobile phone camera-based video scanning of paper documents on video frames captured at a resolution of 1280x720 pixels captured by a Nokia Lumia 920 smartphone. A custom application is developed for the capture phase of video scanning of paper documents. Some screenshots of this application are presented in Fig. 3. The video scanning application runs only in landscape mode to force the user to hold the smartphone with two hands and thus ensuring a stable orientation of the smartphone during capture.

During these preliminary experimentations the mosaic image reconstruction was performed on a laptop computer. The video scanning of paper documents was performed in an office environment with normal lighting conditions. The phone was kept parallel to the document plane by following the visual feedback on the orientation of the mobile phone camera, and the swipe was performed slowly and carefully.

The experimentation dataset comprises fifteen A4 sized pages of scientific research papers; containing mostly printed textual content (in English). Some document pages contain also tables and mathematical equations. The document pages were printed on A4 pages and were digitized in three different modes, as given below:

1) image scanned by a classic scanner at 300dpi grayscale 2) photo taken by the smartphone at a resolution of 1280x720 pixels 3) video scanned by the smartphone at a resolution of 1280x720 pixels

The image scanned by a classic scanner at 300dpi serves as reference for evaluating the quality of reconstructed mosaics. We use the Levenshtein distance [START_REF] Levenshtein | Binary codes capable of correcting deletions, insertions and reversals[END_REF] as metric for comparing the OCR results of the images with ground truth. The Levenshtein distance between two string sequences is the edit distance between them i.e. the minimum number of single character edits (insert, delete, replace) required to change one string sequence into the other.

Table I presents a comparison of the Google Drive OCR results on documents for the three digization modes. The OCR results on the reconstructed mosaic from video scanning of the paper document at 1280x720 pixels, are better than the results on the single image captured at 1280x720 pixels. However they are not as good as those of the classical scanner scanned image. One important reason for this is that our method does not perform any camera calibration i.e. any wide angle lens correction on the captured video frames. The wide angle lens noise is thus inherited by the mosaic image and it eventually effects the OCR results. A second reason is that when the document pages were placed on the table for video scanning there can be small curvature at the corners whereas in case of a scanner this curvatures are flattened by closing the scanner lid. Our method does not perform any preprocessing of the document page. Some examples of reconstructed mosaic images in grayscale from mobile phone camera-based video scanning of A4 sized paper documents are presented in Fig. 5.

Apart from the advantages of portability and liberty of scanning the documents of different sizes, the video scanning of paper documents is interesting as it allows to obtain the mosaic image at a higher resolution than the resolution of the camera used for capturing the frames. A mobile phone camera with a resolution of 1280x720 pixels can take a single image of an A4 sized page at 98 dpi as given by:

Fig. 1 .

 1 Fig. 1. Block diagram for mobile phone camera-based video scanning of paper documents.

  (a) Raw accelerometer readings of video scanning of a document. (b) Calibrated & smoothed accelerometer readings of video scanning of a document. (c) Raw accelerometer readings associated with video frames. (d) Calibrated & smoothed accelerometer readings associated with video frames.

Fig. 2 .

 2 Fig. 2. Accelerometer data recorded for the X-axis of the mobile phone sensor during video scanning of a paper document.

  (a) Misaligned triangles resulting in perspective distortion. (b) Well-aligned triangles avoiding perspective distortion.

Fig. 3 .

 3 Fig. 3. Screenshots of the interface during video scanning of a paper document.

Fig. 4 .

 4 Fig. 4. Windows phone accelerometer sensor coordinate system.

  dpi = 1280 × 720 8.27 × 11.69 = 98 (4) where, 8.27 × 11.69 is the size of an A4 page in inches.

Fig. 5 .

 5 Fig. 5. Reconstructed mosaic images for A4 sized paper documents from video scanning at 1280x720 pixels.

TABLE I .

 I

			EXPERIMENTAL RESULTS
			Levenshtein distance between ground truth
			and Google Drive OCR results
		# chars	Classical scanner	Single image at	Mosaic from video
	Page	in page	image at 300 dpi	1280x720 pixels	scan at 1280x720 pixels
	01	2569	36	1330	37
	02	2311	44	808	57
	03	1854	10	456	56
	04	2353	10	441	125
	05	2438	26	885	69
	06	2495	22	1061	41
	07	2171	20	678	542
	08	2524	72	623	99
	09	1422	296	785	354
	10	2482	235	1417	487
	11	2085	9	431	22
	12	3286	382	3230	1442
	13	4299	61	3786	111
	14	3638	107	2112	1451
	15	3924	300	3149	606
	Mean	2657	109	1413	367

ACKNOWLEDGMENT

The piXL project is supported by the "Fonds national pour la Société Numérique" of the French State by means of the "Programme d'Investissements d'Avenir", and referenced under PIA-FSN2-PIXL. For more details and resources, visit http://valconum.fr/index.php/les-projets/pixl.

Whereas when the same mobile phone camera is used for video scanning of A4 sized pages, the reconstructed mosaic images as shown in Fig. 5 are at 159 dpi, 174 dpi and 163 dpi respectively.

Table II show a comparison of computation times for the reconstructed mosaics from video scanning at a resolution of 1280x720 pixels of A4 sized pages. Computation times are shown for mosaic construction with and without optimized registration using accelerometer sensor data. Using the optimiztion a mean speed up of 23 seconds could be realized.

IV. CONCLUSION

We have presented our ongoing research on the mobile phone camera based video scanning of paper documents. Our method employs the gyroscope sensor of the phone for providing a visual feedback to the user for avoiding perspective distortion, and the accelerometer sensor of the phone for optimizing the keypoint feature descriptor-based image mosaicing technique. Our preliminary experimentation shows that the video scanning of documents not only allows to reconstruct the full page mosaic image of a document page from its mobile phone camera-based video scanning, but also reconstructs the full page mosaic image at a better resolution than the physical limits of the camera lens. The work is in progress and we are working on the detailed experimental evaluation of the method along with an implementation on the smartphone platform.

Our ongoing research focus is on employing the gyroscope data for correcting perspective distortion of the frames in addition to the visual feedback. A second direction of ongoing research is to use super-resolution techniques for improving the quality of the mosaic image. In near future we will explore the use of the ambient light sensor for incorporating the lighting conditions of the video scan environment in mosaic reconstruction. In medium term we have planned to include a preprocessing step in our system for rectifying various geometric noises from mobile phone camera captured document image frames.