Rafael Linares

Kévin Berger

Melika Hinaje

Bruno Douine

Jean Lévêque

Design of a Vector Magnet Generating up to 3 T with 3 Axis Orientation

Keywords: LTS coils, Helmholtz coil, superconducting applications, vector magnet

The Low Temperature Superconducting (LTS) coils are designed to produce strong magnetic fields at the working temperature of liquid Helium. Some kinds of them have been created to produce strong and homogeneous magnetic fields, as in Magnetic Resonance Imaging (MRI) systems. However, the critical and mechanical parameters of superconducting material limit the performances of the coils. Most of the time, LTS coils are built to produce a homogeneous magnetic flux density in one single direction. Consequently, it is necessary to rotate the subject or the sample, when a study is needed depending on the magnetic field orientation. This is not suitable for special kinds of experiments involving High Temperature Superconducting (HTS) systems or biomedical applications in which the subject must not be rotated. In order to satisfy these requirements, we have designed a new coil setup. The design was made using a 3D electromagnetic software in order to calculate the performances of the coils system, to optimize the magnetic field homogeneity and to find a cost effective solution. The first results reported in this paper concern the design of a magnet of three Helmholtz coils, each of which generates a unidirectional magnetic field (x, y and z). Each Helmholtz coil can generate 1.73 T, so the flux density magnitude in the center of system could reach 3 T. The price of the all system, including the cryostat, is around 45 k€.

I. INTRODUCTION

ESEARCH and development in energy management problems is booming. Since superconductivity has the advantage of saving energy with its low electric losses [START_REF] Hulm | Superconducting Solid Solution Alloys of the Transition Elements[END_REF], it is desirable to make superconducting applications accessible. Superconductors are mainly used to generate high magnetic fields. Some superconducting magnets are used in material and biomedical research [START_REF] Kitaguchi | Development of a Bi-2223 HTS Magnet for 3T MRI System for Human Brains[END_REF]. In these areas, it is sometimes necessary to have a strong magnetic field that can be oriented without the need to modify the sample position.

Thus, a special magnet in which the magnetic flux density can be set with high homogeneity is needed [START_REF] Furuse | Design and Fabrication of HTS Coils for a Vector Magnet[END_REF]. This magnet is composed by several coils which can individually provide a magnetic flux density in an expected direction. This set of coils can produce two and/or three dimensional magnetic fields and is called a vector magnet.

In previous work [START_REF] Kirschvink | Uniform magnetic fields and double-wrapped coil systems: Improved techniques for the design of bioelectromagnetic experiments[END_REF], the study of various coil geometries was carried out, and the goal was to establish important This work is supported by the Region of Lorraine. Rafael-Antonio Linares-Lamus, Kévin Berger, Melika Hinaje, Bruno Douine and Jean Lévêque are with the Groupe de Recherche en Electrotechnique et Electronique de Nancy, Faculté des Sciences et Technologies, Université de Lorraine, 54506 Vandoeuvre-lès-Nancy, France (e-mail: name.surname@univ-lorraine.fr). relationships to produce homogeneous magnetic fields, e.g. the Merrit coil system which consists of four concentric square coils, that provides a large homogeneous space.

The main objective of our work is to find a geometry being able to combine homogeneity, strong field and low cost. In this article, we propose a design of vector magnet with the capability of producing a strong magnetic field in any direction with high homogeneity. In the designed vector magnet, the strength of the magnetic flux density norm |B| in the working space will be 3 T. In order to ensure this, the Helmholtz coil must provide the three components (Bx, By, Bz) of the magnetic flux density with an individual strength of 1.73 T. Subsequently, we plan to feed the coils of this magnet with fuel cells in order to study the stability of such system [START_REF] Hinaje | Using Fuel Cell as a Power Supply for Superconducting Coil[END_REF]. This magnet would probably be the first three axis vector magnet being able to characterize HTS samples with a large working space.

II. MATERIALS AND METHODS

A. Description

We want to characterize High Temperature Superconductor (HTS) samples in a high and directional magnetic field. The desired strength of the magnetic flux density was defined to 3 T, and the working space, where the sample will be placed, is a sphere of 100 mm in diameter. In order to ensure the proper orientation, the magnetic flux density B must have three components, Bx, By and Bz, which can be controlled independently. The access of the sample into the magnet is designed as a tube that passes through the coils and allows two ways to access at the center. Then, various types of coils arrangements were compared for uniformity, cost and ease of construction. Some authors have already reported some interesting results on coils configurations so-called: "Merrit", "Lee-Whiting" [START_REF] Kirschvink | Uniform magnetic fields and double-wrapped coil systems: Improved techniques for the design of bioelectromagnetic experiments[END_REF], "Helmholtz coil" [START_REF] Bronaugh | Helmholtz coils for calibration of probes and sensors: limits of magnetic field accuracy and uniformity[END_REF], [START_REF] Johansson | Helmholtz coils for characterization of magnetic sensors[END_REF], and "sector coil" [START_REF] Rossi | Electromagnetic design of superconducting dipoles based on sector coils[END_REF]. From these studies, we have concluded that the arrangement of coils that best suits to our constraints is the Helmholtz coil. Its construction is less complex and quicker than the others. Moreover, its magnetic field homogeneity and strength are enough for the required tasks.

Helmholtz coil consists of two circular coils of radius R, placed along a common axis, and separated by a distance equal to their radius. It can produce magnetic fields with homogeneity greater than 90 % in the desired working space. This meets the conditions that we have set for the proposed magnet. A Helmholtz coil will be placed along each axis of a Cartesian coordinate system. The whole system will consist of six circular coils, two per axis, shown in Fig. 1.

The value of the radius of each Helmholtz coil is calculated R in order to facilitate the assembly of the coils, and also to fulfill the conditions of homogeneity. The radius Ry of the smallest Helmholtz coil, which is placed along the y-axis, is taken as reference. The other radii are then calculated according to Ry as follow: Rx = 1.310 Ry for the Helmholtz coil following the x-axis and Rz = 1.644 Ry for the z-axis. The Helmholtz coils are serially coupled to each other. The parameters of the vector magnet were shown in Table I.

The cross-section of each coil is square shape. Indeed, this leads to a biggest space of homogeneous magnetic flux density compared to circular shape or rectangular shape. Numerical calculations have shown the influence of the cross section shape of a coil on the homogeneity of B [START_REF] Rossi | Electromagnetic design of superconducting dipoles based on sector coils[END_REF]. Moreover, the cross-sectional area of each coil depends on its number of turns:

( ) 2 , ,
, , wire / 2 / ,

x y z x y z S N d F π = (1) 3/ 2 , ,
, , , , 0 5 with 0.75 and 4

x y z x y z x y z B R F N I µ   = =     (2)
The coil system is placed in the air and does not contain iron in its construction. Indeed, the magnetic saturation of the iron can affect the homogeneity of B in the working space and can produce undesirable induced currents when B varies. In order to avoid these perturbations, a resistant structure that can hold the vector magnet in an accurate position has to be built: assembly pieces in epoxy resin were chosen. The influence of the thermal contraction of the epoxy resin structure has been neglected at this time. However, it will be evaluated at the final manufacturing stage. Hence it is important to ensure the clamping of the coils as it directly affects the homogeneity. The mechanic constraints are calculated to design the supporting structures.

The formula ( 3) and ( 4) represent the homogeneity H0 calculated with the maximum, minimum and average value of B inside of the working volume.

( ) ( )

0 max min 1 100 [%]   - =  - ×     H B B B ( 3
)
The homogeneity can be also calculated for each component in order to assure its quality using

( ) ( ) , , , , , , 0 ,
, max min 1 100 [%]   -   = - ×     x y z x y z x y z x y z B B H B (4)
The superconducting NbTi wire SC56S53 of Supercon Inc.

[9] was chosen for its high proportion of superconducting (SC) material (Cu:SC, 0.9:1). The diameter dwire of the wire is 0.54 mm. With 56 NbTi filaments, this wire can carry a high current density. Conversely, its thermal stability is less than that of wire with higher ratios, as the MR16 wire (Cu:SC, 11:1) dedicated for MRI imaging magnets requiring more stability. Moreover, to prevent the superconductors from transiting to the normal state, the magnet must have a load characteristic that does not exceed the characteristic of the wire. In other words, the current that flows in the magnet must never be higher than the critical current Ic of the SC wire. Thus, for a maximal value of the magnetic flux density Bmax over the wire, there is a maximal value of current I to be used. In order to avoid an accidental and devastating transition of the magnet, which can be very fast when using LTS wires, a 20 % margin from the Ic (|B|) curve is proposed as a criterion in Fig. 2. It is important to define the load line of the coil system. For this purpose, a parametric study has been conducted for different values of I and Ry.

B. Simulations

To ensure the right configuration of the magnetic calculation made by the electromagnetic calculation software, we have used an analytical calculation to verify the accuracy of the obtained results. This verification process was carried out using analytical equations for a Helmholtz coil. Then we have created a Comsol model for the entire magnet. Comsol is a 3D software for electromagnetic calculations using finite element method. The model was designed with the electromagnetic module and configured with the "multi-turn coil" function. This function allows to manage a large number of turns and ensures uniformity in the spatial distribution of the turns. The coils are placed in the air. A very fine mesh was used throughout the magnet, in order to minimize the potential errors that a bigger mesh would have caused. Unfortunately, this high level of mesh leads to a larger calculation time. Indeed, when the three Helmholtz coils are fed with current, there is no plane of symmetry, and this leads to a completely 3D problem, i.e. the magnetic flux density has three components that vary with the three space coordinates. In order to reduce the computational domain and time, we have considered each Helmholtz coil as an axisymmetric problem with the midplane as symmetry plane, and applying the principle of superposition, only one eighth of the whole geometry has been studied. Fig. 3 shows the magnetic flux density norm |B| in a 1/8 model of the vector magnet and the location of Bmax.

Within the constraints mentioned in Section II.A, we found solutions that meet with the specifications of the vector magnet. Then, the model parameters I and Ry will allow the implementation of the parametric study under the critical values Bc and Ic. The whole geometry of the model coil is built and automatically adjusted according to the value of the reference radius Ry. The method used for building the geometry prevents the coil overlapping each other, and ensures their separation.

The worst operating condition is obtained when the 3 Helmholtz coils are fed with the same current I. Thus, to find Bmax over the superconducting wire, it is necessary to make this study with all the coils serially connected. This value of Bmax is used to estimate the critical current Ic, and defines then the operating current that avoid the normal conducting transition of the magnet. In Fig. 2, with the value of Bmax about 8 T obtained during the simulation of the system with Ry = 15 cm, the value of the maximal operating current I is 110 A, which is 80 % of the critical current Ic.

In order to analyze the homogeneity of the magnetic flux density, a Comsol post processing calculation of ( 3) and ( 4) was performed.

III. RESULTS

Simulations with Comsol Multiphysics helped us for the determination of the critical values for the proper operation of the desired vector magnet. The relevant values to be determined are the maximal magnetic flux density Bmax over the superconducting wire, the homogeneity H0, Bx, By, Bz and |B| in the working space. In Fig. 4, the curves corresponding to Bx, By and Bz at Ry = 15 cm and I = 110 A are shown. The curve labeled Bx is the x-component of the magnetic flux density produced by the Helmholtz coil following the x-axis. By and Bz are labeled in the same way, while |B| is the magnetic flux density norm, produced by the whole magnet, plotted along the y-axis.

The maximum magnetic flux density Bmax over the vector magnet has been calculated using Comsol. The results for Bmax from the parametric study are shown in Table II. The values of the magnetic flux density norm |B| at the center of the vector magnet are presented in Table III. At a given current, e.g. 100 A, |B| does not change if the radius Ry is modified.

Table IV represents the total length of wire required to construct the vector magnet, its weight and the corresponding cost, assuming 0.6 €/m (June 2015).

The variation of the homogeneity which is computed in the entire working space, denoted as H0 vs. Ry, is plotted in Fig. 5.

IV. DISCUSSION

The simulation of the vector magnet has been accomplished with some constraints. In order to facilitate the direction control of B, the model ensures that all the components of B have the same value in the center of the working space, as shown in Fig. 4. Moreover, the magnetic flux density produced by each Helmholtz coil, fed with the same current, has to be of the same magnitude in the center of the vector magnet.

The sizes of vector magnet that are the most interesting for us correspond to a value of radius Ry set between 12 and 15 cm. For higher values, it would require a greater length of wire which not helps to improve significantly the homogeneity but increases quickly the cost, beyond our budget. Indeed, the only benefit of a bigger coil would be a better homogeneity in a given volume. As regards the cryostat of the vector magnet, a thickness of about 30 mm has been selected for the stainless steel wall of the cryostat and the thermally insulating layers inside. Radius Ry of 15 cm enables us to gain a sphere of 10 cm for a working space and a surplus of 2 cm between the cryostat and wire of the magnet to place the mechanical structure. That is why a radius Ry of 15 cm has been chosen.

The values of Bmax shown in Table II must be compared with the 80% Ic(|B|)) line of the wire in Fig. 2. If these values of Bmax exceed the 80% Ic(|B|)) line, the stability of the vector magnet is not guaranteed, and the values of Bmax and |B| are reported in a cell with a white background in Table II and Table III, respectively. As an example, at a current of 100 A, Bmax must not exceed 8.17 T in order to prevent the superconducting wire from transiting to the normal state. Since Bmax is always lower than 8.17 T at 100 A, all the values reported on the corresponding line of Table II have a gray background. Therefore, the parameters with a gray background in Table II and Table III are the only values that will ensure the proper operation of the vector magnet in the superconducting state.

As shown in Table III, a bigger Ry does not increase |B| at the center; however the Bmax value is increasing with increasing Ry. A way to change this tendency consist of increasing the cross section of the coils, i.e. the number of turns. Indeed, Bmax decreases if the cross section of the coils increases. It was unnecessary to increase the cross section of the coils in our case since a design of vector magnet that satisfies all our requirements has been found. The range of values from 100 A to 120 A for I, and from 10 cm to 19 cm for Ry are shown because the vector magnet can produce a magnetic field norm around the desired one of 3 T in this range.

Fig. 5 shows the variation of the homogeneity H0 vs. Ry. One can obtained a homogeneity of 98.87% for Ry = 15 cm. The curve shows the behavior of H0 as a function of the magnet size. If a larger homogeneity is required, the size of the magnet must be increased, leading irrevocably to an increase of its cost. The price of the cryostat for such a system is around 20 k€. By adding 16.1 k€ for the NbTi wire and approximately 9 k€ for the raw material (epoxy resin…) and the machining, the total cost of the designed vector magnet grows up to 45 k€. The values of weight represented in Table IV are useful to design the thermal performances of the cryostat.

V. CONCLUSION

The objective of this work was to study and to design a superconducting vector magnet which is able to produce a high, directional and uniform magnetic flux density to perform the characterization of HTS samples. The results demonstrated that a coil system consisting of Helmholtz coils can be adapted to these requirements with a total cost of around 45 k€. Due to the lack of iron in the structure, the behavior of the vector magnet is completely linear. The effects of an external magnetic screen in order to reduce the magnetic pollution around the magnetic system, i.e. the 5 gauss curve, remain to be examined. The presence of the three Helmholtz coils significantly reduces the performance of the wire unlike a single magnet. Since the design of the vector magnet is mostly influenced by Bmax, we can say that the Bmax value is the key variable in this design. The mechanical problems have not been addressed in the present article because work is ongoing. The uniformity of the magnetic flux density has a nonlinear behavior in the defined working zone. The configuration composed of three Helmholtz coils has a very good homogeneity. It fits well to the applications in electrical engineering for which 95 % of homogeneity is sufficient. Thus, we have designed an important device that will be very useful for the research in the superconducting field. This system is currently under construction and will be tested in the first semester of 2016.

Fig. 1 . 2 FFig. 2 .

 122 Fig. 1. Proposed vector magnet: system of 3 Helmholtz coils. The sphere located in the center represents the working space. TABLE I VECTOR MAGNET PARAMETERS Name Description Units Ic, Bc Critical current and magnetic flux density A, T I, B Operating current and magnetic flux density A, T Bx,y,z Magnetic flux density following each axis of the Helmholtz coil T Nx,y,z Number of turns of each Helmholtz coil tr Rx,y,z Helmholtz coil radius for each axis m H0 Homogeneity in the working volume % Bmax Maximal magnetic flux density over the coils T dwire Diameter of the superconducting wire mm Sx,y,z Cross-sectional area of each Helmholtz coil mm 2 F Filing factor: ratio of wire volume/cross section T = Tesla, A = Ampere, tr = turns, m = meters.

Fig. 3 .

 3 Fig. 3. Picture of the magnetic flux density norm |B| in a 1/8 model of the magnet and the location of Bmax.

Fig. 4 .

 4 Fig. 4. Values of Bx(x), By(y), Bz(z), and |B(y)| obtained with Ry = 15 cm and I = 110 A. As an example, the curve labeled Bx(x) is the x-component of the magnetic field produced by the Helmholtz coil following the x-axis.

3

  Axis length, l (cm) Magnetic flux density, B (T)|B (y)| B x (x) B y (y) B z (z)Location of Bmax

  .00 3.00 3.00 3.00 3.00 3.00 3.00 120 3.27 3.27 3.27 3.27 3.27 3.27 3.27 3.27 3.27 3.27 Magnetic flux density norm |B| at the center of the vector magnet for different values of Ry and I.

					TABLE III
			|B| AT THE CENTER OF THE VECTOR MAGNET
						Ry (cm)
	I (A)	10	11	12	13	14	15	16	17	18	19
	100 2.72 2.72 2.72 2.72 2.72 2.72 2.72 2.72 2.72 2.72
	110 3.00 3.00 3.00 3

TABLE IV TOTAL

 IV LENGTH, WEIGHT AND COST OF WIRE .4 17.1 20.1 23.3 26.8 30.5 34.4 38.6 43.0 m 12.3 14.9 17.7 20.8 24.1 27.7 31.5 35.5 39.8 44.4 Cost 7.1 8.6 10.3 12.1 14.0 16.1 18.3 20.6 23.2 25.8 Total length, weight and cost of the wire SC56S53 required to build the vector magnet for different values of Ry. The values are presented for cost and size comparison. Total length L = km, Weight m = kg, Cost = k€.Fig. 5. Variation of the homogeneity H0 vs. Ry.Radius of the Helmholtz coil following the y-axis, R y (cm)

	99 100			98,5	98,9	99,1	99,3 99,4 99,5 99,6 99,6 99,7
	Homogeneity, H 0 (%)	96 97 98	96,3	97,4	98,1
	Ry (cm) 10 L 11.9 1494,7 11 12 13 14 15 16 17 18 19 94 95 10 11 12 13 14 15 16 17 18 19 20 21 22


Nesrine Ouannes
email: nesrineouannes@gmail.com

Noureddine Djedi
email: n.djedi@univ-biskra.dz

Herve Luga
email: herve.luga@univ-tls1.fr

Yves Duthen
email: yves.duthen@univ-tls1.fr

Modeling a bacterial ecosystem through chemotaxis simulation of a single cell

Keywords: Artificial life, Bacterial chemotaxis, Single cell, Virtual ecosystem

Introduction

Artificial life research has made progress in the study of adaptive behavior through computational models of artificial organisms. Remarkably simple chemical reactions can perform movements toward attractants, and are in principle capable of modulating the behavior of artificial organisms. We will demonstrate in this work that a simple process of evolved bacterial chemotaxis can explain the emergence of more complicated behaviors seen in bacterial population dynamics. In terms of bacterial chemotaxis, this issue can be explored by looking at the influence of signaling network parameters on the spatiotemporal dynamics of bacteria that migrate toward chemical attractants and away from repellents. Cell chemotaxis is one of the simplest behaviors known, and it likely is one of the first behaviors to have existed in the history of life on earth. In bacterial chemotaxis such as that of Escherichia Coli (E. Coli), a bacterium typically tends to swim in a random walk, punctuated with periods of straight swimming sprints (or runs) interrupted by brief tumbles that causes the bacterium to change direction when the concentration of attractant or repellant is uniform or undetectable. In response to attractant gradient, this random walk becomes biased, and the bacteria tumble less frequently when encountering increasing attractant concentrations (i.e., they swim longer runs), and tumble more frequently when the attractant concentration decreases [START_REF] Adler | Chemotaxis in bacteria[END_REF].

Phosphorylation cascade in a chemotactic pathway was simulated by a lot of works, in [START_REF] Mello | Perfect and nearperfect adaptation in a model of bacterial chemotaxis[END_REF] a theoretical analysis of a full system of ODEs included phosphorylation Abstract We present in this paper an artificial life ecosystem in which bacteria are evolved to perform chemotaxis. In this system, surviving bacteria have to overcome the problems of detecting resources (or sensing the environment), modulating their motion to generate a foraging behavior, and communicating with their kin to produce more sophisticated behaviors. A cell's chemotactic pathway is modulated by a hybrid approach that uses an algebraic model for the receptor clusters activity, an ordinary differential equation for the adaptation dynamics, and a metabolic model that converts nutrients into biomass. The results show some analysis of the motion obtained from some bacteria and their effects on the evolved population behavior. The evolutionary process improves the bacteria's ability to react to their environment, enhancing their growth and allowing them to better survive. As future work, we propose to investigate the effect of emergent bacterial communication as new species arise, and to explore the dynamics of colonies. cascade, and other works focusing on chemotaxis bacteria [START_REF] Rao | Design and diversity in bacterial chemotaxis: a comparative study in Escherichia coli and Bacillus subtilis[END_REF][START_REF] De Gennes | Chemotaxis: the role of internal delays[END_REF]. More recently, a bacterial chemoreceptor model was developed in [START_REF] Nishikawa | Nonadaptive fluctuation in an adaptive sensory system: bacterial chemoreceptor[END_REF]. The present work presents a bacterial ecosystem where individual E. Coli's chemotactic pathway is simulated by modeling signal processing by mixed chemoreceptor clusters using a rapid equilibrium (algebraic) model, adaptation through methylation simulated by ordinary differential equations (ODEs), and E. Coli's characteristic flagellar motor-induced runs and tumbles [START_REF] Vladimirov | Dependence of bacterial chemotaxis on gradient shape and adaptation rate[END_REF]. The metabolism is modeled as a set of intracellular chemical reactions. The goal of this study is to demonstrate the importance of allowing to recycle dead cell matter released in the environment.

Bacterial chemotaxis

The chemotaxis process consists of three stages: chemoreception, signaling, and adaptation. Methyl-accepting chemotaxis proteins are located along the cell surface. These proteins act as chemoreceptors and bind with chemicals in the environment. If a nutrient attractant is detected outside of the cell, through MCP, the level of production of protein CheA decreases because the receptors state shifts to the off state. It has been shown that the activity of the receptor cluster depends on the local ligand concentration and the methylation level according to the MWC (Monod-Wyman-Changeux signal processing) model [START_REF] Endres | Chemosensing in escherichia coli: two regimes of two-state receptors[END_REF]8]. CheA binds with phosphate in the cell. The phosphate group is transferred from the active CheA to the response regulator CheY. The concentration of CheYp modulates the motor and its behavior makes the cell run or tumble.

MWC model

We applied the MWC model for a mixed receptor cluster [START_REF] Endres | Chemosensing in escherichia coli: two regimes of two-state receptors[END_REF]8]. The probability A that receptor cluster is active is dependent on ligand concentration and the methylation state of the receptors and calculated as where F = F on -F off , and F on/off is the free energy of the cluster to be on/off as a whole. Hence, the average activity per receptor in the cluster is A. The total free-energy difference in the mean-field approximation is F = n r f r (m), which is just the sum of the individual free-energy differences between the receptor on and off states.

where [S] is the ligand concentration, and K on/off r is the dissociation constant for the ligand in the on and off state, respectively.

(1)

A = 1/(1 + e F ), (2) f r (m) = f on r -f off r = ǫ r (m) + log 1 + [S ]/K on r 1 + [S ]/K off r ,

Adaptation model

Adaptation is modeled according to the mean-field approximation of the assistance-neighborhood (AN) model [START_REF] Endres | Chemosensing in escherichia coli: two regimes of two-state receptors[END_REF][START_REF] Hansen | Chemotaxis in escherichia coli: a molecular model for robust precise adaptation[END_REF]. It is assumed that both enzymes work at saturation ([CheR]= 0.16, [CheB]= 0.28) [START_REF] Li | Cellular stoichiometry of the components of the chemotaxis signaling complex[END_REF]:

The average methylation level evolves in time as

The parameter k indicates the adaptation rate relative to the wild-type adaptation rate V that is the rate of receptor methylation (see Eq. 3) [START_REF] Vladimirov | Dependence of bacterial chemotaxis on gradient shape and adaptation rate[END_REF].

Kinase activity

Both ligand binding and receptor methylation affect the activity of CheA. CheA kinase activity [START_REF] Vladimirov | Dependence of bacterial chemotaxis on gradient shape and adaptation rate[END_REF] is calculated as

(varrint into [0,1])
where A is the probability that receptor cluster is active, CheY tot is the total CheY concentration that is equal to 9.7 according to [START_REF] Li | Cellular stoichiometry of the components of the chemotaxis signaling complex[END_REF], and K A =5 and K Y =100 are the rate constants according to [START_REF] Vladimirov | Dependence of bacterial chemotaxis on gradient shape and adaptation rate[END_REF].

CheY phosphorylation

The concentration of CheYp is obtained as a function of active CheA from the steady-state equation [START_REF] Kollmann | Design principles of a bacterial signalling network[END_REF].

where CheY tot is the total CheY concentration, and CheZ is the total CheZ concentration, CheA is the active [CheA],

k y = 100 µM -1 s -1 , k Z = 30/[CheZ]s -1 , Y = 0.
1 are the rate constants according to [START_REF] Kollmann | Design principles of a bacterial signalling network[END_REF], [START_REF] Sourjik | Binding of the E. coli response regulator chey to its target measured in vivo by fluorescence resonance energy transfer[END_REF].

The CCW motor bias

The CCW motor bias depends on CheYp concentration in the following form [START_REF] Cluzel | An ultrasensitive bacterial motor revealed by monitoring signaling proteins in single cells[END_REF].

Bacterial metabolism

The metabolism is responsible for essential reproduction, growth, and development cycles, which are affected by gene expression and bacterial movement. An organism's genome contains instructions that encode the ability to metabolize one or more substrates present in the environment. Metabolism of a food either accelerates or

(3) dm/dt = a(1 -A) [CheR] -bA [CheB] (4) m(t + t) = m(t) + kV t (5) CheA = CheA tot AK A /(AK A + K Y CheY tot ), (6) CheY = CheY tot K Y A/(K Y CheA + K Z CheZ + g y ), (7) mb = mb 0 CheY(1 -mb 0 ) + mb 0
decelerates a bacterium's replication rate by a factor that is positive (nutrient) or negative (toxin), respectively. In this model, every bacterium is represented by a genome that encodes basic properties such as motion, energy absorption, toxin removal, and waste production. These properties are adjusted at each time step as the genome is expressed, and mutations are applied after each division. Forrest and Jones' simulation [START_REF] Forrest | Modelling adaptive systems with Echo. Mechanisms of adaptation[END_REF] allows for simple material cycling through absorption while agent is created and release agents die and their bodies deteriorate. We adopted a similar approach where the material within dying bacteria is released in the medium and available as a source of energy for surrounding bacteria. Metabolism is calculated as the organism's total energy (sum of constant basal energy and energy foraged via metabolized nutrients) and subtracted from the energy spent to tumble or run. This metabolic model allows bacteria to stabilize their energy consumption and eventually reach a splitting threshold. The bacterial energy cycle is implemented as follows:

-M t is the total metabolic expenditure (stepwise energy depletion rate); -M 0 is the basal metabolic level (initial level at birth, equal to 25); -A is the bacteria nutrient absorption rate (genetically encoded and described in the next section); -M F is the food source energetic value (2 units); -M W is the energetic value extracted from waste con- sumption (+1 unit); -M T is the metabolic cost of consuming toxins (-2); -M M is the metabolic cost of movement (-1 unit); -mb 0 is the tumble frequency as a result of bacterial chemotaxis; -M S is the metabolic cost of cell division (Mt/2).

Genetic representation

In E. Coli chemotaxis, flagellar rotation is generated by an intracellular system of moderate complexity modulated by sensing the medium in which the bacterium lives. This system exhibits properties such as sensing, adaptation, memory, and motor modulation. To adjust these properties in order to simulate bacterial population behaviors, we use a genome that encodes factors involved in various levels of the metabolic pathway responsible for chemotaxis. Each gene in the genome is composed of two different types of encoding as is shown in Fig. 1. The first type is a binary encoding that is used in four genomic loci. It is used to encode variable levels of the bacterium's ability to detect nutrients or toxins with the same receptor. The signals ( 8)

M t = (M 0 + A(M F + M W ) + M T + mb 0 M M + M S ) t
collected by these receptors constitute the environmental information fed to the chemical pathway generating chemotaxis. The values in M and T increase with the detected nutrient concentration. A is the bacterium's nutrient assimilation rate (consumption capacity), which affects its metabolism. The gene (AC) encodes a bacterium's sensitivity to detect autoinducer molecules from its own bacterial strain, allowing it to communicate with its own kind by secreting autoinducer molecules. The second type of encoding uses real numbers and is employed in the remainder of the genome. The following loci constitute the building blocks of the bacteria's chemotactic pathway: Cluster activity (CI), kinase activity, methylation level (Me), CheY phosphorylation (Y), and motor bias (Mb), which are employed in Eqs. (1-7).

Experimental results

The set of experiments presented here have been designed to answer the question of whether an agent-based simulation evolving simple bacterial chemotaxis can explain the emergence of more complicated behaviors at the population level. The work presented here makes a distinction between three types of chemical substances: (1) nutrients diffusing from multiple points in the environment, (2) nutrients released from dead bacteria (or waste), and (3) toxic compounds. Although the environment has been discretized, bacteria are free to move in the continuous two-dimensional space by translating their location. Each bacterium has one cell, with identical size, shape, and chemotactic network controlling its movements. The reader is invited to see [START_REF] Ouannes | Modeling a bacterial ecosystem through chemotaxis simulation of a single cell[END_REF] for the remaining parameters used in the setup of the chemotactic network. All runs presented in Fig. 2 (left) show a fast population growth in the first twenty simulation cycles, fueled by abundant nutrients and low levels of toxicity. This increase is followed by a high reproduction regime (or split operations) compensated by as many death events, which leads to a constant population size for about 200 generations. From generation 200 to 300, the sharp drop in population size is due to two factors: a high death rate caused by the depletion of food resources is compound with bad MCP and toxin avoidance abilities. From generation 300, and at every 300 cycles, the growth rate often increases and is the result of the combined effect of individually improved genomic traits via mutation. The number of species varies greatly during an experiment, indicating that a high population turnover is taking place due to bacterial death and division. It is important to keep in mind that the speed of the whole population dynamics is not solely due to food abundance in the environment or to population size, but is for the most part due to evolved skills. We monitored the changes in total population energy for all 30 runs [Fig. 2, middle], where resources are being consumed, resulting in a decrease in the population's collective energy at the beginning of the run. Since cell divisions happen at very high frequency, this process accounts for most of the energy depletion in the first 300 generations. The simulation parameters have been adjusted so that within thirty generations, after most of the population died of starvation, a few survivors remain. In iteration 300, when new nutrient resources are added to the environment, the bacteria feed on it and multiply until there is no more nutrients in the environment, at which point their energy decreases again but more slowly than before. After thousands of time steps, the surviving bacteria are more efficient energy consumers and are able to keep a more stable energy level. This effectiveness is due to improved nutrient detection sensitivity (MCP), and to a greater extent to improved consumption ability. We also present some examples of trajectories from a few simulated swimming bacteria in Fig. 2-right. Each bacterium responds to the presence of nutrient sources with long runs and short tumbles, or by random walks in the absence of detectable traces of nutrients (Steady State). Several characteristic angles formed between two periods of long straight swimming (a 'run') marked by a period of abrupt change in orientation (a tumble) can be seen in the figure. These successive angles illustrate how bacteria can progress from their initial positions toward a preferred, nutrient-rich zone. Two characteristic evolutionary phases have been observed in this system: in the beginning (20) of the first 300 cycles, the bacteria first execute long runs from a nutrient source to another, but switch their oscillation bias to transition to a random walk when all resources are depleted.

Discussion

The results shown above demonstrate that a simple simulation of single-celled organisms with biological mechanisms such as simple chemical reactions allows us to model more complicated emergent behaviors at the population level. We observed in the results a healthy growth rate (or increasing bacterial count) for several hundreds of epochs, as long as the resources are abundant in the environment, and the population's collective lifespan increases because the evolved bacteria consume less energy, thanks to their optimized metabolism and their increased ability at gathering resources. This system is thus favorable to the emergence (or adaptation) of efficient sensory apparatus to detect food and avoid toxins, allowing individuals to better avoid deadly situations and improve reproduction allowing the population to thrive. When bacteria are moving, consuming, and splitting, their chemotaxis network is optimized in order to control their movement.

Figure 3 explains how the chemotaxis network responds to intracellular and environmental changes. Four internal states have been isolated in this response: Steady State, Fully Inactive, Adapted, and FullyActive State. For each of these states, the protein concentrations have been measured and compared to [START_REF] Vladimirov | Dependence of bacterial chemotaxis on gradient shape and adaptation rate[END_REF]. The Steady State is characterized by the bacteria performing a random walk and exploring the environment with initial parameter values of CheA =0.0164, methylation = 1.92, CheY = 1.92, and motor bias = 0.65. When the bacteria detect food sources, they enter a consumption regime where transmembrane Fig. 2 Left-hand the 'Growth rate' runs, which we have replicated 30 times with quantitatively the same results, representing the optimal values of the whole of the bacteria for 5,000 steps. Middle the energy of the evolved population of bacteria for 30 runs at 5,000 genera-tions. Inside the figure, a zoom in the same runs for 500 generations. Right-hand path realized by some bacteria in 2D space for the first receptors sense the changes of attractant concentrations and become inactive. As the attractant binds to the receptors, it inhibits autophosphorylation activity of CheA. The CheY phosphorylated by the groups received from CheA (CheYp) diffuses to the flagellar motors and induces changes in motor rotation, causing the bacterium to run (move forward). The increase of attractant concentration switches the receptor to the 'off' state (i.e., Fully Inactive State), resulting in an initial fast decreases of kinase activ- ity (CheA) (to 0.002) and hence CheY level, and causes longer runs (i.e., mb 0.75). The decrease of ChA activity is followed by a slow CheR-dependent adaptation. In the Adapted State, both on and off states are equiprobable, and both enzymes CheR and CheB are working for methylation and demethylation simultaneously. In this regime, methylation increases receptor ability to stimulate CheA activity. A removal of attractants shifts the system to the on state (or Fully Active State) that activates CheA autophosphorylation (0.047) and hence the downstream CheY phosphorylation.

Conclusion and future work

Our model was designed to model growth and behavior of bacterial ecosystems by simulating a group of bacterial cells in a stepwise fashion. To analyze and visualize the evolved behaviors, we present data that characterize bacterial positions, energy, and state in the cellular reproduction cycle. These results demonstrate that bacteria are able to evolve through mutation. The constructed model of chemotactic E. coli employed a hybrid model for simulating chemical pathways, with a mix of algebraic, ODE, and stochastic components instead of a fully stochastic model with an evolutionary algorithm to evolve a population of bacteria.

In future work, we will build upon the existing chemotaxis network to generate more diversified behaviors and capable bacteria that can emerge as distinct species with different behaviors to eventually form colonies. This platform will allow us to test our model with different environmental and non-environmental conditions.

Fig. 1

 1 Fig. 1 The bacterium's genome

Fig. 3

 3 Fig. 3 Graph of variations in concentrations of proteins and enzymes used inside bacteria during the chemotaxis process


Olivier Cadot

Cyril Touzé

Arezki Boudaoud

Linear versus nonlinear response of a forced wave turbulence system

Linear versus nonlinear response of a forced wave turbulence system

Olivier Cadot and Cyril Touzé Unité de Recherche en Mécanique (UME), ENSTA-ParisTech, Chemin de la Hunière, 91761 Palaiseau Cedex, France

Arezki Boudaoud

Laboratoire de Physique Statistique, UMR 8550 du CNRS/ENS/Paris 6/Paris 7, 24 rue Lhomond, 75231 Paris Cedex 5, France A vibrating plate is set into a chaotic state of wave turbulence by a forcing having periodic and random components. Both components are weighted in order to explore continuously intermediate forcing from the periodic to the random one, but keeping constant its rms value. The transverse velocity of the plate is measured at the application point of the force. It is found that whatever the detail of the forcing is, the velocity spectra exhibit a universal cascade for frequencies larger than the forcing frequency range. In contrast, the velocity spectra strongly depend on the nature of the forcing within the range of forcing frequencies. The coherence function is used to extract the contribution of the velocity fluctuations that display a linear relationship with the forcing. The nonlinear contribution to the velocity fluctuations is found to be almost constant, about 55% of the total velocity fluctuations whatever the nature of the forcing from random to periodic. On the other hand, the nonlinear contribution to the fluctuations of the injected power depends on the nature of the forcing; it is significantly larger for the periodic forcing ͑60%͒ and decreases continuously as the randomness is increased, reaching a value of 40% for the pure random forcing. For all the cases of intermediate forcing from random to periodic, a simple model of the velocity response recovers in a fairly good agreement the probability density function of the injected power. The consequence of the existence of a linear-response component is discussed in the context of the fluctuation-dissipation theorem validation in experiments of out-of-equilibrium systems.

I. INTRODUCTION

Out-of-equilibrium dissipative systems with a large number of degrees of freedom still display open challenging questions on their statistical distributions of energy and energy fluxes. In the past decade, theoretical advances have emerged by taking into account the chaotic character of the microscopic time evolution ͓1͔, leading to important results such as the fluctuation theorem ͑FT͒ ͑or Galavotti-Cohen theorem͒ ͓2-4͔ or the pairing theorem by Dettmann and Morris ͓5͔. However, the extension of these results to turbulent large-scale and dissipative systems is still a matter of debate, and experimental validations do not also provide a clear picture. The first experimental validations of the FT were provided with systems having very few degrees of freedom ͓6-10͔ and with forcing being either continuous or periodic but never random. Second, as mentioned in ͓11,12͔, the apparent verification of the FT might be hindered by a too small range of explored fluctuation amplitudes or averaging times.

More recent experiments have achieved a correct statistical range, showing that the FT theorem is not fully verified; see, e.g., ͓13͔ for the fluctuations of power injection in randomly granular gases, ͓14͔ for gravitocapillary waves, ͓15͔ for forced turbulence on a free surface, and ͓16͔ for a randomly forced plate. Interestingly, all these experimental violations of the FT were found with a forcing having a random component to drive the system out of equilibrium. In the same context, the theoretical results obtained by Farago ͓17,18͔ for a Langevin linear model with external random forcing give analytical results for the departure from the FT predictions, so that the experimental results provided in ͓12͔ appear as a direct application. More generally, the nature of the forcing ͓19-21͔ is often overlooked in the investigation of turbulent systems because it is known to control the injected power statistics. For instance, in ͓20,21͔, two different statistics of injected power were achieved by a temporal modification of the forcing, whereas the local turbulent statistical properties were identical in both cases. The discussion above raises the question of the effect of the randomness of the forcing on the dynamics of such a system. We note, however, that this question was recently addressed for a thermal system ͓22͔.

In the present paper, a vibrating plate is set into a chaotic state of wave turbulence ͑WT͒ by a forcing having both periodic and random components. The vibrating plate provides a perfect setting to study energy fluxes and statistical distributions of global quantities. The WT regime is rather well established ͓23-25͔ and the persistence of waves has been experimentally verified ͓26͔. The plate displays a random superimposition of out-of-plane bending waves of broadband frequencies. The first attempt for verifying the FT has been reported in ͓16͔, showing a good agreement when the external forcing is periodic, whereas the random forcing leads to an important departure from the theoretical FT prediction. In that context, the natural question arising is that of the influence of the nature of the forcing on the response of the turbulent system and the statistical distribution of global variables.

The paper is organized as follows. Section II provides the details of the experimental setup, allowing a forcing that spans the continuum between the random forcing and the sinusoidal forcing. Section III provides a first decomposition for the velocity field into a turbulent part and a forcingcorrelated part, which is tested by predicting the probability density function ͑PDF͒ of the injected power. In Sec. IV, this decomposition is refined using the coherence function in the frequency space, which gives access to the different contributions to the fluctuations of the velocity and the inject power. In Sec. V, the two approaches are compared and the relevance of the results to the FT is discussed.

II. EXPERIMENT AND FORCING

The experimental setup is exactly the same as in ͓16͔. Briefly, a plate is locally forced by a magnet-coil system. The vibrating plate is a steel plate suspended by each of its corners to a rigid frame. The plate was chosen for its very high modal density, obtained by large dimensions, 2 m ϫ 1 m, for a thickness of h = 0.5 mm. Material properties were estimated as Young's modulus E = 200 GPa, Poisson's ratio = 0.3, and mass per unit volume = 7800 kg/ m 3 . The forcing device consists of a coil and a permanent magnet simply magnetized on the steel plate. It has been shown in ͓27͔ that in this configuration, the force F acting on the magnet is proportional to the current I͑t͒ circulating in the coil F͑t͒ = KI͑t͒. The current is measured by inserting an Ohmic resistance of 0.12 ⍀ in series with the coil. The proportionality constant found is K = 0.456 N / A. The normal velocity v at the application point of the forcing is measured with a laser vibrometer from Polytec ͑model OFV 056͒, as in ͓16͔.

The normal velocity v and the coil intensity I are simultaneously acquired at the sampling frequency of 5000 Hz. A force signal generated by a PC is amplified with a QSC audio RMX 2450 professional power amplifier which supplies the coil. The force signal is built from a two component function F ␣ Ј͑t͒ ͑where ␣ ͓0,1͔ represents the proportion of the first ͑periodic͒ component and prime symbol xЈ denotes reduced by standard deviation x ͒ as

F ␣ Ј͑t͒ = ␣ ␤ SЈ͑t͒ + ͑1 -␣͒ ␤ Ј͑t͒. ͑1͒
The first component is a sine function SЈ͑t͒ = ͱ 2sin͑2f 0 t͒ and the other, a band limited random function Ј͑t͒. The normalization coefficient ␤ = ͱ ␣ 2 + ͑1-␣͒ 2 ensures that F ␣ Ј = 1. The two time series SЈ͑t͒ and Ј͑t͒ are generated from a computer using LABVIEW. The sampling rate is 5000 Hz and the time duration 600 s. The frequency of the periodic part is f 0 = 75 Hz. The random part Ј͑t͒ has been generated by filtering a Gaussian white noise with a bandpass ͑5-75 Hz͒ filter. In order to have the current in the coil to be proportional to this force, the voltage signal delivered by the PC has been compensated by the transfer function between the amplifier input and the current circulating in the coil. This transfer function has been measured when applying a white noise as a voltage signal. In Fig. 1, we show the power spectra of the current ͑proportional to the force applied on the plate͒ for the two limiting types of forcing, F 0 Ј for pure random force ͑with ␣ =0͒ and F 1 Ј for periodic force ͑with ␣ =1͒. The random forcing is characterized by a very sharp bandwidth. For intermediate forcing ␣, this two spectra are combined following Eq. ͑1͒. We performed 11 experiments with ␣ = ͓0,1͔ by steps of 0.1. The magnitude of the force F has been adjusted to reach a wave turbulence state and re-mains fairly constant for all the experiments. For each experiment, we measured simultaneously the current circulating in the coil and the velocity of the forcing point with the laser vibrometer. For a given ␣, time series of a signal represents 2 400 000 points ͑approximately 600 s long͒. Auto and cross power spectra are time averaged over windows of 5 s ͑20 000 points͒, which gives a frequency resolution of 0.2 Hz. This time averaging is denoted by "͗ ¯͘" all the way through the paper.

III. CORRELATIONS IN THE REAL SPACE

A. Properties of injected power

We consider here the injected power

I͑t͒ = v͑t͒F ␣ ͑t͒. ͑2͒
The mean power ͗I͑t͒͘ and its standard deviation I are plotted versus ␣ in Fig. 2. They are normalized using the standard deviations of the force F and of the velocity v . The mean injected ͑or dissipated͒ power decreases slightly as the randomness of the force increases. The dependence upon the forcing detail is even more pronounced for the fluctuations I of the injected power. For the pure random forcing the fluctuations are 33% larger than for the pure periodic forcing. Our objective is to understand this strong dependence on the nature of the forcing and its consequences on the measured statistical quantities.

B. Decomposition of the velocity

As the power input is a measurement of the correlation between the velocity v and the force F ␣ , we first follow ͓16͔ and seek a decomposition of the velocity as the sum of two parts-one proportional to the forcing and one uncorrelated with the forcing. We write v = u + ᐉF ␣ . The condition that ͗uF͘ = 0 is achieved, say ͗vF ␣ ͘ = ᐉ͗F ␣ 2 ͘, determines the value of ᐉ and the velocity decomposition reads

v = u + r v F F ␣ , ͑3͒
where r is the correlation coefficient between v and F ␣ , and stands for the variance of the quantity . At this stage, Eq. ͑3͒ involves the hypothesis that the correlated response is proportional to the force, while any linear operator for ᐉ would achieve ͗uF͘ = 0. We will generalize to a linear operator in Sec. IV B.

This decomposition extracts the part of the velocity, the nonlinear feedback u, which does not contribute to the average input power ͗I͘ = ͗vF ␣ ͘. The correlation coefficient r = ͗I͘ / ͑ v F ͒ is shown in Fig. 2. It has a weak dependence on the forcing detail, with a maximum for the pure periodic forcing ͑r = 0.678͒ and a minimum for the pure random forcing ͑r = 0.589͒. Now we make the hypothesis that u and F ␣ are statistically independent variables in order to describe the probability density function of the resulting injected power. We further assume that u is Gaussian as it superimposes fluctuations at many scales. In Sec. III C, we investigate whether these hypotheses are consistent with experiments.

C. Probability distributions and predictions for injected power

Using the hypotheses introduced above, we aim at predicting the PDF of the injected power shown.

First, the PDFs of the force and the velocity are displayed in Figs. 3͑a͒ and3͑b͒. We can check in Fig. 3͑a͒ that as prescribed in Eq. ͑1͒, the force statistics change continuously from Gaussian to periodic as ␣ increases from 0 to 1. In contrast, the velocity response statistics in Fig. 3͑b͒ are not very sensitive to the forcing detail. They are rather Gaussian, with a slight negative skewness that increases as the randomness of the forcing increases ͑as ␣ → 0͒. The PDFs of the normalized injected power I ͗I͘ are plotted in Fig. 3͑c͒. We can observe a significant evolution from exponential tails when the forcing is random ͑␣ =0͒ to Gaussian tails when the forcing is periodic ͑␣ =1͒. These two limiting cases were previously studied in ͓16͔. We can also notice that the PDFs are not sensitive to ␣ when randomness of the forcing is strong as they are all superimposed for ␣ in the range of 0-0.4. We hereafter generalize the work in ͓16͔ where two models for the power PDFs were derived for the limiting cases ␣ = 0 and ␣ = 1. From Eqs. ͑2͒ and ͑3͒, the injected power reads

I = uF ␣ + r v F F ␣ 2 . ͑4͒
Then, introducing

q = r ͱ 1 -r 2 ,

͑5͒

the PDF of I can be computed as

f I ͑I ˜= I/͗I͒͘ = ͵ f F ␣ ͑F ␣ Ј͒f u ͑uЈ͒␦ͩI ˜-1 q F ␣ ЈuЈ -F ␣ Ј 2 ͪdF ␣ ЈduЈ,

͑6͒

where f F ␣ . and f u are the PDFs of the forcing and of the nonlinear velocity feedback, respectively. We assume this velocity feedback to be Gaussian, as suggested by the nearly Gaussian character of the total velocity ͓Fig. 3͑b͔͒,

f u ͑uЈ = u/ u ͒ = 1 ͱ 2 expͩ- uЈ 2 2 ͪ.

͑7͒

From force definition ͑1͒, its PDF can be computed as

f F ␣ ͑F ␣ Ј͒ = ͵ f S ͑SЈ͒f ͑Ј͒␦ͩF ␣ Ј -␣ ␤ SЈ - 1 -␣ ␤ ЈͪdSЈdЈ,

͑8͒

where

f ͑Ј͒ = 1 ͱ 2 expͩ- Ј 2 2 ͪ ͑9͒
is the PDF of the Gaussian noise and

f S ͑SЈ͒ = 1 ͱ 2 -SЈ 2 ͑10͒
is the PDF of the sinusoidal force. Therefore, Eq. ͑8͒ becomes

f F ␣ ͑F ␣ Ј͒ = ␤ ͑1 -␣͒ ͱ 2 3 ͵ -ͱ 2 ͱ 2 1 ͱ 2 -SЈ 2 ϫexp ͫ -1 2 ͩ ␤FЈ -␣SЈ ␣ -1 ͪ 2 ͬ dSЈ, ͑11͒
while replacing Eq. ͑7͒ in Eq. ͑6͒ yields the PDF of the injected power as

f I ͑I ˜= I/͗I͒͘ = ͱ 2 ͵ 0 ϱ qf F ␣ ͑F ␣ Ј͒ F ␣ Ј ϫexp ͫ -1 2 ͩ ␤F ␣ Ј -␣SЈ ␣ -1 ͪ 2 ͬ dF ␣ Ј . ͑12͒
Both integrals ͑11͒ and ͑12͒ are evaluated numerically.

1.E-08

1.E-06

1.E-04

1.E-02

1.E+00 1.E+01 1.E+02 1.E+03

1.E-08

1.E-06

1.E-04

1.E-02 It is worth noticing that this model only depends on the correlation coefficient r ͑measured experimentally͒ and ␣ ͑imposed by the forcing͒. Model ͑12͒ is plotted for three different values of ␣ = 0.9; 0.7; 0 in Fig. 4. It is found to match the PDF obtained experimentally when varying ␣ in the range ͓0,1͔. It is surprising to find such a good agreement while the PDFs of the velocity are not exactly Gaussian in the experiment as assumed in the model ͓see Fig. 3͑b͔͒; the injected power is not sensitive to the shape of the PDF of the velocity response but is dominated by the statistics of the imposed force. The case ␣ = 0 in Eq. ͑12͒ is identical to the model used in ͓14,16͔.

1.E+00 1.E+01 1.E+02 1.E+03 f (Hz) (a) (b) f 0 f m f 0 f (Hz)
The success of these predictions shows that it is fair to consider the velocity as a sum of two contributions: one linearly correlated with the forcing and one that is independent from the forcing. However, within this framework, the extraction of either contribution in the time domain would require a better knowledge of the linear response, which is at the moment simply approximated by a proportional model. For instance, a simple phase shift of the linear response in Eq. ͑3͒ would lead to different time series of u͑t͒, while the PDF results would remain unchanged. In Sec. IV we aim at identifying these two contributions in the frequency space.

IV. CORRELATIONS IN THE FREQUENCY SPACE

A. Velocity spectrum

The autopower spectra of the velocity ͉͗v ˜͑f͉͒ 2 ͘ are shown in Fig. 5. Notably, for frequencies lower than f 0 , the amplitude of the spectra depends on ␣, while for frequencies larger than f 0 , all the spectra are superimposed. It demonstrates that the high-frequency part is insensitive to the nature of the forcing. This observation suggests to decompose the energy of the fluctuations v 2 into three parts, each corresponding to different frequency domains,

v 2 Ϸ R 2 + P 2 + C 2 . ͑13͒
The first energy R 2 = ͐ f m f 0 -␦f 0 ͉͗v ˜͑f͉͒ 2 ͘df corresponds to the frequency domain of the random forcing: f m = 5 Hz to f 0 -␦f 0 = 74.8 Hz ͑excluding the peak at f 0 ͒. The second energy corresponds to the narrow domain of the periodic forcing: P 2 = ͐ f 0 -␦f 0 f 0 +␦f 0 ͉͗v ˜͑f͉͒ 2 ͘df with ␦f 0 = 0.2 Hz. The last energy C 2 = ͐ f 0 +␦f 0 +ϱ ͉͗v ˜͑f͉͒ 2 ͘df corresponds to the turbulent cascade domain. Looking at Fig. 5͑a͒, we observe that the spectra have similar shapes within the forcing frequency domain, which can be checked by a rescaling using the energy in the random frequency domain, R 2 ; in Fig. 5͑b͒, all power spectra superimpose satisfactorily within the frequency range ͓5 Hz, 75 Hz͔. The three contributions ͓Eq. ͑13͔͒ to the total energy are shown in Fig. 6 for ␣ ͓0,1͔. The sum is very close to v 2 because contributions for frequencies lower than f m are negligible. As could be expected from the shape of spectra for ͗f͘75 Hz in Fig. 5͑a͒, the contribution of the cascade domain is fairly independent of ␣.

B. Decomposition of the velocity using the coherence function

Here, we refine the decomposition of Eq. ͑3͒. We seek a decomposition of the velocity v in the form of the sum of a turbulent feedback u͑t͒ and a linear response

v͑t͒ = u͑t͒ + LF ␣ ͑t͒. ͑14͒
To do so, we define the linear operator L so that ͗u ˜͑f͒F ˜␣͑f͒ ء ͘ = 0, ͑15͒ where ء stands for the complex conjugate and ˜͑f͒ for the Fourier transform of the function evaluated at the frequency f. We then make the hypothesis that this operator is independent of the forcing. The self-consistent predictions made hereafter will give ground to this hypothesis.

A more direct way to compute L is to use the coherence function ␥͑f͒ that measures the correlation between the excitation F͑t͒ and the response v͑t͒ in the frequency space,

␥͑f͒ = ͗v ˜F ˜␣ ء ͘ ͱ ͉͗v ˜͑f͉͒ 2 ͉͗͘F ˜␣͑f͉͒ 2 ͘ . ͑16͒
By taking the average of the product of the Fourier transform in Eq. ͑14͒ by F ˜␣ ء and using ͗u ˜F ˜␣ ء ͘ = 0, we obtain

␥͑f͒ = L f ͱ ͉͗F ˜␣͑f͉͒ 2 ͘ ͉͗v ˜͑f͉͒ 2 ͘ , ͑17͒
where L f is the Fourier representation of the linear operator L defined as LF ˜␣ = L f F ˜␣. The modulus of the coherence, ͉␥͑f͉͒, is plotted in Fig. 7 for the frequency domain of the forcing, where the coherence has a meaning. We find that the coherence has a well-defined plateau with a magnitude that increases as ␣ decreases. The pure random forcing is then the most coherent forcing: the linear correlation between the response and the forcing is the strongest. When the periodic forcing is present ͑␣ 0͒, the coherence is always 1 at frequency f 0 . In the following, we use this coherence function to disentangle the different contributions to the fluctuations of velocity and power input.

C. Coherent and incoherent contributions to the velocity fluctuations

The definition of the operator L implies that the energy of the fluctuations of the velocity v is the sum of a nonlinear incoherent contribution and linear coherent contributions,

v 2 = u 2 + ͗͑LF ␣ ͒ 2 ͘. ͑18͒ The linear contribution L 2 = ͗͑LF ␣ ͒ 2 ͘ = ͗͑L f F ˜␣͒ 2 ͘ can be com- puted directly using Eq. ͑17͒, L 2 = ͵ 0 ϱ ͉␥͑f͉͒ 2 ͉͗v ˜͑f͉͒ 2 ͘df . ͑19͒
In order to probe the effects of the nature of the forcing, we further decompose the total linear contribution L 2 into two parts, L 2 = RL 2 + PL 2 , where

RL 2 = ͵ f m f 0 -␦f 0 ͉␥͑f͉͒ 2 ͉͗v ˜͑f͉͒ 2 ͘df , ͑20͒
corresponding to the frequency range of the random forcing, and

PL 2 = ͵ f 0 -␦f 0 f 0 +␦f 0 ͉␥͑f͉͒ 2 ͉͗v ˜͑f͉͒ 2 ͘df , ͑21͒
corresponding to the frequency of the periodic forcing. The nonlinear contributions are deduced by subtraction

RNL 2 = R 2 -RL 2 ͑22͒
for the frequency domain of the random forcing and

PNL 2 = P 2 -PL 2 ͑23͒
around the frequency of the periodic forcing. The total nonlinear contribution to the fluctuations energy of the velocity response in Eq. ͑18͒ is the sum

u 2 = NL 2 = PNL 2 + RNL 2 .
Figure 8 shows the results for the different contributions obtained with Eqs. ͑20͒-͑23͒. In the frequency range of the random forcing ͓see Fig. 8͑a͔͒, the fluctuations of the linear response RL 2 decrease as ␣ increases and tend to be zero for ␣ → 1 because the coherence function vanishes. The nonlinear response RNL 2 has a weaker dependence on the forcing detail, decreasing as ␣ increases. The decrease of RL 2 and RNL 2 with ␣ could be expected because the random part the force F ␣ decreases. In the frequency range of the periodic forcing ͓see Fig. 8͑b͔͒, the energy of fluctuations of the nonlinear response is almost zero ͑the coherence is close to 1 at f 0 ͒, while the energy of the linear response increases with ␣ following the weight of the periodic part of the force F ␣ .

It is worth noticing that the coherence function is exploitable only in the frequency domain of the forcing, say ͓5 Hz, 75 Hz͔. In the other domains, the measured coherence is not strictly equal to zero as it should be in the absence of forcing. This is due to the weak current induced in the coil by the magnet motion according to the vibration of the plate. Hence, the resonant frequencies of the plate are also present in the force signal but at a very low level. This can be seen in Fig. 1 on the spectrum of the periodic force F ␣=1 : the background noise below 10 -7 over the complete frequency range corresponds to the resonances of the plate that are observable on all the velocity spectra in Fig. 5.

D. Predicting the coherent contributions to velocity fluctuations

To push our analysis further, we now compute the linear response L in Eq. ͑18͒ using form ͑1͒ of the force F ␣ .

The total energy of the fluctuations of velocity response ͑18͒ is written in the Fourier space as

͵ 0 ϱ ͉͗v ˜͑f͉͒ 2 ͘df = ͵ 0 ϱ ͉͗u ˜͑f͉͒ 2 ͘ + ͉͗LF ␣ ͑f͉͒ 2 ͘df , ͑24͒
which using Eq. ͑1͒ becomes, for the power spectrum of v͑t͒,

͵ 0 ϱ ͉͗v ˜͑f͉͒͘ 2 df = ͵ 0 ϱ ͉͗u ˜͑f͉͒͘ 2 df + ͵ 0 ϱ F 2 ␣ 2 ␤ 2 ͉LSЈ ͑f͉͒ 2 df + ͵ 0 ϱ F 2 ͑1 -␣͒ 2 ␤ 2 ͉LЈ ͑f͉͒ 2 df . ͑25͒
Actually, the terms on the right-hand side belong to different frequency domains. The first term, which accounts for the turbulent feedback, spans all frequencies; the second term vanishes except around the frequency f 0 , while the third term is restricted to the bandwidth of the random forcing. Splitting these three terms into the frequency domains of the random forcing, of the periodic forcing, and of the cascade, respectively, we obtain

R 2 = ͵ f m f 0 -␦f 0 ͉͗u ˜͑f͉͒͘ 2 df + F 2 ͑1 -␣͒ 2 ␤ 2 ͵ f m f 0 -␦f 0 ͉LЈ ͑f͉͒ 2 df , ͑26a͒ P 2 = ͵ f 0 -␦f 0 f 0 +␦f 0 ͉͗u ˜͑f͉͒͘ 2 df + F 2 ␣ 2 ␤ 2 ͵ f 0 -␦f 0 f 0 +␦f 0 ͉LSЈ ͑f͉͒ 2 df , ͑26b͒ C 2 = ͵ f 0 +␦f 0 +ϱ ͉͗u ˜͑f͉͒͘ 2 df . ͑26c͒
These three contributions can be rewritten as sums of a linear part and a nonlinear part,

R 2 = a R ͑1 -␣͒ 2 ␤ 2 v 2 + RNL 2 , ͑27a͒ P 2 = a P ␣ 2 ␤ 2 v 2 + PNL 2 , ͑27b͒ C 2 = CNL 2 , ͑27c͒
where

a R = F 2 v 2 ͵ f m f 0 -␦f 0 ͉LЈ ͑f͉͒ 2 df ͑28a͒ and a P = F 2 v 2 ͵ f 0 -␦f 0 f 0 +␦f 0 ͉LSЈ ͑f͉͒ 2 df . ͑28b͒
The linear parts of the energy fluctuations, say RL 2 / v 2 = a R ͑1-␣͒ 2 / ␤ 2 for the random frequency domain and PL 2 / v 2 = a P ␣ 2 / ␤ 2 for the periodic frequency domain, are in good agreement with the measurements displayed in Fig. 8͑a͒. For each case, a R and a P are determined by a best fit ͑see caption of Fig. 8͒.

This agreement indicates the self-consistency of the decomposition of the velocity into a linear part and a nonlinear part ͓Eq. ͑14͔͒. We thus obtain a reliable measurement of the corresponding contributions to the energy of velocity fluctuations, v 2 = NL 2 + L 2 , as a function of the detail of the forcing ␣, where

NL 2 = RNL 2 + PNL 2 + C 2 , L 2 = v 2 -NL 2 . ͑29͒
Both contributions are plotted in Fig. 8͑c͒. It can be seen that the forcing detail ␣ does not affect much the relative weight of these contributions. The nonlinear response is always found to be slightly larger than the linear response, about 55% of the total-energy fluctuations.

E. Coherent and incoherent contributions to the fluctuations of the injected power

We now turn to the fluctuations of the injected power,

I͑t͒ = v͑t͒F͑t͒. ͑30͒
Similarly to the velocity, we can extract from the energy of the injected power fluctuation the contributions due to linear and nonlinear responses.

Using decomposition ͑14͒ of the velocity response, the injected power is rewritten as

I = uF ␣ + LF ␣ F ␣ , ͑31͒
which leads to both nonlinear and linear contributions,

I 2 = u 2 F ␣ 2 + ͗͑LF ␣ F ␣ -͗LF ␣ F ␣ ͒͘ 2 ͘. ͑32͒
The nonlinear contribution which is the first right-hand side term in Eq. ͑32͒ is calculated directly from the nonlinear contribution of the velocity response,

INL 2 = u 2 F ␣ 2 = NL 2 F ␣ 2 ,
while the linear contributions is simply deduced from

IL 2 = I 2 -INL 2 .
Figure 9 shows both linear and nonlinear contributions to the total fluctuation energy of the injected power. The nonlinear contribution INL 2 remains fairly constant whatever the forcing detail ␣ is. In contrast, the linear one, IL 2 , is very sensitive to the forcing detail and increases with the randomness of the forcing. Hence, for the pure periodic forcing ͑␣ =1͒, the dominant response in the power fluctuations is nonlinear, while for the random forcing ͑␣ =0͒ the dominant response is linear.

V. DISCUSSION

A. Relation between the two decompositions of the velocity

In the decomposition of the velocity following Eq. ͑14͒, the Fourier representation of the linear operator L could be computed using the coherence function. The first decomposition ͓Eq. ͑3͔͒ amounts to stating that

L = r v F ͑33͒
is a constant ͑r is the correlation coefficient between the force F and the velocity͒, which is the simplest approximation possible. We will refer to this decomposition as the proportional model. Using Eq. ͑33͒, we can compute the linear response in Eq. ͑26a͒ for the random frequency domain and in Eq. ͑26b͒ for the periodic frequency domain. The proportional model yields a P = a R = r 2 . The correlation coefficient r = ͗I͘ / ͑ F v ͒ was shown in Fig. 2. For the pure periodic forcing, the measured value of r gives a P = r 2 = 0.46, which matches exactly the value given by the best fit in Fig. 8͑c͒. When randomness is introduced, the best fit gives a R = 0.43 which is different from r 2 as predicted just above. Nevertheless, a R remains close to r 2 . All this means that the Fourier representation of the linear operator is actually close to a constant.

B. Linear response vs nonlinear response

The fluctuations of the response of the plate have been decomposed into a linear ͑coherent͒ component and a nonlinear ͑incoherent turbulent feedback͒ component.

͑a͒ For the velocity response, the ratio between these two contributions is fairly constant; the nonlinear contribution, NL , is always slightly larger than the linear one, L ͑see white symbols in Fig. 10͒. The proportional model gives a good estimation of this ratio, especially on the side of the periodic forcing ͑close to ␣ =1͒.

͑b͒ Concerning the injected power fluctuations, the ratio of nonlinear to linear fluctuations varies significantly with the forcing detail. The energy of fluctuations due to the turbulent feedback, INL , is clearly dominant for a periodic forcing while it is smaller than the linear part of fluctuations, IL , for a random forcing. Actually we can even observe in Fig. 9 that the variations of the total power fluctuation are completely ascribed to the linear part of the fluctuations.

This discrepancy between the velocity and energy input could appear for other global quantities characterizing a turbulent system because global quantities often reside in the frequency domain of the forcing and then mix linear and nonlinear responses. While the nonlinear response could be universal, the linear one is certainly not. For instance, a global quantity such as the injected energy on long times has been extensively studied in the framework of the Gallavotti-Cohen theorem ͑or fluctuation-dissipation theorem͒. One might wonder about the range of applicability of this theorem. An answer is given by the theoretical work of Farago ͓17,18͔ on the linear case of a particle submitted to viscous damping and a random Gaussian external force. He found that the conclusions of the theorem are not met. Experimentally, when a wave turbulence system is randomly forced Falcon et al. ͓14͔ and Cadot et al. ͓16͔ showed that results are very similar to those of Farago ͓17,18͔. In ͓16͔, the same system was forced periodically, and the conclusions of the theorem were met. The present work gives the explanation that the random forcing is dominated by the linear response that might invalidate the application of the theorem. One could therefore speculate that the FT only concerns the nonlinear part of the power fluctuations.

VI. CONCLUSION

The wave turbulence state of a vibrating plate was studied when the driving spans the continuum between a periodic force and a random force, but keeping the same rms fluctuations. Using the coherence function, we were able to disentangle the linear contribution to the velocity fluctuations from the nonlinear ones. Furthermore, a simple additive model consisting of a linear proportional part and a turbulent feedback allowed us to recover the statistics of the injected power.

The main experimental finding is that, when the forcing is random, the linear response dominates the energy of the fluctuations of the injected power, while in contrast, the nonlinear response dominates for the periodic forcing. This result may help understanding the observation ͓13-16͔ of the invalidation of the fluctuation theorem ͓2-4͔ for randomly forced systems, while it seems to be applicable to deterministically forced ones ͓16͔.

FIG. 1 .FIG. 2 .

 12 FIG.1. ͑Color online͒ Power spectrum of the force for the two limiting types of forcing: the periodic forcing ␣ =1 ͑black͒ and the random forcing ␣ =0 ͑red͒.

FIG

  FIG. 4. ͑Color online͒ Probability density function of the injected power divided by the mean. Continuous lines are the functions f I ͑I ˜͒ in Eq. ͑12͒ provided by model ͑3͒ for which there is no adjustable parameter.

  FIG.3. ͑Color online͒ Probability density function of ͑a͒ the force, ͑b͒ the velocity divided by their standard deviation, and ͑c͒ the injected power divided by its mean value versus ␣ spanning ͓0,1͔ in steps of 0.1. In ͑a͒ and ͑b͒, the thick dashed lines are the Gaussian statistics.

  FIG. 5. ͑Color online͒ Power spectra of the dimensionless velocity at the application point of the force for ␣ spanning ͓0,1͔ in steps of 0.1. In ͑a͒, the spectra are presented as ͉͗v ˜͑f͉͒ 2 ͘ / v 2 . In ͑b͒, the spectra are presented as ͉͗v ˜͑f͉͒ 2 ͘ / R 2 ͑see text͒. The frequencies f m = 5 Hz and f 0 = 75 Hz stand for the bounds of the forcing in the frequency space as shown in Fig. 1.

  FIG.8. Decomposition of the velocity fluctuations using the coherence function. ͑a͒ Linear ͑ RL 2 ͒ and nonlinear ͑ RNL 2 ͒ contributions to energy of the fluctuations in the frequency domain of the random forcing R 2 . ͑b͒ Linear ͑ PL 2 ͒ and nonlinear ͑ PNL 2 ͒ contributions to energy of the fluctuations in the frequency domain of the periodic forcing P 2 . Dashed lines are a R ͓͑1-␣͒ 2 / ␤ 2 ͔ with a R = 0.43 in ͑a͒ and a P ͑␣ 2 / ␤ 2 ͒ with a P = 0.46 in ͑b͒, see text. ͑c͒ Linear and nonlinear contributions to the energy fluctuation of the velocity response. All contributions are drawn as a function of the forcing nature as defined by ␣.

  FIG. 9. Linear IL and nonlinear INL contributions to the injected power fluctuations I , normalized using the standard deviations of the velocity v and force F .

  FIG. 10. Ratio of the nonlinear to linear response for the velocity and the injected power. Proportional model ͑3͒ for the linear response ͑dashed lines͒ is compared to the experimental data using the coherence decomposition ͑symbols͒.


Marwa Chafii
email: marwa.chafii@supelec.fr

Lamarana Diallo

Jacques Palicot

Faouzi Bader

Rémi Gribonval

M Lamarana Diallo

Adaptive Tone Reservation for better BER Performance in a Frequency Selective Fading Channel

Keywords: Tone Reservation (TR), Peak-to-Average Power Ratio (PAPR), Orthogonal Frequency Division Multiplexing (OFDM)

published or not. The documents may come

I. INTRODUCTION

Multicarrier modulation (MCM) systems are widely used in wireless and wireline communication applications; for example, discrete multi-tone modulation is used in the asymmetric digital subscriber line (ADSL) standard, and orthogonal frequency division multiplexing schemes (OFDM) are used in terrestrial digital television, WiFi, and in the long term evolution (LTE). MCM systems have distinct advantages over single-carrier modulation systems, especially in a frequency selective fading channel, but they suffer from high fluctuations of the signal envelope. The efficiency of the high power amplifier in the transmission side, is then decreased because of these large variations. The peak-to-average power ratio (PAPR) is the random variable usually used to measure the variations of the MCM signals.

There are several PAPR reduction techniques such as partial transmit sequences approaches, selective mapping approaches [START_REF] Bauml | Reducing the Peak-to-Average Power Ratio of Multicarrier Modulation by Selected Mapping[END_REF], clipping and filtering techniques [START_REF] Li | Effects of Clipping and Filtering on the Performance of OFDM[END_REF], and tone reservation (TR) techniques [START_REF] Tellado-Mourelo | Peak to Average Power Reduction for Multicarrier Modulation[END_REF]. The latter family is one of the most popular schemes and it consists in using a set of reserved tones to design a peak cancelling signal for PAPR reduction. As TR works on reserved tones, no additional signal processing is required at the receiver side to extract the data information. Besides, the computation of the peak cancelling signal can be achieved thanks to an efficient convex optimization algorithm such as the second order cone programming [START_REF] Zabre | SOCP Approach for OFDM Peak-to-Average Power Ratio Reduction in the Signal Adding Context[END_REF], or a gradient project algorithm [START_REF] Tellado | Efficient Algorithms for Reducing PAR in Multicarrier Systems[END_REF]. Due to these reasons, TR is quite popular for practical implementations and therefore seems to be a promising technique for commercial standards such as Digital Video Broadcasting-Terrestrial (DVB-T2). In the literature, there are several selection methods of the reserved tones in the TR reduction technique. The choice of the reserved tones can improves the TR performance as well as the performance of the system.

In this paper, we propose a new method to select the reserved tones in the TR method, in order to improve the bit error rate (BER) performance of the system, in the presence of a frequency selective fading channel. The selection of the reserved tones is adapted to the channel characteristics. Under the assumption that the channel state is perfectly known at the transmitter, the reserved tones are the ones that are the most affected by the channel. For the classical selection methods, the reserved tones are used to reduce the PAPR, their utility is then limited to this function. Our new selection technique make the reserved tones also useful for improving the BER performance of the scheme.

II. ADAPTIVE TONE RESERVATION

A. Notations and definition

Throughout this paper, an OFDM symbol x(t) of duration T u is used and expressed as follows:

x(t) = M -1 m=0 X m e j2πmF t , 0 ≤ t ≤ T u , (1)
where M is the total number of carriers, F = 1 Tu is the intercarrier space, mF the m-th frequency and X m the input symbol carried out by the m-th carrier. After oversampling the signal by a factor L, x l such that l = 0, . . . , LM -1, are the discrete time domain samples at the instant lT e where T e = Tu LM and x = [x 0 , . . . , x LM -1 ] is the discrete OFDM symbol (vector). We will use capital letters to denote Fourier transform.

B. Adaptive Tone Reservation (ATR) principle

Instead of reserving the tones randomly, or according to a predefined set known at the transmitter and the receiver as in DVB-T2, they can be selected differently in order to improve the BER performance of the system. Our proposed method consists in selecting the carriers most affected by the frequency selective fading, and in using them to reduce the PAPR, rather than carrying useful information. Since these carriers will be discarded in the receiver side, the useful information will be carried in the less affected carriers, which helps improving the BER performance.

Thus, every frequency selective fading channel h, the M R tones that correspond to the smallest |H(f )| are reserved for the PAPR reduction process, according to Algorithm 1. We assume that there is a feedback of the channel state information from the receiver to the transmitter, and the channel is stationary during this period of time.

In order to understand this new selection solution, we consider in the first place a specific channel. Table I defines the frequency selective fading channel of Proakis C [START_REF] John | Digital Communications[END_REF]. The considered number of carriers, M , is equal to 1024. Figure 1 shows the amplitude frequency response of the channel. Based on the observation of this figure, it can be seen that the channel has a spectral null around carrier indices 396 and 630. The frequency coefficients around these indices are very small, and the relative noise around them is then very high.

By applying the new proposed selection method, the M R reserved carriers (M R is equal to 12 is this case) are selected as described in Algorithm 1, and they are: Let R be the set that contains the indices of the M R reserved tones. In the ATR method this set is obtained thanks to Algorithm 1, and in the TR approach this set is predefined as in the DVB-T2 standard. With this set of reserved tones, the main goal of the TR, and also the ATR approach is to find the effective peak cancelling signal in order to reduce the PAPR. Several analyses are proposed in the literature for this purpose [START_REF] Tellado-Mourelo | Peak to Average Power Reduction for Multicarrier Modulation[END_REF], [START_REF] Zabre | SOCP Approach for OFDM Peak-to-Average Power Ratio Reduction in the Signal Adding Context[END_REF], [START_REF] Tellado | Efficient Algorithms for Reducing PAR in Multicarrier Systems[END_REF]. In this paper, the TR method that uses the gradient project algorithm will be used to compute the signal to be added [START_REF] Tellado | Efficient Algorithms for Reducing PAR in Multicarrier Systems[END_REF]. This method consists in maximizing the Signal to Clipping noise power Ratio (SCR). Let

Algorithm 1 Get Reserved Tones

f (r) = Ae j arg(r) , if |r| ≥ A r, else, (3)
denotes the clipping function, with A a predefined threshold. Adding corresponds to replacing x by x + F R C where F R is the sub matrix of the inverse discrete Fourier transform matrix containing only the column indexed in the set of the reserved tones. The SCR to be maximized is defined as:

SCR = x 2 2 x + F R C -f (x + F R C) 2 2
.

Maximizing ( 4) is equivalent to minimizing its denominator by solving the following optimization problem min

C x + F R C -f (x + F R C) 2 2 . (5)
Starting with the initial condition c (0) = 0 MR , the peak cancelling signal can be computed iteratively thanks to the gradient project algorithm featuring the following equation [START_REF] Tellado-Mourelo | Peak to Average Power Reduction for Multicarrier Modulation[END_REF]

c (i+1) = c (i) -µ |x (i) l |>A x l + c (i)
l -Ae j arg(x l +c (i) l ) P l , [START_REF] John | Digital Communications[END_REF] where P l represents the l-th line of matrix P = F R (F R ) H . Therefore, by adding the useful signal x in the two terms of (6), we obtain the reduced PAPR signal after the (i + 1)iterations as the following

x (i+1) = x (i) -µ |x (i) l |>A x l + c (i) l -Ae j arg(x l +c (i) l ) P l , (7)
where x (i+1) = x + c (i+1) . Note that, in the TR approach, because of the set of reserved tone positions is predefined (is the same for each x), the matrix P = F R (F R ) H can be precomputed and stored. Therefore, the numerical complexity of the TR method is O(M L). In the ATR method, since the set of reserved is not predefined (it depends on the channel), the ATR method has the same complexity as the matrix

P = F R (F R ) H .

III. RESULTS

Simulations in this Section have been performed for OFDM with a number of carriers equal to M = 1024, and a number of reserved tones equal to M R = 12 using 4QAM, 16QAM and 64QAM for BER performance, and 4QAM for PAPR analysis. The inter-carrier spacing is considered equal to ∆F = 5kHz.

A. BER Performance

We evaluate the BER performance (in a frequency selective fading channel) of the proposed adaptive tone reservation (ATR) method based on adaptive tone selection, and compare it with a classical TR method which consists in reserving the first M R tones. The channel delay τ and the average path gain G av that characterize the channel considered in these simulations are defined in Table II, which represents the extended typical urban (ETU) model used in the LTE multipath channel standard [START_REF][END_REF]. The gradient project algorithm is used in the TR simulations, with a number iterations of N = 20, stepsize µ = 0.25, and threshold A corresponding to 6.5 dB. The BER and the PAPR simulations are performed using n 1 = 10 4 and n 2 = 10 6 iteration loops respectively. Figure 2 gives the BER performance of the OFDM system based on the proposed ATR with adaptive tone selection, and compares it with the TR method with specific tone selection. As depicted in Figure 2, for a BER of 10 -3 , the proposed method gives a gain of 5 dB in terms of signal-to-noise ratio (SNR) for the different constellations used in the simulation (4QAM, 16QAM, 64QAM). This can be explained by the fact that the carriers affected by the deep fading of the channel do not carry useful information, and they are discarded in the receiver side. Therefore, they do not affect the BER of the whole system.

B. PAPR Performance

The proposed selection technique improves the BER in frequency selective fading channel, without reducing the PAPR performance compared with the classical method. In order to confirm this claim, we simulate the PAPR performance using the ATR reduction method based on the adaptive selection method and we compare it with a specific selection method which consists in reserving the first M R tones. The PAPR performance is evaluated using the complementary cumulative distribution function which is the probability that the PAPR exceeds a certain value γ. Figure 3 confirms that there is no degradation of the PAPR when the new selection method for ATR is used. In this paper, we have proposed a new selection method for adaptive tone reservation. It has been shown through simulations that this new proposed selection technique allows 5 dB gain in terms of signal-to-noise ratio for a bit error rate of 10 -3 , for different constellations, in frequency-selective Rayleigh fading channel for uncoded OFDM system, without any degradation of the PAPR performance compared with a classical selection method.

The future work is about evaluating the gain in terms of BER performance for Coded-OFDM.

Figure 1 :

 1 Figure 1: Amplitude Response of Proakis C Channel.

  Require: H, M R positions ← [ ] Sort the elements of the channel H in ascending order and return their indices [-, indexes] ← sort(H) Select the M R first indexes that correspond to the M R smallest values of H for i ← 1 : M R do Return the positions of the reserved tones which are the most affected by the channel selectivity positions ← [positions, indexes(i)] end for

Figure 2 :

 2 Figure 2: BER Performance of the Uncoded OFDM based on the Adaptive Tone Reservation for different constellations.

Figure 3 :

 3 Figure 3: PAPR Performance Comparison between OFDM based on the Adaptive Tone Reservation and OFDM based on Classical Tone Reservation.

Table I :

 I Coefficients of Channel C of Proakis.

	h0	h1	h2	h3	h4
	0.027	0.460	0.688	0.460	0.027

Table II :

 II Channel Delay and Power Profile of ETU Channel.

ACKNOWLEDGMENT

This work has received a French state support granted to the CominLabs excellence laboratory and managed by the National Research Agency in the "Investing for the Future" program under reference Nb. ANR-10-LABX-07-01.