Adrien Lasbouygues

Benoit Ropars

Robin Passama

David Andreu

Lionel Lapierre

Atoms Based Control of Mobile Robots with Hardware-In-the-Loop validation

Mobile robots are nowadays used in a wide variety of missions especially for exploring environments and allowing scientists who study these environments to gather data about them. However the monolithic design of control laws is often a hindrance to adapting control solutions from a robotic application to another. We thus propose a new control description paradigm strongly focused on modularity as well as integrating constraints coming both from the system capabilities and from the concerns of control engineering (such as stability). In this paper, we present the main ideas behind our approach and illustrate it through Hardware-In-the-Loop simulation.

I. INTRODUCTION

The first mobile robots were developed in the late 60's. Their development has never stopped since leading to more and more autonomous robots capable of performing missions in more and more complex environments. Thus they began interesting other sectors of the scientific community as a tool to explore environments too dangerous or too costly for human exploration. This is especially true in our application context which is underwater robotics (for karstic exploration, archeology or biology) in which we exploit applicationdependent models of the environment.

However a broader use of robots in scientific applications is largely impeded both by the time and cost of developing a control for a given robotic application (i.e. the resulting interaction between a robot, its environment and the task performed) and by the difficulty to integrate knowledge from a specialist of the environment to this application.

The main underlying reason behind these issues is the way control laws are designed by control engineers. Indeed as explained in [1], control engineers focus on solving the control problem at hand often neglecting other considerations. Thus the resulting control laws are very monolithic, merging knowledge together while strongly hiding relationships between them as well as design assumptions.

Moreover control engineers rarely take into account the fact that the control will then be implemented on a Software and Hardware architecture which brings its own constraints and limitations (for instance sensors update rates or computation power). The monolithic nature of control equations also strongly limits software engineers leeway when implementing the control thus not allowing them to

The authors would like to thank the French Labex NUMEV for supporting this research. 1 Adrien Lasbouygues, Benoit Ropars, Robin Passama, David Andreu and Lionel Lapierre are with LIRMM, University of Montpellier, Montpellier, France. E-mails: {adrien.lasbouygues, ropars, robin.passama, andreu, lapierre}@lirmm.fr 2 Benoit Ropars is with Ciscrea, Toulon, France. E-mail: bro@ciscrea.fr use the full scope of their proposals in term of software modularity and architectural solutions.

We thus believe that a strong focus should be put on modularity concerns from the very beginning of the control design and up to the implementation of this control. This is the guideline of the approach presented in this paper.

We are currently working on the modified version of the Jack ROV which was presented in [START_REF] Ropars | Thrusters dead-zones compensation for the actuation system of an underwater vehicle[END_REF]. This paper is organized as follows. In section II, we discuss related works as well as our positioning towards them. In section III, we describe the basic entity of our approach. In section IV, we describe how to compose these entities together to describe a control law. Section V presents an insight on how Constraints are studied and show their importance through a rather simple heading angle control example which is validated using Hardware-In-the-Loop (HIL) simulation. Finally this paper ends with a conclusion and an opening towards future work.

II. RELATED WORKS AND POSITIONING

A. Modularity in control engineering

If modularity is a recent concern for control engineers, a few approaches already tried to take it into account.

Interesting approaches are the ones based on hybrid control such as [START_REF] Toibero | Stable switching control of wheeled mobile robots[END_REF]. Indeed modularity is achieved by partitioning a complex control problem into simpler subproblems and to develop an individual solution for each of these subproblems. The complete solution is obtained by switching between the subproblem solutions according to the situation faced by the robot. However implementation on a software architecture is not taken into account in this approach.

Another interesting approach is Motion Description Language (MDL) introduced in [START_REF] Brockett | On the computer control of movement[END_REF] which was then extended in works over MDLe [START_REF] Hristu-Varsakelis | The MDLe Engine: a software tool for hybrid motion control[END_REF]. Focused on hybrid control, these approaches describe a control as a sequence of kinematic state machines with a lifetime T that can be aborted by an interrupt. Such an entity is called Atom or Modal Segment and if we use the same name for our basic entities (since in both cases they are "minimal" entities composed to form more complex ones), our definition of an Atom is very different as exposed in section III.

Another noteworthy work is the Modelica modelling language. Though Modelica primary concern is not control law design, it presents relevant concepts that could be extended to control design. Modelica is a language aimed at simplifiying the exchange of models and model libraries used for modelling complex systems [START_REF] Elmqvist | Modelica -the new object-oriented modeling language[END_REF]. It builds on the benefits of object-oriented methodology and applies its principles for general-purpose physical modelling thus allowing to model complex systems by assembling components [START_REF] Elmqvist | Modelica -the new object-oriented modeling language[END_REF].

Modelica was then given the possibility of describing controllers through the addition of several libraries such as the ones presented in [START_REF] Baur | Modelica libraries for linear control systems[END_REF] and [START_REF] Bonvini | A modelica library for industrial control systems[END_REF]. But Modelica approach is not focused on control description and so does not encompass issues related to its real-time implementation.

B. Modularity and software engineering

On the other hand, modularity has long been a major concern for software engineers especially when designing robot control architectures. It has led to the development of several modular architectures such as CLARAty [START_REF] Volpe | CLARAty: Coupled layer architecture for robotic autonomy[END_REF]. One can also cite GenoM [START_REF] Ceballos | GenoM as a robotics framework for planetary rover surface operations[END_REF], an approach based on the generation of modules in order to integrate control algorithms to allow their composition and implementation with no precise information about their content. One can also refer to robotic middlewares such as ROS [START_REF] Quigley | ROS: an open-source robot operating system[END_REF]. All these works aim at developing modular control softwares and control modularity and requirements are not explicitely adressed.

In [1], a layered architecture was also proposed to improve modularity on sensing and actuating steps. Thus a Sensor Fusion layer and an Actuator Fusion layer are added to separate concerns between the control design and the physical system drivers. Moreover requirements and capabilities of the system are expressed through the interface of the components which encompasses data's unit, value, valid range as well as measurement uncertainties. However modularity of the control layer is not taken into account.

In [START_REF] Derler | Cyber-physical system design contracts[END_REF], design contracts based on temporal constraints and delays are proposed to organize the interaction between control and software engineers. These contracts reify the necessary tradeoff between control and software engineers when implementing a control. However the design process in itself is not detailed.

C. Positioning

As we saw in previous subsections, no work tried to propose a comprehensive approach encompassing the whole control design process from its description to its implementation while putting a strong and effective emphasis on modularity. Moreover, while most approaches such as [1] propose a separation of concerns between control engineering and software engineering, we are considering the integration of both concerns throughout the design process aiming at a better integration of a control on the software architecture used to implement it.

Figure 1 presents an overview of the proposed methodology. On the one hand of the process, the control engineers describe the control from a composition of basic entities called Atoms which encapsulate knowledge coming from various sources. We will give a more precise description of this part in sections III and IV.

On the other hand of the process, the software engineers will implement the control on the Software and Hardware architecture of the implementation target using the projection of Atoms as software entities using a dedicated API. This part will not be specifically adressed in this paper.

However the frontiers between these different domains are rather porous since expertise of one domain specialist might also be required in the other one.

Thus we introduce between them a "Bridge" which, taking advantage of the precise control description allowed by Atoms, will allow us to study our control properties to check if the control can be implemented on our target system. And if it can, this study allows us to produce information that will help the software engineer in his task. This part will be presented in more details in section V. As we want to describe the control of a robot in a modular way, we have to determine precisely which criteria our modular approach should satisfy in order to design the entities which will be the base of our control design paradigm in the most suitable way.

A. Knowledge in control engineering

First of all we should define what we mean by the term knowledge.

We define as knowledge any relevant information allowing the robot to perform its mission online as well as any data that could be gathered to allow offline post mission processing.

This knowledge can take heterogenous forms (parameters, models, these models being maps, equations, ...). Moreover knowledge comes from various sources. For instance control engineers bring control equations, environment specialists provide their models of the environment, robot manufacturer knows the robot's hardware specifications. Thus we would like that control engineers who design the control equations and software engineers who then implement them to be able to handle this knowledge without being familiar with it.

Finally knowledge is only the expression of the current set of information we know about a specific point (this is particularly true for environment knowledge).

Thus knowledge is evolving as we gather data through experiments or when theoretical breakthroughs are achieved. We thus would like to be capable of making models evolve without having to redesign our whole control laws.

We should then design the control no longer as a monolithic block but as a composition of entities. These entities shall encapsulate knowledge to allow it to be manipulated easily without requiring information about their precise content and allow to modify their content without affecting the surrounding entities. We are thus required to design a self-sufficient interface to prevent the user from making any assumptions that could lead to composition errors.

Finally these entities should also convey the information that would be relevant for a software engineer to implement this control within his software architecture on a given hardware architecture while respecting the requirements from the control engineer.

B. Atoms: the basic entity

Based on previous statements, we created an entity named Atom to support our new control description paradigm. These Atoms are mainly inspired from the object paradigm and are in their structure similar to Modelica objects. Moreover we paid a particular attention on their interface as it will be discussed in the next subsection.

They are defined as: Definition 1: An Atom is a minimal and indivisible composable entity encapsulating a piece of knowledge. An Atom, A, is represented as a 5-tuple: A = (Name, Phy, Int, IntPar,Constraints)

The Inerface parameters, IntPar, will be defined in the next subsection and Constraints will be presented in section III-D.

The interface, Int, is what allows an Atom to be connected to other Atoms. It is a 3-tuple (Ne, Pr, IS):

• The Needs (Ne): the set of knowledge required by the Atom. • The Products (Pr): the set of knowledge produced by the Atom. • The Internal Storages (IS): the set of knowledge stored by Atoms (for data historization or integration of local states for instance). Any knowledge item belonging to one of these sets is called an interface element. It is a pair (Name, Datatype). The Datatype of an interface element i will be noted as D(i) and will be defined more precisely in the next subsection.

The Physics, Phy, is the heart of the Atom. It is an application defined as:

Phy : D(Ne 1 ) × ... × D(Ne Dim(Ne) ) × D(IS 1 ) × ... × D(IS Dim(IS) ) → D(Pr 1 ) × ... × D(Pr Dim(Pr) ) × D(IS 1 ) × ... × D(IS Dim(IS) )
However there exists a specific kind of Atoms called External Knowledge Atoms characterized by an undefined Physics. They are used to explicit the elements that are not directly part of the control composition but either produce knowledge used by the control or use knowledge produced by the control (for instance sensors, actuators, network links, files to store data). Being implementation dependent, their Physics is not defined and they will be projected as software entities in a different way Atoms are.

We must note that the Atom paradigm, on the opposite of other object-oriented approaches such as Modelica, doesn't integrate inheritance mechanisms.

C. Interface characterization : the Datatypes

As mentioned previously, the focus on the interface definition is crucial for allowing an efficient use of Atoms. Thus the Datatype of the interface elements should convey all the required information to use the Atom.

The Datatype is a 3-tuple (Type, Frame, Axis). It contains the Type of the exchanged data which can be either a physical quantity (i.e. Force, Speed, Acceleration and so on) or an arithmetic type (i.e. Boolean, Unsigned Integer and so on). While in other approaches such as [1] or [START_REF] Elmqvist | Modelica -the new object-oriented modeling language[END_REF], interfaces express the unit of the data, in our approach, the Type is an object that encapsulates the data and offers different getters and setters (in the way they are defined in object methodology) for the different units associated with the quantity. Thus the element is specified in terms of quantity and no longer in terms of unit.

But can we consider that knowing the Type is sufficient ? • The control is performed in sonar frame.

• The sonar measurements are converted to robot frame in the sensor's driver. • The sonar measurements are converted to robot frame in the control Atom. • The sonar is mounted in such a way that sonar and robot frames coincidate and thus no conversion is required. Without more information we can't sort out which one is the right assumption. And when we have to perform this kind of assumptions at numerous locations when reusing a control description, there is little chance our implementation will work and finding which error(s) were made on the interpretation of the control would be a very difficult task. This is why in addition to the type, Atom Datatypes also encompass a Frame field used to indicate in which frame is expressed the data. In the previous example, we indicate that the Product of the sonar is expressed in the PROXIMET ER frame while the Need of the control is expressed in the ROBOT frame. Thus direct connexion between them is impossible and we are forced to add an Atom in charge of performing a frame shift from the PROXIMET ER frame to the ROBOT frame as shown in Figure 3, explicitely indicating the presence of a frame shift.

Similarly we add another field called Axis that is in charge of indicating to which axis or plane of the Frame is related the data. With this approach, we are capable of detecting and preventing misconnections between Atoms during the design phase. To understand their importance let us consider the Atom presented in Figure 4 which is in charge of controlling the heading angle (Yaw) of a robot (either with a static reference or for path following).

Its Physics contains the following equation:

acc = acc des + K1 * de + K2 * e
where e is the Yaw error, de its derivative, acc des the desired acceleration profile and acc the resulting acceleration to be applied to the robot's dynamic model. where the Axis is set to z to indicate that the angle revolves around the z axis of the WORLD frame, NOFRAME and NOAXIS indicate that this Need is not related to a particular Frame and Axis and FrequencySquared is the quantity name we have associated with the unit Hz 2 .

To end the description of interface, we must introduce a final notion which is Interface Parameters. Indeed the previous equation used to control the Yaw angle of our robot could as well be used to control its Roll or its Pitch, the only difference being the Axis related to the interface. We could make three Atoms for those three situations but that would obviously be a waste of time (and also a quite error-prone process). Thus we introduce a parameter that will allow us to indicate that the Axis may vary depending on how the Atom is used. Hence in the previous interface, z can be replaced by an Interface Parameter.

We can of course also define Interface Parameters on the Frames. However defining Interface Parameters on Types is not allowed in our approach to avoid the design of models which would be described in a too abstract way.

D. Temporal constraints

We need to specify a given number of additional information to be capable of efficiently conveying knowledge from control description to its implementation.

These information, that we call Constraints, should allow us to encompass both the requirements of the control engineer and what the system can do, as described in [START_REF] Derler | Cyber-physical system design contracts[END_REF].

We have chosen to focus on temporal constraints since they are particularly relevant in our context. Indeed control engineers study stability of their control in continuous time through tools such as Lyapunov functions. However since control is executed in a discrete-time environment, periods of execution become detrimental to maintain the stability of the controlled system (though current tools to determine these periods are very limited). On the other hand, temporality is instrumental for the software engineers to know with which mechanisms of their architecture and middleware they should wrap knowledge (for instance consider the difference between synchronous and asynchrounous execution). Moreover it is crucial for them to know how to set up the architecture and especially at which period the synchronous components should run.

Thus temporal constraints are applied both on the Needs of the Atom and on its Physics. A Constraint is a pair (value, properties) with the set of available properties depending on the current Constraint value.

On the Physics, Constraints represent how the Physics should be run temporally and thus how Products values are updated. There exists three possible Constraint values:

• Constant: the Physics is run only once (initialization).

• Sporadic: we don't know when the Physics will be recomputed. It is associated with a property called C s phy . This property refers to an Atom describing the conditions under which the physics should run (an "important enough" change in a input value for instance). This woud allow us, for instance, to describe controllers based on Event-triggered control [START_REF] Heemels | An introduction to eventtriggered and self-triggered control[END_REF].

• Periodic: the Physics should run on a periodic basis. It is associated with two properties called t comp max representing the worst-case computation time of the Physics (this value of course depends on the implementation target) and T phy which is the allowed or possible set of intervals of periods for the Atom. For instance, in the case of a control Atom, T phy will correspond to the set of periods for which we can guarantee that the controlled system remains stable. On the other hand, for a sensor, T phy represents its set of update rates.

Constraints related to Atom Needs describe how the value of the Needs should be updated through time to ensure that computations performed in Physics always produce the expected result and/or preserve system stability. There are 4 kinds of constraints:

• Constant: the Need's value is set once and for all during initialization. • Sporadic: we don't know when the Need's value will be updated. • Coupled: the Need has to be updated every time the Physics is runned. This induces a temporal coupling between the two Atoms linked together. • Periodic: the Need is updated periodically (indepentently of Physics Constraints). It is associated with a property T Need specifying the allowed update periods. The Constraint on a Need n, noted Ctr(n) and the Constraint on a Product (and hence on the associated Physics) p, noted Ctr(p), are said compatible if they are both Constant, both Sporadic, both Periodic or Coupled with Periodic. This is noted Ctr(n) ≡ Ctr(p).

IV. A CONTROL LAW AS A COMPOSITION OF ATOMS

In the previous section, we introduced a minimal entity called Atom. To describe a control, we now need to connect them together.

A. The Molecules

Since Atoms are basic entities, it would be useful to allow to have an entity allowing us to use groups of Atoms without having to describe by hand their association every time we need them. This is the role of Molecules defined as:

B. Composition of atoms

The process of connecting Atoms and/or Molecules together through their interface is called Composition. It results in creating a Composition of Atoms which is the atomic description of a control law. However in order to use an Atom or a Molecule in a Composition, we have to set up how it will be used. This means we have to valuate its Interface Parameters and its Constraints (though the properties associated with the Constraints don't have to be since they may depend on the implementation target). This process will be referred to as Instantiation and the resulting elements as Instances.

A connection between a Product of an Instance and a Need of another one is called a Link. Let n be a Need and p a Product. The Link between them is noted L(n, p) and we have : f an-in(n) = 1 and f an-out(p) = ∞.

V. THE ATOMS THROUGH EXAMPLE

In order to allow explanations and figures to remain simple to explain and to display, we will focus on a quite simple application: a heading angle control algorithm.

Figure 5 presents the legend that will be used in next Figures. Moreover, the Links between Atoms have been merged for the sake of readability of the figures. Such a graph is presented in Figure 7 to Figure 9. The "Bridge" presented in Figure 1 then consists in manipulating this graph to structure the control and extract relevant information for its implementation as shown in Figure 6. Our robot, the Jack, is equipped with an Inertial Measurement Unit (IMU), a Loch Doppler and a Depth sensor. From these sensors, we can reconstruct our robot's state as shown in Figure 7. Moreover, since the way the IMU and Loch Doppler are mounted can influence the measurements, we must perform Frame Shifts in order to express the measurements in Robot's body frame. In Figure 7, ExpertiseIMU and ExpertiseLoch are External Knowledge Atoms used to indicate how they are mounted (information coming from the robot's manufacturer). You must also note three Molecules RobotSpeedsBody, RobotSpeedsWorld (an Instance of the same Molecule as RobotSpeedsBody but with a different Interface Parameter) and RobotPose. They are used as buffers allowing to insulate the computation of the robot's state from the rest of the Composition. This is a design technique ensuring a better modularity in control design since we can change the way they are computed (by introducing sensor fusion for instance or by using a SLAM algorithm to compute robot's pose) without having to reroute the Links referring to Robot's state thus simplifying the changes.

Before the Actuators, the Dynamic Model (Figure 8) is used to turn the desired accelerations into forces and torques expressed in the robot's frame. Desired accelerations coming from the control don't have to be updated at every iteration hence allowing us to introduce Periodic Links (Figure 9) while robot's speeds are required at every computation to determine the damping and coupling terms to compensate them thus inducing a temporal coupling with Robot's state reconstruction. The control (Figure 9) is based on the DynamicAngleControl Atom presented in section III which is used to control the heading angle. The related ExpertiseYawControl allows to set the values of K1 and K2 as well as set to 0 the desired acceleration and speed (since reference is static). The Yaw error is computed by a simple substractor and the desired angle is set by the operator through a Network link whose value is updated on a Sporadic basis. The error derivative is computed by an AngularSpeed substractor. Moreover since the produced acceleration is expressed in the WORLD frame, we have to frame shift it to ROBOT frame in order for it to be used by the Dynamic Model which expects its inputs in ROBOT frame (Figure 9).

Finally we have introduced Separator Atoms between Robot's states and the control which uses them (Figure 9). The Separator Atoms are used to introduce periodic Links thus enforcing time-decoupling between the control Atoms and the ones in charge of reconstructing robot's state. Hence both can run at optimum periods without being constrained by the other part of the Composition.

B. Constraints evaluation

We should now see how a Composition can be structured so that we can merge the Constraints exerted on each Instance into Constraints exerted on their Composition as shown in Figure 6. We begin by checking the control validity. Then in order to simplify the Graph, we merge edges that connect the same nodes with identical Constraints. Finally since our main goal is to determine the period of executions of all periodic Instances, we remove the nodes and edges corresponding to non periodic Instances. The resulting Graph for our example is shown in Figure 10. We then use the periodic Links to separate the Composition that are temporally independent and thus their Contraints will be studied independently.

From this we can determine the Constraints related to each Composition through the Diffusion process (Figure 6).

Let us consider the second one, C2. We have:

t comp max (C2) = 7 ∑ i=1 t comp max (Inst i ) + ∆ com (1)
where ∆ com represents the communication time between the Atoms of the Composition and Inst i the set of Instances included in the Composition.

We can then compute the set of possible execution periods of the composition, T exe , through intersection of allowed period sets of the different Instances:

T exe (C2) = 7 i=1 T phy (Inst i ) ∩ [t comp max (C2) ∞] (2)
We then valuate the different properties. If T exe of any Compositions is an empty set then the control cannot be implemented on the desired target. Otherwise one can choose a period T 2 ⊂ T exe (C2) as the period of execution of all Atoms in the Composition. Similarly we can choose T 1 ⊂ T exe (C1). This approach is for now an offline process but we are considering its embedding in the application to cope with potential changes in the temporality of some Atoms.

C. HIL setup presentation

Our HIL setup is constituted of the BeagleBone Black [START_REF] Coley | BeagleBone Black System Reference Manual -Revision C[END_REF] equipping our robot and which serves as main controller. A Linux kernel 3.8.13 with a Debian distribution and patched with Xenomai 2.6.3 is installed on the Bea-gleBone. It embeds the ContrACT real-time middleware developed at LIRMM [START_REF] Passama | ContrACT: a software environment for developing control architecture[END_REF] that allows to precisely configure software components' timing properties, which reveals very complementary to the Atom approach. The controller is connected to the simulation PC through a network link and both communicate through an UDP protocol ensuring non blocking communications. The simulator contains a model of the robot, environment representation and sensors simulation with noise.

D. Simulation Results

The following results have been obtained with T 1 = 25ms, T 2 = 50ms, K1 = 10.0 Hz, K2 = 10.0 Hz 2 and ψ des = 75 • . Figure 11 presents the heading angle evolution through time. Maximal noise levels on the measurements of ψ and r are noise ψ = 0.1rad and noise r = 0.05rad/s. first seconds of the simulation. There is no direct mapping between Atoms and real-time loops but the precise process is beyond the scope of this paper. C1 is implemented using 4 loops: Sensors and Actuators to replace the drivers of the real components, MDynamic where the Dynamic model is implemented and Navigation to reconstruct robot state. C2 is implemented by 2 loops: Control and Frameshift.

Figure 13 and Figure 15 present the overall execution times associated with the Compositions. They integrate computation and communication times between the real-time loops. We neglected context switching times and these measurements also do not integrate the applicative scheduler's computation times. These times remain inferior to execution periods and the variability observed can be explained by natural computation times variations, calls to non real-time functions with important duration variations (especially the ones used to log data) as well as possible preemptions by the scheduler especially since we use other loops to communicate with the simulator and to store data. They also illustrate how hard it is to have a precise estimation of t comp max values for Atoms and Compositions hence underlining the importance of HIL simulation in the design process.

Figure 12 and Figure 14 represent the evolution of delays of the real-time loops start times compared to the theoretical ones after the initial scheduling. For the execution i of a control loop, the delay is computed as:

delay i = t start the i -t start real i (3)
These delays are influenced by computation durations which impact start times of loops. Moreover since the BeagleBone possesses only a single-core processor, potential delays on a loop may also impact other loops start times. However the delays remain very limited and, if some appear, the scheduler manages deadlines to compensate for them. The overall delays remain within admissible bounds due to system imprecision as shown in Table I. The loop MDynamic is the most affected one by delays because it is the most computation time demanding loop and is very susceptible to be impacted by other loops (since it needs data from most of them). Thus it is the hardest one to manage for the scheduler explaining delays superiors to the other loops.

VI. CONCLUSION AND FUTURE WORK

We have presented a new control design paradigm based on composable entities named Atoms and which benefits from object-oriented principles in order to design modular control laws.

We have illustrated it through an example described with this paradigm and shown how a control can be represented as a Knowledge Association Graph which is then manipulated to determine the Constraints applied to our Compositions, Constraints which are then used by software engineer to parameterize the real-time software architecture.

Further work is currently being done on defining a new entity called Alternative used to integrate concepts from switching and hybrid control to our control description methodology. Assessing the impact of the decomposition of monolithic control in terms of performance is also an important issue. Finally we would like to push forward the notion of Knowledge Association Graph to try an see if any further control properties could be studied from it as well as which kinds of Constraints it would be relevant to add to the already existing temporal constraints.

Fig. 1 .

 1 Fig. 1. The "Bridge" between control design and software implementation III. ENCAPSULATING KNOWLEDGE: TOWARDS THE ATOMS

Fig. 2 .

 2 Fig. 2. The simple example of a control Atom using sonar measurements Let us consider the Figure 2 example where measurements performed by a sonar are used by a control Atom. How can we interpret the control designer's intentions ? With no other information the following assumptions are possible:

Fig. 3 .

 3 Fig. 3. With the Atom Datatype, we are forced to express explicitely the shift from PROXIMETER to ROBOT frame.

Fig. 4 .

 4 Fig. 4. Interface of the Atom containing the Yaw control equation It is clear that e, de, acc des and acc should relate to the same Frame and Axis, if not the control's behavior will not be the expected one. The Atom's Interface can then be described as : Ne = {(e, [Angle,WORLD, z]), (de, [AngularSpeed,WORLD, z]), (K1, [Frequency, NOFRAME, NOAXIS]), (K2, [FrequencySquared, NOFRAME, NOAXIS]), (acc des, [AngularAcceleration,WORLD, z])} Pr = {(acc, [AngularAcceleration,WORLD, z])}

Definition 2 :

 2 A Molecule is a composable entity whose Physics is made of a composition of Atoms or Molecules. A Molecule M is a 5-tuple: M = (Name, Phy, Int, IntPar,Constraints) with Phy = (E, l, L) where E is the set of entities (Atoms or Molecules) used in the Physics, l the set connections of these entities with the interface of the Molecule and L the set of connections between the entities.

Finally

  a valid Composition is defined as: Definition 3: A Composition made of a set I of Instances and a set L of Links, noted C = (I, L), is valid if and only if L i is valid for all i. Definition 4: A Link L(n, p) is valid if and only if D(n) = D(p) and Ctr(n) ≡ Ctr(p).

Fig. 5 .

 5 Fig. 5. Legend used on Figures 7 to 9

A

  . Control Description Thanks to the fine grain description allowed by atoms we can represent a control description as a directed multigraph known as Knowledge Association Graph. Definition 1: A Knowledge Association Graph G C = (V, E) is the directed multigraph associated with the Composition C = (I, L) so that: V = I and E = L

Fig. 6 .

 6 Fig. 6. The Different steps linking control design and implementation

Fig. 7 .

 7 Fig. 7. Our example described as a Knowledge Association Graph part 1

Fig. 8 .Fig. 9 .

 89 Fig. 8. Our example described as a Knowledge Association Graph part 2

Fig. 10 .

 10 Fig. 10. The Graph with only periodic Atoms and separated Compositions

Fig. 11 .

 11 Fig. 11. Heading Angle vs setpoint with sensor noise Figures 12 to 15 contain the execution time measurements performed on the real-time loops in which Compositions C1 and C2 have been implemented. We focused on the 20

Fig. 12 .

 12 Fig. 12. C1 delays to theoretical start times

Fig. 13 .

 13 Fig. 13. C1 execution durations

  REFERENCES[1] P. Ulbrich, F. Franzmann, C. Harkort, M. Hoffmann, T. Klaus, A. Rebhan, and W. Schröder-Preikschat, "Taking control: Modular

Fig. 15 .

 15 Fig. 15. C2 execution durations

TABLE I AVERAGE

 I DELAYS OF LOOPS AND DESRIED PERIODS.

	RT Loop	Comp.	Av. delay(ms) Period(ms) Error(%)
	Sensors	C1	-0.3583	25	1.433
	Actuators	C1	-0.1535	25	0.614
	Navigation	C1	-0.2808	25	1.123
	MDynamic	C1	-1.0358	25	4.143
	Control	C2	-0.4878	50	0.976
	Frameshift	C2	-0.0173	50	0.035


The quest of vertical flight: the development of helicopter Thierry M. FAURE Maître de Conférences (HDR), Centre de Recherche de l'Armée de l'air (CReA), École de l'Air, 13 661, Salon Air, France Vendredi 4 septembre 2015

The idea of vertical flight

The dream of flight, first developed by the observation of nature, was inspired from birds and insects.

However, the whirling fall and glide of seeds from trees like the sycamore, allowing a wide dissemination of the specie, is an alternative way of flight. That is probably the origin of conception of the Chinese top (reported by Ge Hong circa 400 BC), which is a small toy consisting in two blades attached to a stick, which rotation between the hands of an operator, generates lift and its vertical flight [START_REF] Leishman | Principles of Helicopter Aerodynamics[END_REF]. The idea of using blade rotation, to fly a machine carrying a man, was developed by Leonardo da Vinci at the end of the XV th Century. His thinking was to adapt the Archimedes' screw, which is a device for raising water from the ground level, to the vertical translation of a machine in the air (Figure 1). If no light and highstrength material existed at that time and no engine was available to provide enough power to fly that machine, its design was wrong in its mechanical concept. Contrary to the seeds falling from a tree, or the Chinese top, the shaft-driven rotor generates a torque which is transmitted to the nacelle making it swirl and leads to a loss of control. Thus, the lack of a proper counter-reacting torque system is the main drawback of that design. In addition, the rotor shape is far from an optimum use in the air. In Russia, Mikhail Lomonosov realized in 1754 a model inspired by the Chinese top to raise a charge from the ground [START_REF] Leishman | Principles of Helicopter Aerodynamics[END_REF]. That model consisted in two rotors powered by a wound-up spring device; the charge was inside a box in translation between two vertical poles. That model did not require counterreacting torque system since it was linked to the floor. In France in 1784, Launoy and Bienvenu developped an evolution of the Chinese top including a power source, a string wounded around the rotor shaft and tensioned by a crossbow [START_REF] Ross | Flying Windmills[END_REF]. For the first time, torque reaction was cancelled using two counter-rotating rotors (Figure 2). That device was successfully demonstrated before the French Academy of Sciences and renewed the interest of rotor vertical flight in the scientific community.

First designs and models

During the XIX th Century, the cancellation of torque reaction was well understood and solved by counter-rotating rotors. Many vertical motionned rotor-driven machines were designed and sometimes tested. They were generally moved by spring or steam power, but the lack of engines developing sufficient power-to-weight ratios did not lead these machines to success. Following the flights of a hotair balloon, by Pilatre de Rozier and the Marquis d'Arlandes, and of a hydrogen balloon, by Charles and Robert, in 1783, a key issue was to control the trajectory. That period was then marked by the struggle between the "lighter than air" and "heavier than air" flight supporters, the latter being dominated by the defendant of vertical rotor propulsion. It is reported [START_REF] Ross | Flying Windmills[END_REF] that in England William Henry Phillips, the inventor of the extinguisher, developed and flew in 1842 a miniature model of vertical rotor-driven machine powered by steam ejected at the tip of the blades. Because that rotor was not driven by shaft rotation, a torque reaction system was useless. After making a four-bladed version of Launoy and Bienvenu's model in 1796 [START_REF] Ackroyd | Sir George Cayley: The Invention of the Aeroplane near Scarborough at the Time of Trafalgar[END_REF], Sir George Cayley developed the concept of Convertiplane in 1843 [START_REF] Cayley | Retrospect of the progress of aerial navigation, and demonstration of the principles by which it must be governed[END_REF] (Figure 3). That kind of girodyne was intended to take-off vertically with two pairs of counter-rotating rotors, then to fly horizontally by the conversion of the rotors into lifting circular wings, the horizontal propulsion being addressed by two vertical rear propellers. That machine was never built. In France, Gustave de Ponton d'Amécourt (1861) made a steamed engine model called for the first time "helicopter" [START_REF] De Ponton D'amécourt | Collections de mémoires sur la locomotion aérienne sans ballons[END_REF]. To save weight, that coaxial counter-rotating rotors machine was the first to use aluminum, which was a new light material at that time (Figure 4). A rise of different helicopter designs and models followed, and many thought that heavier-than-air flight would be realized very soon. The literature with Jules Verne's novel "The Clipper of the Clouds" (« Robur le Conquérant ») published in 1886 [START_REF] Verne | Robur le Conquérant[END_REF] inspired that standpoint. In 1878, the Italian Enrico Forlanini flew during 20 s a miniature model, using steam power, which reached the height of 13 m [START_REF] Leishman | Principles of Helicopter Aerodynamics[END_REF]. Up to 1880, all the helicopter models employed short blades. Thomas Alva Edison developed small helicopter models and was the first to found the need for a large rotor diameter with low blade area to give good hovering efficiency. Unfortunately, the engines tested did not provide a sufficient power-to-weight ratio.

Theoretical developments

During this period, in the field of theoretical approach, the Scottish William Rankine in 1865 [START_REF] Rankine | On the Mechanical Principles of the Action of Propellers[END_REF] and the British William Froude in 1878 [START_REF] Froude | On the Elementary Relation Between Pitch, Slip and Propulsive Efficiency[END_REF] and his son Robert Edmund Froude in 1889 [START_REF] Froude | On the Part Played in Propulsion by Differences of Fluid Pressure[END_REF], developed the propeller theory. First applied to naval propulsion, this theory was then adapted to vertical flight dynamics from the beginning of the XX th Century. Many improvements followed, with Drewiecky's concept of a hybrid momentum-blade element (1900) [START_REF] Drzewiecky | Des hélices aériennes : théorie générale des propulseurs hélicoïdaux et méthode de calcul de ces propulseurs pour l'air[END_REF] and other approaches including tip vortices dynamics by Joukowski (1906) [START_REF] Joukowski | On annexed bounded vortices[END_REF]. A practical helicopter would not limit its operating range to vertical translation, but has to combine vertical and horizontal motions. The forward flight of a helicopter is an important issue which produces a non-uniform velocity distribution on the rotor and creates a rolling moment. The way to suppress that rolling moment is the principle of flapping propeller blades, suggested in 1904 by Colonel Charles Renard. An operational system of flapping blades was patented in 1908 by Louis Charles Bréguet [START_REF] Leishman | Principles of Helicopter Aerodynamics[END_REF]. He built a whirling arm test rig in 1905 and found the way of producing a forward thrust by tilting the rotor blades. The effects of forward flight speed on the rotor flow were also analyzed by Joukowski in 1909 [START_REF] Joukowski | Théorie tourbillonnaire de l'hélice propulsive[END_REF]. He proved that the asymmetry of velocity and pressure on the main rotor generates asymmetry in forces and moments on the propeller in forward flight. But, unless Bréguet, Joukowski did not propose a technical solution to address this issue. The concept of cyclic pitch was one way for attaining full horizontal control and was patented in 1906 by the Italian Gaetano Arturo Crocco. He recognized that, in order to operate properly in forward flight, a helicopter needed a change of the blades pitch, to take into account the dissymmetry in the aerodynamic loads between the advancing blade into the relative wind, and the blade retreating away from the wind. Other torque-reaction systems were also suggested, as the synchropter, which consisted in two counterrotating side-by-side intermeshing rotors, and was patented by Max Bourcart in 1902 [START_REF] Bourcart | Machine à voler[END_REF].

First attempts of piloted vertical flights

The better performance of the internal combustion engine obtained at that time is linked to the development of airplanes, but that gain was not large enough to see a successful helicopter flight.

However, that new type of engine and the increase of their power-to-weight ratio provided to helicopter pioneers, a way of building machines able to carry a man. France was the main workshop of the research on helicopter. It is said that Paul Cornu, a bicycle maker, was the first one to make a vertical flight of 30 cm from the ground during 20 s on 13 November 1907 (Figure 5). However, the engine power was limited to 24 hp, and would not be enough to sustain hovering flight out of ground effect. The frame of the helicopter consisted in tandem rotors to cancel torque reaction. Note the short rotor blade shape, which was not the best aerodynamic design. Because of its absence of stability system, it was secured by ropes tighten by operators on the ground. This is probably the reason why it is not mentioned as the first helicopter flight. Again, the machine took advantage of ground effect. Nevertheless, the Bréguet brothers had a good knowledge of aerodynamics as they used rotors with long blades and low area.

Stability control and horizontal flight

The theoretical improvements on rotor wing dynamics were rapidly followed by the development of innovative technologies. In Danemark, Jen C. Ellehammer built a co-axial rotors helicopter in 1914 which made many hops into the air. Each rotor was fitted with six short blades and a lower disk covered with fabric, to serve as a parachute in case of an engine failure. The main innovation of that helicopter is the first use of a rotor cyclic pitch mechanism. In the United States, a Russian émigré expert in rotor-wing theories [START_REF] De Bothezat | The General Theory of Blade Screws[END_REF], Georges de Bothezat, built for the U.S. Army a quad-rotor helicopter (Figure 9). It was fitted with collective and cyclic blade pitch variations and four smaller rotors to help control. That machine flew successfully at a low altitude and a low forward speed in 1922, but the performances were not sufficient and the Army stopped its funding in 1924. Raúl Pateras Pescara, from Argentina, worked in France and Spain and focused on the design of the rotor, the key element of the helicopter, which he transformed into a complex mechanical system. He developed many co-axial rotors helicopter. Each rotor was made with five biplane blades with a cyclic pitch and control mechanism; yaw was controlled by differential collective pitch. The biplane blade design of the rotor needed bracing wire and generated a large drag, and the first models were underpowered [START_REF] Leishman | Principles of Helicopter Aerodynamics[END_REF]. A new version with a more powerful engine was the Pescara 2F (Figure 10 Again in France, Étienne OEhmichen, funded by the Peugeot car company, developed a quad-rotor helicopter with eight additional rotors for horizontal direction control (Figure 11). That machine was the first one with flight capabilities not limited to vertical and straight line motion, but demonstrated some maneuverability [START_REF] Oemichen | My Experiments with Helicopters[END_REF]. In May 1924 OEhmichen was awarded a prize of the FAI for a closed-loop circuit of one kilometer in 7 min 40 s.

The autogiro

An important milestone in the elaboration of an operational helicopter is the autogiro. That gyroplane consists in a traditional airplane with a propeller and fixed wing, with the addition of a non-powered horizontal rotor above the fuselage, freely rotating on its shaft by the action of the forward speed. It was unable of vertical take-off and landing, but take-off and landing distances were largely shorten in comparison with a traditional airplane. That b) design was developed by the Spanish Juan de la Cierva in order to provide additional lift by the autorotation of the horizontal rotor, to prevent an eventual engine failure. The first model C-4 flew in 1923 (Figure 12-a). De la Cierva elaborated and realized new technical solutions as the first practical application of flapping hinges as a means of equalizing the lift on the parts of the rotor in forward flight [START_REF] De La Cierva | The Development of the Autogiro[END_REF]. De la Cierva was invited to Great Britain in 1925 by the Weir Company and a first autogiro model, the C-19, was produced from 1932 to 1936 (Figure 12b). In later models, de la Cierva added a lag hinge to each blade to cancel stresses caused by centrifugal forces, completing the articulated rotor hub design. The autogiro was very popular in the 1930s and was produced by many companies. The Austrian-born British engineer Raul Hafner developed the cyclic pitch control on its autogiro AR-3 in 1935, which provided a means of increasing collective pitch and also tilting the rotor disk without tilting the rotor shaft with a control stick as in de la Cierva's direct control system [START_REF] Leishman | Principles of Helicopter Aerodynamics[END_REF]. Ironically, although motivated by air safety, Juan de la Cierva died in an airplane crash in 1936. The autogiro development was then stopped and the helicopter was preferred.

The first reliable helicopters

With the maturation of the key technologies in rotor control, like the cyclic and collective pitch control, the flapping and lagging hinges and anti-torque reaction rotors, the emergence of a useful helicopter was possible in the years 1930s. Many rotor configurations emerged at that time (Figure 13), and they are still in use by present helicopter manufacturers. In Belgium, the Russian émigré Nicolas Florine was the first to fly a tandem rotor helicopter (Figure 13a). That rotor configuration is more compact than the quad-rotor machines previously tested by Bréguet, de Bothezat and OEhmichen. In 1933, it reached the altitude of 5 m during a flight of 9 min (Figure 14). Anton Flettner built an autogiro in 1935, and then developed the synchropter Fl-265 in 1939 (Figure 17).

This type of helicopter possesses two counter-rotating side-by-side intermeshing rotors (Figure 13-d).

The Fl-265 was much easier to control than the Fa-223, and was the first helicopter to demonstrate transition into autorotation and then back again into powered flight. An improved version, the Fl-282 Kolibri (Figure 18), was produced.

Conclusion

The development of the helicopter is a step-by-step story resulting from original design ideas, fundamental aerodynamics theories and innovative technical solutions. That course benefitted the internal combustion engine improvements, because of the large power-to-weight ratios needed, and the a) b)

emergence of light and high-strength materials, but the key of a practical helicopter relied in its control, addressed by the elaboration of the complex mechanical system of the rotor head. That long time development of knowledge and technology, in the first decades of the XX th Century, led many pioneers interested in helicopter development, like Bréguet and Sikorsky, to bring a contribution, to turn their mind to other aeronautical fields before achieving a successful helicopter. It is worth noticing that the funding provided by their own aircraft company, together with their skills in aeronautics, were important elements of success.

Figure 1 :

 1 Figure 1 : Aerial screw.

Figure 2 :

 2 Figure 2 : Launoy and Bienvenu's invention (1784).

Figure 3 :

 3 Figure 3 : George Cayley's Convertiplane (1843) [4].

Figure 4 :

 4 Figure 4 : Helicopter of de Ponton d'Amécourt (1861).

Figure 5 :

 5 Figure 5 : Helicopter of Cornu (1907). Previously the same year, brothers Louis Charles and Jacques Bréguet tested the Gyroplane No. 1, which was under development since 1905 and based on the work of Professor Charles Richet. The machine consisted in four double four-bladed rotors arranged at the end of each branch of a cross, this geometry is called quad-rotor, and was driven by a 40 hp engine (Figure 6). It is reported that the Gyroplane No. 1 lifted off and reached 60 cm during 1 min on 14 August and 29 September 1907.

Figure 6 :

 6 Figure 6 : Gyroplane No. 1 (1907).In Russia, Igor Ivanovitch Sikorsky (1910) built vertical flight models. After the failure of the S-1 helicopter, the S-2 was fitted with a more powerful engine. It consisted in two three-bladed counterrotating rotors and made only small hops from the ground. Discouraged, Sikorsky turned his mind to fixed-wing airplanes.

Figure 7 :

 7 Figure 7 : Sikorsky S-2 (1910).

Figure 8 :

 8 Figure 8 : Ellehammer's helicopter (1914).

Figure 9 :

 9 Figure 9 : de Bothezat's helicopter (1924).

  ) which established on 18 April 1924 a world record of flight distance in straight line of 736 m with an altitude of 1.8 m in 4 min 11 s recorded by the Fédération Aéronautique Internationale (FAI), but the helicopter control was limited.

Figure 10 :

 10 Figure 10 : Pescara 2R (1921) and 2F (1923).

Figure 11 :

 11 Figure 11 : OEhmichen's helicopter No. 2 (1924).

Figure 12 :

 12 Figure 12 : de la Cierva's autogiros: a) C-4 (1923), b) C-19 (1932).

Figure 13 :

 13 Figure 13 : Rotor configurations: a) tandem rotors, b) co-axial rotors, c) lateral side-by-side rotors, d) synchropter, e) main rotor with tail anti-torque rotor.

Figure 14 :

 14 Figure 14 : Florine's tandem rotor helicopter (1930). After establishing a noted aircraft company, Louis Charles Bréguet came back to the helicopter development with René Dorand. They improved the coaxial counter-rotating rotors concept (Figure 13b) with a cyclic pitch control. The yaw angle control was obtained by applying a differential torque on a rotor. In 1936 the helicopter established a 62 min flight of 44 km (Figure 15).

Figure 15 :

 15 Figure 15 : Bréguet-Dorand (1936).

Figure 16 :

 16 Figure 16 : Focke's helicopters: a) Fa-61 (1936), b) Fa-223 (1940).

Figure 17 :Figure 18 :

 1718 Figure 17 : Flettner Fl-265 (1939).

Figure 19 :

 19 Figure 19 : Igor Sikorsky testing his helicopters a) VS-300 (1939), b) VS-300A (1940).

Figure 20 :

 20 Figure 20 : Sikorsky R-4 (1941).


L Marioni

J Alves

F Bay

E Hachem

EFFECT OF M-EMS ON IN-MOULD TRANSIENT FLOW DURING CONTINUOUS CASTING

Electromagnetic braking (EMB) and electromagnetic stirring (EMS) are widely used to increase the efficiency in continuous casting process of steel. In particular, in-mould applications (M-EMS) allow decreasing the free surface fluctuation, controlling the velocity field and decreasing the turbulence of the flow. Since laboratory-scale tests are not fully representative of the process and industrial measurements are both expensive and difficult to carry out, numerical simulation is a strong tool to study and optimize these electromagnetic applications in steel industry. These simulations cannot fully model the process because the problem results to be strongly multiphysical and the simulation of all the phenomena involved would lead to huge computational costs, which is one of the main limits in the current situation. For this reason, this work aims at starting searching a coupling algorithm which could guarantee both computational efficiency and accuracy of the final results.

INTRODUCTION

In continuous casting processes the liquid metal is injected into the mould. The final steel shell is obtained after the solidification which starts in the mould and continues in the strands. Electromagnetic devices such as stirrers and brakers are well known technologies used to improve both the quality of the final product and the casting speed. The main defects of the final shell in terms of micro-structure and surface cracking can be directly related to in-mould phenomena: temperature variation, velocity and pressure of liquid steel, free-surface behaviour and slag entrainment are some of the main causes of defects in the final product. With the purpose of controlling the process and preventing the final product from defects, the process has been enhanced with electromagnetic devices such as stirrers (EMS) and brakers (EMB). The main difference is that the stirrers work under the supply of AC current and produce dynamics magnetic fields. Brakers are permanent magnets or circuits fed by DC current, thus they produce constant magnetic fields. Despite the differences they are based on the same idea that the superposition of a magnetic field to the metal flow will generate Lorentz forces which can drive the flow in accordance with the process' design. The physical phenomena occurring in the mould is a multiphysics problem which includes liquid flow, multiphase analyses, electromagnetic computation, heat transfer and solidification processes where each of these physics depends on the others. Numerical simulations are one of the main tools used to design the process because it is too complex to be solved analytically. Also, experimental data is very expensive and difficult to obtain because of the high temperature of the casting process. The main drawback of these simulations is the high computational effort required to simulate the highly interconnected phenomena such as electromagnetic field, multiphase fluid flow, heat transfer and solidification. The study of an efficient coupling is therefore essential to neglect non-relevant aspects of the modelling and obtain fast and reliable simulation algorithms which can be exploited in the industrial world.

SIMULATION TEST CASE

The main purpose of the current work is to create an efficient coupling strategy able to simulate electromagnetic influence in continuous casting processes. For this reason, a simplified test case has been used. The commercial mould considered in the work has been taken from Singh's work [START_REF] Singh | Large eddy simulations of double-ruler electromagnetic field effect on transient flow during continuous casting[END_REF], but only a quarter of it has been included in the simulation (see figure 1). A simplified electromagnetic stirrer has been positioned away from the mould in the direction of the narrow face and between and under the free surface level. The EMS device is a solenoid fed by a pulsating AC current and produces a horizontal pulsating magnetic field directed from the narrow face in the direction of the nozzle plan. The geometry of the device is simplified since available data on commercial stirrers were not accurate enough to perform a simulation. In the current work we content to simulate a realistic Lorentz force field, since the aim is to study the interaction of the stirrer on the flow more than to obtain specific industrial results. The main data used in the simulation are reported in table I.

MATHEMATICAL MODEL Governing equations

The potential model has been used to model the electromagnetic field (EMF). Thus, the evolution in time of the EMF is describes by the following system:

{ ( ) (1)
where and are the free space and material's relative magnetic permeability respectively. is the magnetic vector potential, is the scalar electric potential, is the medium velocity, and is the electrical conductivity. The convective term can be neglected from the first equation in system (1) because the Reynolds magnetic number is . The same term can be neglected also in the second equation of system (1) because , which leads to . For further justification to these assumptions, we refer to [START_REF] Mureau | Magnetohydrodynamics[END_REF]. These conditions ensure that the EMF field is not affected by the fluid motion in EMS applications, while in EMB the convective term could not be neglected in the second equation of system [START_REF] Singh | Large eddy simulations of double-ruler electromagnetic field effect on transient flow during continuous casting[END_REF]. The electromagnetic volumetric force is computed as: [START_REF] Mureau | Magnetohydrodynamics[END_REF] where the convective term due to fluid motion has been neglected for the aforementioned reasons. Both the time evolving and time averaged Lorentz forces have been considered in two different simulations because the literature does not provide a unique criterion to accept the usage of the averaged field. Such criterion should be dependent on frequency [START_REF] Felten | Numerical modelling of electromagnetically-driven turbulent flows using LES methods[END_REF], on the skin depth phenomenon and on the Stuart number [START_REF] Barna | Numerical Simulations of Mould Electromagnetic Stirring for Round Bloom Strands[END_REF]. The electromagnetic model is solved by using Nedelec edge finite elements [START_REF] Nedelec | Mixed finite elements in[END_REF]. The turbulent fluid flow has been modelled by the incompressible Navier-Stokes equation. The Lorentz forces have been considered as a source term. No buoyancy effect was considered, since it is negligible in comparison to inertia and electromagnetic influences. A variational multiscale (VMS) approach [START_REF] Hughes | The variational multiscale method -a paradigm for computational mechanics[END_REF] has been used to model both the fine and the coarse level of turbulence. Both velocity and pressure fields have been enriched and stabilized according to [START_REF] Hachem | Stabilized finite element method for incompressible flows with high Reynolds number[END_REF]:

{ ∑ ∑ ∑ (3)
with and the momentum and continuity residuals from the modelling of the fine scale, respectively. System (3) considers the fine scale modelling without explicitly tracking it. The stabilizing terms and (whose theoretical derivation is detailed in [START_REF] Hughes | The variational multiscale method -a paradigm for computational mechanics[END_REF]) enable the model to overcome the instability of the classical formulation arising in convection dominated flows and to deal with the spurious pressure oscillations. The coupled system has been simulated through the so called two meshes two solvers approach [START_REF] Barna | Numerical Simulations of Mould Electromagnetic Stirring for Round Bloom Strands[END_REF][START_REF] Fujisaki | In-Mold Electromagnetic Stirring in Continuous Casting[END_REF][START_REF] Thomas | State of the Art in Electromagnetic Flow Control in Continuous Casting of Steel Slabs: Modeling and Plant Validation in 6th International Conference on Electromagnetic Processing of Materials EPM[END_REF]. The electromagnetic problem has been solved with Matelec® over a domain which included the air and the simplified stirrer. The Lorentz forces have thus been imported in Thercast ® which solves the mechanical problem over the fluid domain. This approach allows us to exploit specific solvers and optimized meshes for each involved physic.

Anisotropic mesh adaptation

In [START_REF] Coupez | Metric construction by length distribution tensor and edge based error for anisotropic adaptive meshing[END_REF] a posteriori error estimate based on the length distribution tensor approach and the associated edge based error analysis is proposed. This error estimation algorithm is the base for the anisotropic meshing scheme adopted in the current work developed in [START_REF] Coupez | Solution of high-Reynolds incompressible flow with stabilized finite element and adaptive anisotropic meshing[END_REF] in the case of Navier-Stokes, multi-field re-meshing. The mesh in the fluid domain is computed as output of an error optimization analysis performed along the edges and referred to user defined fields; we propose the mesh to be adapted according to the following multi-field vector:

( ) { (| |) (| |) } (4)
where is the velocity field, is the distance from the boundary, is the volumetric force and are the nodal indexes. Two different strategies have been used to re-mesh:

 Implicit-Static adaptation: the mesh is anistropically adapted to the boundaries and in the weak regions. The obtained mesh is then used as a fixed mesh during the whole computation.  Explicit-Dynamic adaptation: the mesh is dynamically adapted during the computation according to the multicriteria vector (eq. ( 4)) computed from the solution. In figure 2 we can see how the mesh is anisotropically adapted to the flow in the nozzle area.

RESULTS

Flow simulation

The flow developed in the mould is the classical double-roll flow. In figure 3 the magnitude and the streamlines on the average velocity are plotted. The average velocity has been computed in the time interval when a nominally steady state is reached. In figure 3 the abnormal width of the injection flow can be noticed, which is an already observed drawback in symmetric simulations [START_REF] Yuan | Study of computational issues in simulation of transient flow ini continuous casting[END_REF]. This average flow is predictive of the macrophenomena occurring in the mould, but lacks information regarding localized turbulence, which has a strong impact on metallurgical properties. The turbulence leads to short-time subflows detected by the non-average solution; this difference can be measured in term of standard deviation, defined as:

√ ∑ ̅ (5)
where is the instantaneous velocity, ̅ is the average velocity, is the number of time-steps and represents the i-th time step. In figure 4 the normalized standard deviation of velocity during the nominal steady state regime is plotted. We see that a large deviation is present in the main flow region and the maximum deviation occurs where the injected flow impacts on the narrow face; the flow in this region is highly turbulent and this turbulence has effect on the solidification of the first layer of solid metal. When the stirrer is turned on (figure 4 right) the standard deviation is decreased in the meniscus zone, which implies the flow to be closer to the steady state. Despite the low level of the magnetic field in the inner zone of the mould, EMS affects the injection flow as well. The flow is smoother and easily turns in the double-roll pattern, while in the non-stirred simulation the unstable flow close to the narrow face led to dynamic effects in the main injection flow. A second effect of the stirrer is to decrease the global velocity of the flow, which is important especially from a simulation point of view. In figure 5 left, the velocity profile along a vertical line 5 cm far from the central axis is plotted. We see that the horizontal velocity related to the main eddy ( is high, especially when compared to the vertical velocity. This means that the simulated flow and its symmetric one (in the second half of the mould) would have a strong impact in this area, which would create turbulence and lead the global flow to be non-symmetric. For this reason, the flow with no EMS should be simulated without any symmetry condition. On the contrary, from figure 5 right we see that the horizontal velocity in the flow with EMS is lower as well as its ratio with the vertical velocity. This result implies that in the stirred case, the geometrical symmetry plane is more likely to behave like a mechanical symmetry plane. The main effects of the stirrer can be detected close to the narrow face, along the symmetry plane and on the injection flow, but not in nozzle area. In figure 5 right the velocity across the exit central line of the nozzle is plotted and no sensible differences between the stirred and non-stirred solution are predicted. The computational effort of the simulation is not highly affected by the weak coupling between the EM and CFD simulations: in figure 6 an analysis of the computing time is presented. The dynamic mesh adaptation is computationally expensive even if the re-meshing has been performed every 14 time steps. As we can see in figure 7, this choice is not good in terms of optimization. The anisotropic adaptation makes the elements to be stretched in the same direction as velocity: this means that the elementary stiffness contribution is the same for every direction and the final system is well conditioned. In-mould flow is far from being steady, so the velocity field can significantly change at every time step, so the mesh, which was optimized for the solution at a certain time-step, does not fit the new solution field. For this reason the resultant system is bad conditioned, except from the time-step when remashing is done, and a large number of iteration is needed to reach the convergence. On the contrary, the statically anisotropic mesh provides a mesh which is generally good for all the time steps and leads to a better conditioned system: the number of required iterations is decreased by 25% with respect to the isotropic mesh (with the same number of elements).

CONCLUSIONS

In this work a non-industrial EMS application was simulated. The flow has been simulated by using a VMS algorithm in order to take into account the small scale turbulence in an implicit way. The adopted anisotropic mesh adaptation algorithm allows the user to use a low number of elements; it also improves the conditioning of the final system, so less iterations are required to solve the problem. A weak coupling scheme is adopted to simulate the coupled electromagnetic-fluid mechanical phenomena. This scheme is computational light, but is able to predict the main features of the simulated application. Finally, a mechanic symmetry condition is considered to be accurate in the stirred flow simulation.

Figure 1 :

 1 Figure 1: Mould geometry (one quarter of the commercial mould)

Figure 2 :

 2 Figure 2: Isotropic starting mesh (left) and anisotropically adapted mesh (right)

Figure 3 :

 3 Figure 3: Averaged velocity magnitude (left) and streamlines in the double-roll area (right)

Figure 4 :Figure 5 :

 45 Figure 4: Velocity standard deviation with no EMS (left) and with EMS (right)

Figure 6 :

 6 Figure 6: Computational time analysis for dynamic re-meshed simulation

Figure 7 :

 7 Figure 7: Average number of iterations to solve the linear system

ACKNOWLEDGEMENTS

The authors give thanks to Transvalor S.A. for funding this research and for the continuous support which allowed the complete exploitation of its software.


L Marioni
email: luca.marioni@mines-paristech.fr

F Bay

E Hachem

Lid-driven cavity highly turbulent flow subjected to high magnetic field Determination of critical time-step for explicit MHD schemes

Acknowledgements:

The flow is highly affected by the magnetic field and its direction.

In some configurations, the CFL condition is not enough to guarantee the numerical convergence.

A new time-step threshold which guarantees convergence in explicit schemes is proposed.

An implicit modelling of the Lorentz force has to be included in the multiscale stabilized finite-elements solver.

The braking effect will be taken into account in the complete simulation of electromagnetic stirring.

Conclusions & Perspectives:

e Determination of the critical time-step:

 Solution:

• Computation of the limit time-step as the time-step which leads to a complete dissipation of the kinetic energy within one time increment when the magnetic field is normal to the velocity.

e Numerical instabilities:

• The Lorentz force is tracked explicitly.

• The excessive size of the timestep leads to completely brake the flow normally to B and to accelerate it in the opposite direction, i.e. . • The CFL condition is not capable anymore to guarantee numerical convergence.

e

Computational model:

Electromagnetic Field:

Lorentz Forces:

Navier-Stokes:

(multiscale stabilized FE [2] ) Anisotropic mesh adaptation [3] .

e  Problem:

 Solution:

 Numerical procedure:

• Highly turbulent flow A representative study case:

• Impose a constant magnetic field to control and brake the flow

• Lid-driven cavity benchmark [1] • Relate the Stuart number (N) to the target braking level



  Electromagnetic braking and stirring in continuous casting process of steel  Coupled simulation in commercial software: Control of fluid flow in CC:

gratefully thank Transvalor for funding this research and for the continuous support which allowed the full exploitation of its software's capabilities.


Vincent Savaux

Faouzi Bader
email: faouzi.bader@supelec.fr

Jacques Palicot
email: jacques.palicot@centralesupelec.fr

CMA-Based Blind Equalization and Phase Recovery in OFDM/OQAM Systems

Keywords: Blind equalization, OFDM/OQAM, CMA, Impulse Noise Cancellation

This paper deals with a constant modulus algorithm-(CMA) based blind equalization method for orthogonal frequency division multiplexing/offset quadrature amplitude modulation (OFDM/OQAM) scheme. The authors propose to combine the CMA with a low-pass filter in order to increase the convergence speed of the equalizer, and to avoid the detection of the received symbols with phase mismatch. Simulations results show that the convergence speed is largely increased, and the proposed method achieves bit error rate (BER) performance very close to that of the perfect channel estimation and equalization.

I. INTRODUCTION

The blind equalization is an interesting solution for the channel inversion (or deconvolution) since it guarantees the signal a maximum spectral efficiency. Indeed, these techniques do not use pilots, the blind equalizer only takes advantage of a few signal features to perform the channel inversion. In [START_REF] Farhang-Boroujeny | Multicarrier Modulation With Blind Detection Capability Using Cosine Modulated Filter Banks[END_REF], B. Farhang-Boroujeny has adapted the blind equalization using the constant modulus algorithm (CMA) proposed by D. N. Godard in [START_REF] Godard | Self-Recovering Equalization and Carrier Tracking in Two-Dimensional Data Communication Systems[END_REF] to the cosine modulated filter bank (CMFB) modulation. In [START_REF] Bölcskei | A Subspace-Based Approach to Blind Channel Identification in Pulse Shaping OFDM/OQAM Systems[END_REF] the authors proposed to use the second order moment of the received signal, while the CMA is used in [START_REF] Farhang-Boroujeny | Multicarrier Modulation With Blind Detection Capability Using Cosine Modulated Filter Banks[END_REF].

In filter bank multicarrier (FBMC), and in particular in orthgonal frequency division multiplexing/offset quadrature amplitude modulation (OFDM/OQAM), only few papers deals with blind equalization. Due the convergence delay of the blind methods and the uncertainty on the phase of the recovered symbols, the blind channel equalization remains a pending challenge not much studied by the scientific community, in particular in FBMC.

In this paper we propose a constant modulus algorithmbased (CMA) blind equalizer for OFDM/OQAM systems. As in [START_REF] Farhang-Boroujeny | Multicarrier Modulation With Blind Detection Capability Using Cosine Modulated Filter Banks[END_REF], we consider a one-tap per carrier equalization, since it is particularly adapted to the multicarrier systems. As mentioned in [START_REF] Lin | Convergence analysis of blind equalizer in a cosine modulated filter bank-based multicarrier communication system[END_REF], this type of technique may lead to local solutions and implies that the symbols are detected with a phase shift of 0 or π. The proposed method consists of filtering the equalizer response in order to smooth the channel frequency response. The principle of the proposed technique is simple, since it is similar to a low-pass filter as used in channel estimation [START_REF] Schoukens | The Interpolated Fast Fourier Transform: A Comparative Study[END_REF]. The addition of a differential phase coding allows to suppress any uncertainty on the phase of the recovered symbols.

*CentraleSupélec is the former institution of Vincent Savaux, actually with b<>com research institute at: vincent.savaux@b-com.com.

Furthermore, the proposed techniques largely increases the convergence speed of the CMA. The idea behind CMA applied in OFDM/OQAM is to iteratively update the algorithm by comparing the real part of the output of the equalizer with a given real constant that depends on the constellation size. However, it is worth noting that CMA is equivalent to Sato's algorithm [START_REF] Sato | A Method of Self-Recovering Equalization for Multilevel Amplitude-Modulation Systems[END_REF], and the constant norm algorithm (CNA) [START_REF] Goupil | New Algorithms for Blind Equalization: The Constant Norm Algorithm Family[END_REF], due to the use of real transmitted OQAM symbols. Therefore, this limits the choice of the cost function to those that are adapted to real constellation, whereas numerous cost functions adapted to different complex constellations have been proposed in the scientific literature [START_REF] Goupil | New Algorithms for Blind Equalization: The Constant Norm Algorithm Family[END_REF]- [START_REF] Barbarossa | Blind Equalization Using Cost Function Matched to the Signal Constellation[END_REF].

The rest of the paper is organized as follows: Section II presents the system model and the CMA, and we propose our equalizer with IN cancellation in Section III. Simulations shows the performance of the method in terms of mean square error (MSE) and BER in Section IV, and we draw conclusions in Section V.

II. SYSTEM MODEL

We consider the transmission of OFDM/OQAM symbols over a frequency selective channel, with perfect time and frequency synchronization. At the frequency-time position (m, n), the output of the analysis filter bank (AFB) is written as

y m,n = H m,n x m,n + j (p,q)∈Ω H p,q x p,q < g > m,n p,q Im,n +w m,n , (1)
where x m,n is the real symbol transmitted at frequencytime position (m, n), H m,n is the complex channel frequency response (CFR). A detailed description of the synthesis filter bank (SFB) and the AFB is available in [START_REF] Kofidis | Preamblebased channel estimation in OFDM/OQAM systems: A review[END_REF]. w m,n and I m,n are the complex additive white Gaussian noise (AWGN) and the interference terms whose variances are σ 2 and σ 2 I respectively. The term j < g > m,n p,q is called the intrinsic interference due to the prototype filter g m,n (t) that is given by

g m,n [i] = g(i -n M 2 )e 2jπm M (i- L f -1 2 ) e jφm,n , (2)
where the phase term is φ m,n = (π/2)(m + n) + mnπ as defined in [START_REF] Kofidis | Preamblebased channel estimation in OFDM/OQAM systems: A review[END_REF]. It should be noted that the real symbols x m,n are transmitted with a period M 2 twice shorter than that of OFDM signal where symbols are complex. The used filter has a length L f = KM , where K is the overlapping factor and M is the total number of carriers. Thus, even in absence of the multipath channel, there will be some intercarrier and intersymbol purely imaginary interferences. In fact, g m,n is orthogonal in the real field, which allows to perfectly recover the real OQAM symbols. However, it can be seen in ( 1) that the CFR H p,q induces some distortions from the timefrequency positions (p, q) that belongs to the set Ω 1 and have to be estimated and corrected.

In this paper, the signal is blindly equalized, by using the simple one-tap per-carrier CMA equalizer. To summarize, the aim of the CMA is to tends toward the optimal equalizer coefficient F opt m,n thanks to the update algorithm, which can be written as

F m,n+1 = F m,n -µy * m,n sign(a m,n )(|a m,n | -γ). ( 3
)
where F m,n is the equalizer coefficient, and a m,n = Re(F m,n y m,n ), µ is the step-size parameter, and γ is a realvalued constant which is defined as γ = E{|x| 4 } E{|x| 2 } (see [START_REF] Godard | Self-Recovering Equalization and Carrier Tracking in Two-Dimensional Data Communication Systems[END_REF] for more details). The CMA has a low complexity, but it suffers from two drawbacks: i) it requires a long convergence time, and ii) it is phase blind, i.e. the real equalized symbols a m,n can be detected with a phase equal to 0 or π. In the following, we propose a simple method that reduces the delay of convergence of CMA, and avoids the phase incertainty.

III. CMA BLIND EQUALIZATION WITH IMPULSE NOISE CANCELLATION

It has been shown in [START_REF] Farhang-Boroujeny | Multicarrier Modulation With Blind Detection Capability Using Cosine Modulated Filter Banks[END_REF], [START_REF] Lin | Convergence analysis of blind equalizer in a cosine modulated filter bank-based multicarrier communication system[END_REF] that the CMA applied to the real OQAM symbols converges toward the optimal coefficient value |F opt m,n | with a phase shift equal to 0 or π. The solution that consists of using a differential phase coding as in [START_REF] Farhang-Boroujeny | Multicarrier Modulation With Blind Detection Capability Using Cosine Modulated Filter Banks[END_REF] do not ensure a low bit error rate (BER) value when a large number of received symbols are detected with an opposite phase. In this section, we propose to consider the issue of convergence toward F opt m,n or -F opt m,n from an original point of view, which is described just below.

A. Equalization with Impulsive Noise Cancellation

Fig. 1 shows the real part of a frequency response of the invert of the equalizer coefficients 1/F opt m,n (using CMA) compared with the exact channel frequency response H m,n , knowing that |F opt m,n H m,n | = 1. It can be seen that the frequency positions m=11, 29, and 92 correspond to the opposite minima -F opt m,n . The optimal solutions with opposite phase depicted in Fig. 1 can be seen as an impulsive noise (IN) in the channel frequency response, which can be filtered in a simple way as described in [START_REF] Vaseghi | Advanced Digital Signal Processing and Noise Reduction[END_REF]. The IN cancellation can be seen as a low-pass filtering, but performed in the time domain. This technique is also usually used in channel estimation techniques in order to filter the noise or to perform an interpolation [START_REF] Schoukens | The Interpolated Fast Fourier Transform: A Comparative Study[END_REF]. In 1 Ω is the set of positions p, q that surround the position m, n the considered blind equalization context, the IN cancellation should be seen as a smoothing stage applied to the equalizer coefficients F m,n , and allows full blind equalization as it does not require any pilot for the phase estimation. The proposed equalizer regularly switch to the IN cancellation block as depicted in Fig. 2.

The steps into the "IN cancellation" block in Fig. 2 are described as follows:

1) ∀m ∈ 0, M -1 , the channel frequency response (CFR) samples are obtained by Ĥm,n = 1/F m,n , and stored in the vector Ĥn of size M × 1.

2) The estimated channel impulse response (CIR) is achieved by applying an IFFT: ĥn = IF F T ( Ĥn ).

3) The CIR is multiplied by a rectangular window Π defined by

Π m = 1, if m ≤ L Π 0, else , (4)
such as hn = ĥn Π, where the length of Π denoted by L Π must be at least equal to channel length L. 4) The filtered CFR is obtained by the FFT as Hn = F F T ( hn ) and then the equalizer coefficients by F m,n = 1/ Hm,n . Note that in practice, the L value is unknown, so L Π is chosen in order to overestimate the length of the channel. The advantages of the IN cancellation process are multiple: in addition to suppress the opposite solutions that are seen as an impulsive noise, it also avoids the possible convergence of the coefficients toward local minima 2 . Therefore, it offers the possibility of initialization strategies different from the one proposed in [START_REF] Farhang-Boroujeny | Multicarrier Modulation With Blind Detection Capability Using Cosine Modulated Filter Banks[END_REF]. Moreover, as observed in [START_REF] Schoukens | The Interpolated Fast Fourier Transform: A Comparative Study[END_REF] for the channel estimation, the IN cancellation mitigates the noise and the interferences level. In an OFDM/OQAM context in which the interference term may be not negligible, this smoothing stage may appear to be very convenient.

Also, it is worth highlighting that the proposed IN cancellation process using scheme in Fig. 2 may lead to two solutions: i) all the M coefficients F m,n converge toward their optimal solutions without phase shift F opt m,n , ii) all the coefficients F m,n are equal to the opposite expected solutions -F opt m,n . That is the main difference with the CMA in [START_REF] Bölcskei | A Subspace-Based Approach to Blind Channel Identification in Pulse Shaping OFDM/OQAM Systems[END_REF] where each coefficient F m,n may converge to the optimal value or its opposite. In our solution, if all the coefficients F m,n converge to their opposite optimal solution, the differential phase decoding stage leads to a single error 3 which could be easily suppressed by means of a channel code. If no IN cancellation is performed as in [START_REF] Farhang-Boroujeny | Multicarrier Modulation With Blind Detection Capability Using Cosine Modulated Filter Banks[END_REF], a large number of phases may be erroneous because they are detected with a difference equal to π. In the latter case, the differential phase encoding becomes inefficient and therefore the channel coding as well. In the hereafter section, we will focus on the way of applying the proposed IN cancellation stage in Fig. 2.

B. How to apply the IN Cancellation?

There are three possible ways to perform the IN cancellation stage, which are hereafter explained:

• The IN cancellation is used only one time during the convergence process. It could be a convenient solution if only few F m,n coefficients are detected with an opposite phase (see the example in Fig. 1). However, if a large number of phases are detected with a shift equal to π, a single IN cancellation stage is not sufficient to suppress the appearance of the impulse noise and this solution is similar to the CMA approach. Furthermore, this option raises the issue of the choice of the appropriate starting point to proceed the smoothing stage. • The IN cancellation is employed at each iteration. Here, the update algorithm (3) does not have any sense since it is fed by a coefficient F m,n that is not the one obtained at the previous iteration. In other words, this solution is not anymore a blind equalization, and should lead to unsuitable coefficients F m,n . • The operation is performed with a reasonable regularity regarding to the convergence speed of the update algorithm (3). Thus, the process either performs the IN cancellation as long as optimal solutions with opposite phases remain, or performs a noise cancellation due to the nature of the filter. This solution seems to be the way to carry out the proposed CMA with IN cancellation in Fig. 2.

In order to formalize these three solutions, we define N i the time delay between two consecutive smoothing operations. Thus, the first solution corresponds to N i = +∞, the second one to N i = 1 and the third one to a value of N i conveniently chosen. It is worth noticing that simulations revealed that the first two solutions lead to unsatisfactory results. As a consequence, we only focus on the third solution in the next section.

IV. SIMULATIONS RESULTS

In this section, we present simulations results which show the performance of the proposed blind equalizer. Since we are comparing obtained results with the CMA in (3), we will use the same parameters as described in [START_REF] Farhang-Boroujeny | Multicarrier Modulation With Blind Detection Capability Using Cosine Modulated Filter Banks[END_REF] where the CMA is used in FBMC. We consider a constant multipath channel whose CIR is h = [1, -0.2, 0.3, 0.2, 0.1, 0.2, 0.35, -0.2]. The real part of the corresponding CFR is displayed in Fig. 1. The results are presented for a 4-QAM constellation, i.e. for binary OQAM symbols. As a consequence, note that γ = 1. Each OFDM/OQAM symbol is composed of M = 128 carriers. Two initialisation values are tested: a constant initialization F m,0 = 0.01 for any m = 0, 1, .., M -1 such as proposed in [START_REF] Farhang-Boroujeny | Multicarrier Modulation With Blind Detection Capability Using Cosine Modulated Filter Banks[END_REF], and the sub-optimal initialization strategy F s-opt m,0 = E{xm,0} 2 |ym,0| 2 such as proposed in [START_REF] Savaux | Sub-Optimal Initialization for Blind Equalization with Fast Convergence in OFDM/OQAM Modulation[END_REF]. In both cases, the step-size parameters is set at µ = 0.015. In the following, the performance of the data-aided interference approximation method (IAM) is also plotted for reference. A detailed description of the IAM is provided in [START_REF] Kofidis | Preamblebased channel estimation in OFDM/OQAM systems: A review[END_REF]. The performance of the different blind equalizers is assessed by using the mean square error MSE defined as

MSE = E{(|a m,n | -|x m,n |) 2 }. (5)
In Fig. 3, the CMA is compared with the proposed CMA using IN cancellation (denoted by INC in the figures) for N i = 1 and N i = +∞, where it is performed at the 200th and at the 2000th iterations. The results are obtained for SNR=30 dB, and F m,0 = 0.01. Moreover, achieved MSE values with IAM pilot aided pattern is plotted as a reference.

It can be observed that CMA achieves better performance than CMA with IN cancellation when it is performed with N i = 1 and has a similar MSE than CMA using IN cancellation where N i = +∞. Having N i = 1, the MSE has a "chaotic" behavior, but has a lower bound around -14 dB, whereas the CMA achieves a MSE lower than -17 dB. These observations confirm that the gradient update in (3) does not make sense in the case N i = 1. For N i = +∞, one can observe that once the smoothing stage has been performed (at the 200th or 2000th iteration), the trajectories of CMA using IN cancellation are the same as those with CMA. The peaks that appear at iterations 200 and 2000 can be interpreted as a reset of the equalizer coefficient due to the IN cancellation stage, and this indicates that the option N i = 1 is not suitable when too much coefficients F m,n converge to their opposite solutions. More generally, results presented in Fig. 3 show that the IN cancellation with N i = 1 and N i = +∞ does not improve CMA as proposed in [START_REF] Farhang-Boroujeny | Multicarrier Modulation With Blind Detection Capability Using Cosine Modulated Filter Banks[END_REF]. These results are consistent with the a priori remarks we made in Section III-B. Note that none of the presented methods is able to reach the performance of the pilot-aided estimation using IAM. This is mainly due to the facts that: i) IAM is an accurate channel estimation method, and ii) CMA is known to be a simple method, but achieves weak performance.

In Fig. 4 we compare the MSE performance of the proposed CMA in Fig. 2 for a regular IN mitigation stage performed every N i = 500 iterations. Furthermore, the initialization F s-opt m,0 is compared with the constant F m,0 = 0.01 value. It can be observed that CMA with sub-optimal initialization achieves a gain of 2000 iteration at MSE=-15 dB compared with the "classical" CMA such as proposed in [START_REF] Farhang-Boroujeny | Multicarrier Modulation With Blind Detection Capability Using Cosine Modulated Filter Banks[END_REF]. Fig. 4 also highlights that the CMA using IN cancellation with N i = 500 reaches its lower bound equal to MSE=-22 dB at the 3500th iteration, and when F s-opt m,n is used, this MSE floor is reached at the 1000-th iteration. These results show the capability of the proposed CMA-based equalizer with convenient N i value to have a high convergence speed compared with CMA. Furthermore, it must be noticed that this advantage is accompanied with a good performance, as a MSE gain of 6 dB is achieved by the proposed method compared to the CMA at the 4000-th iteration. Also, it can be seen that CMA with IN cancellation achieves a gain of 2 dB of MSE compared with the pilot-aided method.

For proposed CMA approach with IN cancellation, one can notice that peaks appear at the iterations corresponding to the application of IN cancellation stages. This is due to the fact that the IN cancellations gradually suppress the peaks in the equalizer coefficients. Although the proposed equalizer in Fig. 2 with the initialization in F s-opt m,0 has a remarkable convergence speed, it can be observed that the initialization F s-opt m,0 alone with the CMA allows to get a satisfactory fast convergence as well, such as shown in [START_REF] Savaux | Sub-Optimal Initialization for Blind Equalization with Fast Convergence in OFDM/OQAM Modulation[END_REF]. However, the MSE defined in [START_REF] Schoukens | The Interpolated Fast Fourier Transform: A Comparative Study[END_REF] does not take into account the phase shift of 0 or π. In fact it does not affect the absolute values of the expected equalizer coefficients |F m,opt | but plays a key role in the BER performance, as it is shown afterward.

Fig. 5 depicts the MSE-floor values of CMA, and CMA using IN cancellation with sub-optimal initialization value, and IAM, when the steady-state is reached, versus SNRs. Once again, it can be concluded from Fig. 5 that the blind equalizers achieve better performance than the data aided method (at least for the same channel as in [START_REF] Farhang-Boroujeny | Multicarrier Modulation With Blind Detection Capability Using Cosine Modulated Filter Banks[END_REF]) for SNR values <24 dB using CMA and <40 dB using CMA with IN cancellation (N i = 500). Furthermore, the latter has lower MSE values than the "usual" CMA for SNR>10 dB. This result reflects the fact that the IN cancellation stage is also able to smooth the noise of the channel frequency response (step 3) in the algorithm presented in Section III. Fig. 6 shows the BER performance of the proposed blind equalizer versus the SNR, and also that of CMA approach using both the constant initialization and F s-opt m,0 . The BER is obtained by using a differential phase coding and a convolutional coding with code rate 1/2 and a Viterbi decoder. The BER is calculated at the 4000th iteration, by averaging the results on 10 7 random transmitted bits. The data aided estimation and equalization using IAM is also depicted, and the equalizer using perfect CSI is plotted as a reference.

It can be clearly observed that the proposed technique outperforms the blind equalizer using CMA proposed in [START_REF] Farhang-Boroujeny | Multicarrier Modulation With Blind Detection Capability Using Cosine Modulated Filter Banks[END_REF]. This reflects the fact that a large number of coefficients F m,n converge toward their opposite optimal values. Consequently, when CMA is used alone, a large number of symbols are detected with a phase shift equal to π. In that case, the differential phase decoder and the Viterbi decoder are unable to recover the expected correct binary symbols. With the proposed blind equalizer, this problem is solved since all the symbols are detected with the correct phase thanks to both the differential and Viterbi decoders. Moreover, it should be noticed that a SNR gain of 3 dB is achieved by the CMA with IN cancellation compared with the pilot aided estimation and equalization, which is consistent with the MSE result. Furthermore, the CMA with IN cancellation losses less than 0.2 dB compared with the equalizer with perfect CSI. This is an important result since it shows that blind equalization (once the steady-state is reached) is able in our context to outperform data-aided techniques, in terms of both BER performance and spectral efficiency. Fig. 6 also reveals that the use of the sub-optimal initialization has almost no influence on the BER performance. However, it has been shown that largely improves the convergence speed of the equalizer. From the performance analyzed in both Figs. 4 and6, we conclude that the combination of the proposed method with the initialization F s-opt m,0 leads to i) a fast equalizer, ii) a low MSE floor value, and iii) a high BER performance, compared with CMA alone, and with IAM. Furthermore, the proposed with IN cancellation allows to achieve a full blind equalizer and phase recovery, since no pilot is used.

V. CONCLUSION

In this paper, we presented a simple blind equalization technique for OFDM/OQAM systems. The principle is to combine the CMA and a low-pass filtering in order to increase the convergence speed of the blind equalizer, and to avoid the convergence of the coefficients toward wrong phase values. Simulations revealed that the proposed method achieves performance close to the perfect channel estimation and equalization. Further work will deal with the extension of proposed method to larger constellation sizes.

VI. ACKNOWLEDGMENT

Fig. 1 .

 1 Fig. 1. Real part of H and 1/F showing opposite equalizer coefficient values.

Fig. 2 .

 2 Fig. 2. Proposed equalizer structure with impulse noise cancellation.

Fig. 3 .Fig. 4 .

 34 Fig. 3. MSE of CMA and CMA with impulse noise cancellation (INC) versus OFDM/OQAM symbols for N i = 1 and N i = +∞, SNR=30 dB.

Fig. 5 .

 5 Fig. 5. MSE performance versus SNR. Comparison between CMA, CMA using IN cancellation (ICN) with sub-optimal initialization value, and pilotaided method (IAM).

Fig. 6 .

 6 Fig.6. BER versus SNR of the proposed equalizer compared to[START_REF] Farhang-Boroujeny | Multicarrier Modulation With Blind Detection Capability Using Cosine Modulated Filter Banks[END_REF].

In fact, local minima can also be seen as peaks in the estimated CFR

If all the phases are detected with an error equal to π, the difference of phase between to consecutive symbols is the same as if the phase error is 0, excepted for the first symbol.

This work is supported by the project PROFIL (Evolution of the wideband PROfessional Mobile Radio based on the FILter Bank MultiCarrier modulation) funded by the French national research agency (ANR) with grant agreement code: ANR-13-INFR-0007-03.