Yassine Motie

Alexandre Nketsa

Philippe Truillet

An assistance tool to design interoperable components for co-simulation

Keywords: Complex systems, models, FMI, co-simulation, component generation

 encouraging dialogue between disciplines and reducing errors, costs and Development time. In this article, we describe an assistance tool in order to generate black-box components, facilitating this design task for novices.

Introduction

Designing an interactive system is a difficult task. Designing and evaluating a socalled "complex" system is even more so. From the software point of view, this system can be seen as an integrated set of elements interconnected with each other, in order to satisfy in a given environment one or more pre-defined objectives.

In general, the components of this system include both the facilities, the hardware and software equipment, the data, the services, the qualified personnel and the techniques necessary to achieve, provide and maintain its efficiency. A complex system also has many characteristics such as the heterogeneity of its components, their evolution at different time scales and their geographical distribution integrating both digital systems, physical operators and/or human. These complex systems are usually broken down into subsystems, following either a top-down approach (also known as Stepwise Refinement and functional decomposition), or a bottom-up approach (from existing components that need to be reused). This often leads to a failure for these subsystems to speak the same language to each other and to share information effectively and correctly -to be interoperable -which is one of the major problems that complex systems frequently face.

In order to succeed this collaboration at a global level, it is important to opt for an open environment that allows for a continuous dialogue between the different parties.

One way to do this is by applying co-simulation which is defined as the combination of various simulation models, including other tools, for different components of a complex system. Co-simulation provides the abstraction needed for designers to work on their business expertise. The interest regarding co-simulation, is not only triggered by the coupling of the environments but also by the potential efficiency gain of decoupling a large system model. This is exemplified in [START_REF] Hippmann | Efficient simulation of bush and roller chain drives[END_REF] where a model of an engine is split into subsystems, leading to a decrease of the simulation time by an order of magnitude.

FMI [START_REF] Blochwitz | Functional Mockup Interface 2.0: The Standard for Tool independent Exchange of Simulation Models[END_REF] (Functional Mock-up Interface) standard and data mediation were used in a previous work [START_REF] Motie | A co-simulation framework interoperability for Neo-campus project (regular paper)[END_REF] for the structural and semantic interoperability levels respectively, in order to build an interoperable co-simulation framework applied in neOCampus which is a research project supported by the University of XXX [START_REF] Gleizes | neocampus: A demonstrator of connected, innovative, intelligent and sustainable campus[END_REF]. The idea is that tools generate and exchange models that conform to the FMI specification. Such models are called FMUs (Functional Mock-up Units).

The problem is that the transition from the subsystem studied to an FMU component (black box) requires knowledge of the FMI standard, thus constituting an obstacle for our designers, inhibiting thereby their collaborations. Rather than forcing the users to change their behavior to accommodate our framework, we tried, in this work, to optimize the framework's structural interoperability level around how users can, want, or need to use it. Inspired by the user-centered design approach [START_REF] Abras | User-centered design. Bainbridge, W. Encyclopedia of Human-Computer Interaction[END_REF], we studied components generation methods and proposed a prototype-based on the generation tasks to be performed -for partial automation. The idea is to allow the designers to preserve their tools, their favorite languages and their expertise in order to guide them for the co-simulation first step with other heterogeneous simulators.

Related work

The participatory design cycle begins with the analysis of user needs and activities.

The ISO 16982 standard offers, for example, methods of user observation, questionnaires, interviews, or the study of available documents.

To achieve this phase with designers like those of the neOCampus operation, one should understand their needs to co-simulate (need for data, need to validate their simulation,etc.), their working environments, the tools they use, and their levels of expertise with computer (software, FMIs, GUIs,...). \\ Often in a phase following this first, it is useful to put into practice methods of creativity, such as brainstorming [START_REF]Ergonomics of Human-System Interaction-Usability Methods Supporting Human-Centred Design[END_REF], to produce ideas for solutions. There are variants and more specialized methods like the "Group Elicitation Method" [START_REF] Boy | The group elicitation method for participatory design and usability testing[END_REF] which proposes "brainwriting", a written variant of brainstorming. In our case where it's difficult to bring together the different researchers working in neOCampus, we assumed that this phase isn't that relevant at that moment while we collected the main needs of our users and make sure of the importance of the implementation of a black-box components generation tool, and especially to save the learning effort of the FMI standard. We went then directly to the design phase, then replace this brainstorming by adding a debriefing step as shown in (cf. Fig. 1).

For the creation of solutions, the C and D phases of the process, there are many possibilities. The most common is to use low-fidelity prototypes. These are produced as in our case by us designers from the ideas generated collectively. They are used to present to users solutions to evaluate, validate or refute concepts or interactions, and to choose or propose new ideas. For the realization of these prototypes the designer has the choice between several methods. These are often based on the use of visual content.

Rettig [START_REF] Rettig | Prototyping for tiny fingers[END_REF] and Snyder [START_REF] Snyder | Paper prototyping: The fast and easy way to design and refine user interfaces[END_REF] show, for example, the use of paper prototyping, in which manipulative and discussion interfaces are prepared in the form of drawings or collages. The "Wizard of Oz" experiment [START_REF] Kelley | An iterative design methodology for user-friendly natural language office information applications[END_REF] proposes to simulate the interactive functioning of the final prototype. This methodology is often based on such a visual paper mock-up. Serrano [START_REF] Serrano | OpenWizard: une approche pour la création et l'évaluation rapide de prototypes multimodaux[END_REF] proposes with "Open Wizard" a software solution of magician of Oz for multimodal systems. It allows to simulate input modalities but does not allow the simulation of the output modalities. An alternative to the Wizard of Oz is to code low-fidelity prototypes. According to Sefelin [START_REF] Sefelin | Paper prototyping-what is it good for?: a comparison of paper-and computer-based low-fidelity prototyping[END_REF], the results achieved with these prototypes are equivalent to those obtained with paper models. In addition, the interviews conducted at the end of the tests comparing paper models and software prototypes reveal that 22 of the 24 subjects say they prefer working with software prototypes. New technologies like Adobe Flash or MS Silverlight make it easy to create low-fidelity prototypes.

We coded a low-fidelity prototype which addressed two very interesting points: 1. Adequacy between proposed functionalities and user needs 2. Adequacy between the interface and the users Designers working in the neOCampus project mostly use COTS (commercial off-theshelf) simulation software [START_REF] Boer | Distributed simulation and manufacturing: distributed simulation with cots simulation packages[END_REF] to build and test their simulation models. Integrating these models to form a single meta-model is a major issue especially when distributed simulation technologies are not anchored in these software [START_REF] Taylor | Integrating heterogeneous distributed cots discrete-event simulation packages: an emerging standards-based approach, Systems, Man and Cybernetics, Part A: Systems and Humans[END_REF]. We have seen that because of the lack of communication and collaboration between these different designers, their models were built completely in a disconnected way and do not benefit from the exchange of information that can simplify and accelerate their work. Beyond these practices, another problem identified concerns the difficulty of being able to generate the co-simulation components (essentially components implementing the FMI) from the different simulators used. Indeed, the majority of designers are experts in their field but 1/ are not necessarily computer experts and 2/ if they are, not often expert in co-simulation nor in FMI. This time of taking over the FMI technology (time and practice) has led us to propose a mechanism to support the generation of FMU components based on user practices. It's an interface that is intended to be easy to use for novice users, adapted and allows to accompany the process of generation of components ready to connect to our platform (cf. Many efforts have been made aiming to implement and test the FMI standard since its release. [START_REF] Chen | A generic FMU interface for Modelica[END_REF] discusses a generic interface's implementation and technical problems and challenges helping the importation of FMUs into a simulator. [START_REF] Noll | Implementation of modelisar functional mock-up interfaces in SimulationX[END_REF] describes the implementation of FMI in SimulationX. In [START_REF] Elsheikh | Modelica-enabled rapid prototyping of cyber-physical energy systems via the functional mockup interface[END_REF] an integration strategy for rapid pro-totyping for Modelica models into the FMI standard is presented. Since the modelling and simulation step is done separately by each designer, where the model has to be set up and tested against its specifications first. Standalone executables may be launched, traced, and debugged using additional tools like an IDE (Integrated Development Environment). The next step in the FMU generation process is to determine the interface of the simulation model, which is later exposed through the FMI. This consists of the definition of input and output quantities, or states, as well as internal timing (accuracy and precision required by the simulation model) and external timing (simulation step size for data exchange considerations. These informations are gathered with the information about the FMUs architecture in a modelDescription.xml file, which is connected to the software code using functions usually provided by the FMU SDK (software development kit) [START_REF]FMU SDK: free development kit[END_REF]. Some other issues we are facing is filling-up the gap between the semantics of FMI and the semantics of the source formalism of the diverse calculation models (state machines, discrete event, data flow or timed automata) [START_REF] Tripakis | Bridging the semantic gap between heterogeneous modeling formalisms and fmi[END_REF].

The parts of an FMU are put together after being compiled. A zip file is shaped and the ".dll file" is put in the binaries folder for the interrelated platform. The modelDescription.xml file is set in the root folder. It is then a deliberate choice to put the model source files in the source folder. These steps were essential to identify the automatable steps of those where the actions of the user are essential. Finally, based on these identified tasks, we proposed a "medium fidelity" interface, functional on an identified scenario for which we conducted a pre-experiment.

Preliminary study

As mentioned earlier, we found in surveys that a large majority of users were unfamiliar with FMI technology and had trouble "to jump in". That's why we wanted to check that using an interface, to help with the elaboration of FMU components, made sense for the designers (to understand the process and allow to collaborate) and that it was useful (in terms of time saving in the construction of a component for example, ...).

Participants and procedure

Participants: We conducted this pre-experiment with 7 novice participants and FMI experts. 6 participants were adult men (mean age = 20 years old, standard deviation = 5) and 1 woman (age = 27 years old). Participants were recruited from xxxx, xxxx and xxxx laboratories and were familiar with computers.

Equipment:

We have developed the prototype under Java/QT on a laptop running Windows 10 OS (Core I7, 32 GB RAM, 17" screen)

Procedure: After signing the consent form, we exposed the two tasks to be performed for the pre-experimentation. Task ( 1): an open program in the Eclipse framework was provided to the participants. The program was the same for all participants. They were asked to generate an FMU component using the latest version of the JavaFMI library and following the script available via a tutorial we provided them with (library download, class creation, component generation fmu2).

For the second task (2), we asked them to launch and use the graphical interface, helping with the FMU component's generation, that we developed. This interface proposes to guide the user in the generation process from the number of input and output variables, entering names and types of variables, respect of cast and spellchecker, until the download of the generated component, providing a summary of the performed tasks (cf. Fig. 5).

Fig. 4. screenshot of the proposed interface

Five participants (number P1 to P5) performed the task (1) then (2) and the last two ones (number P6 and P7) performed the task (2) before (1) to counterbalance learning effects. Analysis: We recorded and analyzed the processing time of the different participants either by following the prescribed scenario or by using our interface. The number of actions performed was also recorded and compared to the fact that the participants were FMI experts or not, or if they used to work with integrated development environments (and more specifically Eclipse).

Independent Variables: with interface, without interface. Dependent Variables: time (continuous variable), number of actions (discrete variable), FMI expertise: expert or not (categorical variable), experience related to frequency of use of IDE (categorical variable).

The aim is both to verify our initial hypotheses (the main objective of the tool is to allow non-experts to generate FMU components for co-simulation) and to iterate on our solution. We completed the pre-experimentation phase by debriefing with the participants.

Results and discussion

Users gave very mixed opinions on the experience, but most of them told us that following the script and trying to complete the program, to compile it and run it was a daunting task. In fact, one of these participants decided to leave the experiment after 32 minutes and found that the exercise was difficult.

6 users could reach the end of the task (1) but actually only 3 were able to generate a FMU component (two of which were experts FMI and one who frequently used the IDE and java as the main language in everyday work). The other 3 committed code errors that blocked this generation. On the other hand, the generation of the FMU component using our interface (task 2) was successful by all the participants.

We were able to realize that the order of the tasks did not matter much regarding the results in terms of time of completion of the task. Participants are distinguished according to their level of FMI expertise and their use of IDEs (with scale: 1: never, 2: rarely, 3: regularly, 4: all the time).

There is systematically (cf. Fig. 6) a longer realization time (at least a factor of 3) for the execution of the task (1) than for the task (2) regardless of their expertise. The ANOVA analysis also reveals a significant effect of the experience of an IDE on the time of completion of the task (F(1,8)=50.02, p < 0.001). With only two expert users, we cannot say much about the impact of the expertise on the time of completion of the task. Nevertheless, we can note here again a saving of time of realization by using our interface, independently of the degree of expertise of the subject. The number of actions performed (being analyzed) is also sharply down (by a factor of 8).

Conclusion and future works

This preliminary work enabled us to first make an inventory of the practices of the different actors of the neOCampus project. In order to allow the different experts to communicate and collaborate, we understood that it was preferable that they could first and foremost keep their own practices by allowing them to build or improve their own "business" simulator.

We therefore opted for the interoperability of these heterogeneous simulators based on the FMI co-simulation standard overcoming model semantic gaps and offering them a validation platform for co-simulation. In this process, we have targeted one of the most difficult tasks, which is the generation of pluggable components in our platform, from heterogeneous simulations. We have therefore designed an interface facilitating this generation by focusing on automation with an understanding of the generation process. Hence, we pre-evaluated our interface in order to improve it. Our interface will aim the integration of the total control of the co-simulation with potentially visualization tools adapted to each participant in this co-simulation according to its needs.

The experience was enriching and appreciated. This constantly improving interface will not only be a mechanism to facilitate collaboration and simplify the cosimulation of our different systems, but can also be a plus for the FMI community that in all areas and for different uses is brought to generate components and for some languages and simulation environment for which there is currently no documentation or library to do so.

Fig. 1 .

 1 Fig. 1. Our user-centered design approach

  Fig 2).

Fig. 2 .

 2 Fig. 2. Overview of the FMI component generation process

Fig. 3 .

 3 Fig. 3. low-fidelity models implementing the scenario

Fig. 5 .

 5 Fig. 5. Task completion time for tasks (1) and (2)

Table 1 .

 1 results -task completion time

	Participant	Time (1)	Time (2)	IDEuses	FMI expert
	P1	25'	3'	4	oui
	P2	32' -ab.	5'	4	non
	P3	42'	7'	3	non
	P4	65'	14'	1	non
	P5	54'	8'	2	non
	P6	34'	10'	3	non
	P7	30'	13'	4	oui
	(cf. Table 1.) describes the times performed by the 7 participants (participants P6
	and P7 first performed task (2) before (1)).

We first analyzed the information about an FMU stored in the modelDescription.xml file. For example, the latter contains elements like ModelVariables defining all the variables of the FMU that are visible / accessible via the FMU functions. This inventory made it possible to extract the minimum necessary and sufficient variables for a simulation in order to generate as easily as possible the FMU component allowing co-simulation between systems. Then, we built a task model to understand the steps to be performed, in order to build an FMU component, and proposed several models of which (cf. Fig. 3).


Nicolas Odier
email: nicolas.odier@cerfacs.fr

Thierry Poinsot

Florent Duchaine

Laurent Gicquel

Open Archive Toulouse Archive Ouverte INLET AND OUTLET CHARACTERISTICS BOUNDARY CONDITIONS FOR LARGE EDDY SIMULATIONS OF TURBOMACHINERY

Inlet an outlet boundary conditions are essential elements of any CFD predictions and this is even more so for turbomachinery Large Eddy Simulations, either applied to academic or industrial configurations. For compressible solvers, non-reflecting, characteristic inlet boundary condition imposing total pressure, total temperature and flow direction is usually needed, while an outlet relaxation methodology that automatically adapts the outlet static pressure as a function of the desired mass-flow rate rate is used for turbomachinery flow predictions. Establishing such a framework is clearly desirable especially for industrial use of LES. Development and validations remain necessary in such a fully unsteady context as detailled hereafter.

INTRODUCTION

Available computational ressources have enabled since a few decades a growing litterature in the field of unsteady numerical simulations for complex geometries. With this gain of ressources, Large-Eddy Simulation (LES) is becoming a performant tool to predict losses (Scillitoe et al., [START_REF] Scillitoe | Numerical Investigation of Three-Dimensional Separation in an Axial Flow Compressor : The Influence of Freestream Turbulence Intensity and Endwall Boundary Layer State[END_REF], Legget et al. [START_REF] Leggett | Loss Prediction in an Axial Compressor Cascade at Off-Design Incidences With Free Stream Disturbances Using Large Eddy Simulation[END_REF]), acoustic (Pogorelov et al. [START_REF] Pogorelov | Cutcell method based large-eddy simulation of tip-leakage flow Cut-cell method based large-eddy simulation of tip-leakage flow[END_REF], Tyacke et al. [START_REF] Tyacke | Predictive Large Eddy Simulation for Jet AeroacousticsCurrent Approach and Industrial Application[END_REF], Salas and Moreau [START_REF] Salas | Aeroacoustic Simulations of a Simplified High-Lift Device Accounting for Some Installation Effects[END_REF]), or wall heat transfer (Collado Morata [START_REF] Collado Morata | Effects of free-stream turbulence on high pressure turbine blade heat transfer predicted by structured and unstructured les[END_REF]), that are key topics on which industrials need to focus on to increase efficiency and reduce the emitted noise of next generation of engines. Some comprehensive reviews of LES for turbomachinery may be found in Dufour et al. [START_REF] Dufour | Numerical Investigations in Turbomachinery : A State of the Art Large Eddy Simulation Applications[END_REF], Gourdain et al. [START_REF] Gourdain | Large eddy simulation of flows in industrial compressors: a path from 2015 to 2035[END_REF], Tucker et al. [START_REF] Tucker | Progress in Aerospace Sciences Computation of unsteady turbomachinery flows : Part 1 Progress and challenges[END_REF][START_REF] Tucker | Progress in Aerospace Sciences Computation of unsteady turbomachinery flows : Part 2 LES and hybrids[END_REF][START_REF] Tucker | Unsteady Computational Fluid Dynamics in Aeronautics[END_REF]. When considering Large Eddy Simulation and Direct Numerical Simulation (DNS) in a compressor or turbine stage, several points are to be cautiously addressed: among others, mesh and numerical scheme need to be accurate enough to ensure proper initiation and advection of turbulent structure with as little dissipation and dispersion as possible. Wall-treatment must also be taken into account depending on the final application. Depending on the available computational ressource, wall-resolved LES may help transition of the boundary layer or wall heat transfer prediction in academic cases (Segui et al. [START_REF] Segui | Importance of boundary layer transition in a high-pressure turbine cascade using LES[END_REF], Sayadi and Moin [START_REF] Sayadi | Large eddy simulation of controlled transition to turbulence[END_REF]), but remains too expensive for real configurations where hybrid methodologies are to be considered. In this case, many hybrid RANS-LES methodologies have been proposed: Detached-Eddy Simulations (Spalart et al. [START_REF] Spalart | Comments on the feasibility of LES for wings, and on a hybrid RANS/LES approach[END_REF]), Zonal DES (Deck [START_REF] Deck | Zonal-Detached-Eddy Simulation of the Flow Around a High-Lift Configuration[END_REF]) or Zonal LES (Boudet et al. [START_REF] Boudet | Zonal Large-Eddy Simulation of a Fan Tip-Clearance Flow, With Evidence of Vortex Wandering[END_REF]), Wall-modeled LES (Bose and Park [START_REF] Bose | Wall-Modeled Large-Eddy Simulation for Complex Turbulent Flows[END_REF]). To finish, a last crucial point that is often disregarded in favour to the previous points is to impose acoustically non-reflecting inlet and outlet boundary conditions as well as turbulence injection. A usual way to prescribe the inlet boundary condition for turbomachinery applications is the imposition of a total pressure, total temperature, and flow angle. However, imposing such a piece of information while remaining acoustically non reflecting is rarely detailed. Likewise, the outlet boundary condition usually specifies a static pressure, that must comply with the establishement of a radial equilibrium. Finally, this outlet static pressure must be enforced so that the desired mass-flow rate is established within the turbomachinery, given the inlet total pressure and temperature. Determining this outlet static pressure is not an easy task since it is rarely known in advance, and depends on the general work and losses occuring in the considered configuration. Current academic use of LES usually involves iterative procedures to specify these parameters, but with the growing use of unsteady simulations expected in an industrial context (Larsson and Wang [18], Gourdain et al. [START_REF] Gourdain | Large eddy simulation of flows in industrial compressors: a path from 2015 to 2035[END_REF]), such a strategy will no more be acceptable. This specific issue of inlet and outlet boundary condition specification for LES of turbomachinery configurations in an industrial context is addressed in this paper. After a description of the numerical methodology in the second section, the characteristic formulation required for an inlet turbomachinery condition is presented in the third section, based on Odier et al. [START_REF] Odier | A characteristic inlet boundary condition for compressible, turbulent, multispecies turbo-machinery flows[END_REF]. Finally, a relaxation methodology to automatically prescribe the outlet static pressure is proposed in the fourth section, prior to a set of validation test-cases.

NUMERICAL SOLVER

The present study is performed using an unsteady compressible unstructured and massively parallel solver, AVBP, detailed in Schønfeld et al. [START_REF] Schønfeld | Steady and unsteady flows simulations using the hybrid flow solver avbp[END_REF] and Gourdain et al. [21]. This solver, initially devoted to combustion applications, relies on an explicit Lax-Wendroff [START_REF] Lax | Systems of conservation laws[END_REF] 2nd order accurate in space and time numerical scheme, or a TTG4A 3rd order in space, 4th-order in time numerical scheme (Colin and Rudgyard, [START_REF] Colin | Development of High-Order TaylorGalerkin Schemes for LES[END_REF]). To deal with rotorstator interface in this unsteady context, the MISCOG (Multi Instance Solver Coupled through Overlapping Grids) methodology is considered. This methodology, detailed in Wang et al. [START_REF] Wang | An overset grid method for large eddy simulation of turbomachinery stages[END_REF], consists in coupling two LES computations -one for the static domain, the other for the rotating part -, through the use of an overset grid method. This methodology has proven its ability to ensure a 3rd-order in space numerical scheme trough the interface (De Laborderie et al. [START_REF] De Laborderie | Numerical analysis of a high-order unstructured overset grid method for compressible LES of turbomachinery[END_REF]), and has been applied to various turbomachinery configurations: multistage or centrifugal compressor (De Laborderie et al. [START_REF] De Laborderie | Application of an overset grid method to the large eddy simulation of a high-speed multistage axial compressor[END_REF], Dombard et al. [START_REF] Dombard | Large Eddy Simulations in a Transonic Centrifugal Compressor[END_REF], Wang et al. [START_REF] Wang | LES Investigation of Aerodynamics Performance in an Axial Compressor Stage LES Investigation of Aerodynamics Performance in an Axial Compressor Stage[END_REF]), fans (Leonard et al. [START_REF] Leonard | Large Eddy Simulation of a scale-model turbofan for fan noise source diagnostic[END_REF], Odier et al. [START_REF] Odier | Evaluation of Integral Turbulence Scale Through the Fan Stage of A Turbofan Using Hot Wire Anemometry and Large Eddy Simulation[END_REF], Perez-Arroyo et al. [START_REF] Arroyo | Large Eddy Simulation of a Scalemodel Turbofan for Fan Noise Source Large Eddy Simulation of a Scale-model Turbofan for Fan Noise Source Diagnostic[END_REF]), combustor-turbine configuration (Duchaine et al. [START_REF] Duchaine | On the importance of inlet boundary conditions for aerothermal predictions of turbine stages with large eddy simulation[END_REF], turbine stage (Papadogiannis et al. [START_REF] Papadogiannis | Assessment of the Indirect Combustion Noise Generated in a Transonic High-Pressure Turbine Stage[END_REF]), Wang et al. [START_REF] Wang | Noise mechanisms in a transonic high-pressure turbine stage[END_REF]).

INLET CONDITION: P t , T t AND FLOW DIRECTION PRE-SCRIPTION Formulation

As recalled in the introduction, a methodology able to prescribe the inlet target values without acoustic reflection must be considered for an unsteady compressible turbomachinery simulation. The characteristic methodology first described by Thompson [START_REF] Thompson | Time dependent boundary conditions for hyperbolic systems[END_REF] for Euler equations and then extended to Navier-Stokes equations by Poinsot and Lele [START_REF] Poinsot | Boundary conditions for direct simulations of compressible viscous flows[END_REF] is followed here to address this issue. In this case, the methodology consists in rewriting the Navier-Stokes primitives variables in term of 1D characteristics waves, following Eqs. (1):

      L + L - L t 1 L t 2 L S       =           (u n + c) ∂ u n ∂ n + 1 ρc ∂ P s ∂ n (u n -c) -∂ u n ∂ n + 1 ρc ∂ P s ∂ n u n ∂ u t 1 ∂ n u n ∂ u t 2 ∂ n u n ∂ ρ ∂ n -1 c 2 ∂ P s ∂ n           (1)
where n stands for the normal direction to the boundary surface. For a subsonic inlet, L + L t 1 , L t 2 and L S are respectively the ingoing acoustic wave, the two shear waves, and the entropic wave that have to be determined, while L -, the outgoing acoustic wave, is known from the computational domain. Poinsot and Lele [START_REF] Poinsot | Boundary conditions for direct simulations of compressible viscous flows[END_REF] proposed to consider a locally one-dimensional inviscid flow (LODI) to determine the unknown wave amplitudes depending on the desired condition. In this original work, the wave amplitudes are related to account for the imposition of a velocity magnitude, its direction, static pressure and species composition at the inlet. However, wave relations to impose to-

tal pressure P t = P s 1 + γ-1 2 M 2 γ γ-1 , total temperature T t = T s 1 + γ-1 2 M 2 ,
flow direction, and species composition are not provided. Odier et al. [START_REF] Odier | A characteristic inlet boundary condition for compressible, turbulent, multispecies turbo-machinery flows[END_REF] have recently proposed relations dedicated to impose P t , T t and flow direction. Among the different possible relations, the non-reflecting "Pt-Tt-NSCBC-NR" boundary condition is of specific interest for turbomachinery applications. For this condition, derived for a constant adiabatic coefficient γ, the waves are expressed as

L t 1 = - ∂ u t 1 ∂t (2)
L t 2 = - ∂ u t 2 ∂t (3)
Copyright c 2019 by ASME

L k = Y k L S -ρ ∂Y k ∂t (4)
L + = F 1 ∂ T t ∂t + F 2 ∂ P t ∂t + P t rT t • F 2 F 3 + 1 C p • F 1 F 3 + F 1 • T t r ∑ N k=1 r k ∂Y k ∂t F 4 F 2 + F 5 F 1 ( 5
)
L S = ∂ T t ∂t + 1 C p F 3 + T t r ∑ N k=1 r k ∂Y k ∂t -F 5 • L + F 2 (6)
where Y k is the mass fraction of species k, and

F 1 = e c ρ • P t rT t (7)
F 2 = T t ρ - e c ρC p (8)
F 3 = L t 1 u t 1 + L t 2 u t 2 (9)
F 4 = - ρc 2 • P t P s + P t rT t • (γ -1)e c 2c - u n 2 (10)
F 5 = - (γ -1)T t 2c + 1 C p (γ -1)e c 2c - u n 2 (11)
The inlet flow direction is then enforced through the prescription of sin(φ ) and sin(θ ) where φ and θ are the flow angles depicted in Fig. 1. This choice also enforces the two shear waves accordingly since,

∂ u t 1 ∂t = ∂ ∂t → U × (sin(θ )) = -L t 1 (12) ∂ u t 2 ∂t = ∂ ∂t → U × (sin(φ )) = -L t 2 (13) in which → U = u n 2 + u t 1 2 + u t 2 2
. Note that this velocity magnitude is not imposed by the methodology, only the angular components θ and φ are prescribed. A known issue of the LODI non-reflecting boundary condition is its property to "drift" from the target value of interest mainly because of the underlined hypotheses. To address this issue, Rudy and Strikwerda [START_REF] Rudy | A nonreflecting outflow boundary condition for subsonic Navier-Stokes calculations[END_REF] introduced a linear relaxation method writing any quantity of interest temporal derivative ∂ X ∂t dt under the form:

∂ X ∂t dt = -σ X X predicted -X target (14)
where σ X is a relaxation coefficient and comes as a complement to previous relations. In this present study, the relaxation coefficients for every X variables are chosen to be equals, so

σ X = σ P t = σ T t = σ u t 1 = σ u t 2 .

Convergence towards desired values

The convergence of the average quantities of interest towards the prescribed value at an inlet boundary condition is studied in this subsection for an academic test case, to illustrate the effect of the relaxation coefficient. The considered test-case is depicted in Fig. 2, and consists of an inlet boundary condition where total pressure and temperature are imposed, as well as a normal direction to the boundary surface. The outlet characteristic condition imposes a static pressure, and periodic conditions --�--�-----,-------,------,-, ---- convergence towards the prescribed values, as soon as the re laxation coefficients are not zero. lncreasing the relaxation co efficients increases the convergence speed, but is also likely to increase the resulting acoustic refl ectivity. Odier et al. [START_REF] Odier | A characteristic inlet boundary condition for compressible, turbulent, multispecies turbo-machinery flows[END_REF] how ever showed that reflection remains small for the considered con figuration. Enforcement of large relaxation coefficients is rec ommended only for the flow establishement These ax values should then be reduced as much as possible as soon as the flow is established to ensure fully non-reflective inflow. The important message from these tests is that the unphysical establishement of the flow is a strong function of the infl ow and initial guess which need to be carefully addressed for efficient CPU use of such sim ulations.

. •••••••••••••••••-r ••••••••••••• • •••••• 75000
+•••••••••••••••••••••I••••••••••••••••••

.. ••••

••••••• , •••••••••••••• • •••• • •• , • • • •• • •• • •• • •• •• • • •• • ,• •••••••• ••••••••••• . • • •• • •• • •• • •••••••• •• . •• • ••••• • •••• " • "•"
•••• •••• •••• • < •••• •••• •••• •••• -t- •••• •••• •••• • . ..!•••••••••••••••••••• : : . ••••••••• • ••••••••• r •••••••••••••••••

Turbulence injection

A brief insight of the coupling with turbulence injection is given in the present subsection. The proposed boundary condition may be coupled with a synthetic turbulence injection relying on random flow generation based on the methodologies proposed by Kraichan [START_REF] Kraichnan | Diffusion by a Random Velocity Field[END_REF] then extended by Smirnov et al. [START_REF] Smirnov | Random Flow Generation Technique for Large Eddy Simulations and Particle-Dynamics Modeling[END_REF] for non-homogeneous anisotropic fluctuations generation. The general idea consists in generating a fluctuating velocity field from a sample of random Fourier harmonics, thus considering the turbulent flow as a sum of N periodic Gaussian vortices, so shaped that a Passot-Pouquet [START_REF] Passot | Numerical simulation of compressible homogeneous flows in the turbulent regime[END_REF] energy spectrum is reached in the limit N → ∞. The three resulting unsteady velocity components (u 1 , u 2 , u 3 ) are added to the inlet acoustic waves given in Eqs. ( 2), ( 3) and ( 5) governing the mean flow, so that:

L + turb = L + + ∂ u 1 ∂t (16)
L t 1turb = L t 1 + ∂ u 2 ∂t (17)
L t 2turb = L t 2 + ∂ u 3 ∂t (18)
This methodology follows the VFCBC formalism proposed by Guézennec and Poinsot [START_REF] Guézennec | Acoustically Nonreflecting and Reflecting Boundary Conditions for Vorticity Injection in Compressible Solvers[END_REF] for a boundary imposing the velocity and the static temperature. The resulting turbulent characteristics, as well as the acoustic properties with such a turbulence injection are given in Odier et al. [START_REF] Odier | A characteristic inlet boundary condition for compressible, turbulent, multispecies turbo-machinery flows[END_REF].

OUTLET BOUNDARY CONDITION: P s RELAXATION TO-WARD A DESIRED MASS-FLOW RATE

At a subsonic outlet boundary condition, the only wave to be prescribed is the entering acoustic one L -, since all other waves are known and come from the computational domain. A characteristic formulation imposing a static pressure is detailed in Poinsot and Lele [START_REF] Poinsot | Boundary conditions for direct simulations of compressible viscous flows[END_REF], and reads:

L -= ∂ P ∂t (19)
This formulation has been later improved by Lodato et al. [START_REF] Lodato | Threedimensional boundary conditions for direct and large-eddy simulation of compressible viscous flows[END_REF] and Granet et al. [START_REF] Granet | Comparison of Nonreflecting Outlet Boundary Conditions for Compressible Solvers on Unstructured Grids[END_REF] to account for transverse terms, to ensure vortice convection through the outlet boundary. Finally, Koupper et al. [START_REF] Koupper | Compatibility of Characteristic Boundary Conditions with Radial Equilibrium in Turbomachinery Simu-lations[END_REF] have shown that this characteristic formalism is compatible with turbomachinery applications since the radial equilibrium is naturally satisfied.

In practical LES since a total pressure and temperature are specified at the inlet, the static pressure value at the outlet must be prescribed until the desired mass-flow rate is reached within the turbomachinery. Unfortunately, this pressure value is not known a-priori and depends on the overall pressure losses of the configuration, on the work given or received by the flow, and must be determined manually by the user in agreement with the outlet information. A growing use of LES in an industrial context might be expected in the next decade, to complement RANS predictions or 0D analytical models. This manual determination of the outlet static pressure is thus not satisfying in such a context, where several points of a characteristic line are to be simulated. First examples of characteristic line predictions in a complex industrial configuration using LES may be found in De Laborderie et al. [START_REF] De Laborderie | Analysis of a high-pressure multistage axial compressor at off-design conditions with coarse large eddy simulations[END_REF] and Dombard et al. [START_REF] Dombard | Large Eddy Simulations in a Transonic Centrifugal Compressor[END_REF]. A dedicated relaxation methodology that allows an outlet static pressure adaptation is proposed in this section, as well as a best-practise guideline for its use.

Formulation

The proposed formulation consists in modifying the imposed outlet static pressure every T update periods of time, depending on the average mass-flow rate evaluated during this time interval. This relaxation methodology can be expressed by:

P s corrected = Ps + K( Q -Q re f ) (20)
where P s corrected is the static pressure to be imposed at the outlet condition, Ps is the time-averaged static outlet pressure obtained by averaging the simulation over the duration T update , K is a relaxation coefficient expressed in m -1 s -1 , while Q is the timeaveraged mass-flow rate through the outlet condition provided by the simulation during T update . To finish, Q re f is the target massflow rate specified for this specific condition as defined by 0D evaluation for example.

Subsonic turbine-blade channel

The proposed methodology is tested on a linear cooled turbine cascade depicted in Fig. 6, part of the wind tunnel for linear cascade facility at DLR Göttingen (Aillaud et al. [START_REF] Aillaud | on the Use of Periodic Boundary Condition for Large Eddy Simulation of Trailing Edge Cutback Film Cooling With Internal Ribs[END_REF], Rehder [START_REF] Rehder | Design, manufacture and instrumentation of cascade and cooling blades at dlr (european research project aiteb-2)[END_REF][START_REF] Rehder | Investigation of Trailing Edge Cooling Concepts in a High Pressure Turbine CascadeAerodynamic Experiments and Loss Analysis[END_REF]). The mesh consists in 7 millions cells. The inlet condition is treated using the above characteristic formulation with imposed target values P t = 98000 Pa, T t = 293 K and a flow direction to ensure a zero blade angle of attack. An arbitrary initial uniform solution is considered to start the simulation, with P s = 91140 Pa, T s = 287 K in the whole domain, as well as a zero velocity field. Since the experimental velocity fl uctuation is lower that 1 % (Redher [START_REF] Rehder | Investigation of Trailing Edge Cooling Concepts in a High Pressure Turbine CascadeAerodynamic Experiments and Loss Analysis[END_REF]), no inlet turbulence is added. A subsonic case is targeted in this subsection, with a de sired mass-fl ow rate of Q,ef = 0.1 kg/s. The outlet pressure is updated every Ti,pda1e = 25.0 * 10-5 s which corresponds to 0.5 fl ow through blade passage about, and the relaxation coefficient is chosen to K = 2.0 * 10 5 (m.s)- 1 • The value of Kwas obtained after several tests, to find an appropriate value allowing both suff isant outlet pressure amplitude variations, without too large am plitude that would lead to instabilities. Note that the finding of an appropriate K value is the critical part of the proposed methodol ogy.

• i! , ,••••• • • • • • • • • • • • • • t • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • -• ,
. . • •• • • • • • • • • • • • • • • • -• • . • . . • .
. The outlet static pressure evolution and the mass-flow rate through the outlet boundary temporal evolutions are shown in Fig. 7 and Fig. 8 respectively. Both the actual average pressure on the outlet boundary surface (P s in Eq.( 20)) and the prescribed pressure (Psccrrec,,d in Eq.( 20)) are plotted in Fig. 7, showing that the actual mean pressure quickly follows the prescription and that these two pressures evolve towards Ps = 88900 Pa while the mass-flow rate converges toward Q = Q,ef• At time t = 0.0155 s, the relaxation methodology is no longer considered and a fixed static pressure Ps = 88900 Pa is imposed at the outlet, in order to ensure that the desired mass-fl ow rate Q,ef bas actually been reached without forcing the flow into an unphysical solution. The total computational cost is 5800 hCPU using 240 MPI processes. Following this change of outlet conditions no change is ob served, and the final state resulting from it is shown in Fig. 9, depicting a subsonic flow within the inner stator vane channel as targeted.

Turbine blade with choked channel

The maximum mass-flow rate Q,ef = 0.16 kg/s witbin tbe vane is targeted in this subsection to evaluate the proposed methodology in the case of supersonic configuration. The same inlet boundary conditions and initial field as previously described are considered here, the outlet static pressure is updated ev ery Tupdaze = 12.5 * 10-5 s while K = 2.0 * lü5 (m.s)-1 • The mass-flow rate temporal evolution of tbis specific simulation is reported in Fig. 10, showing that the target value is actually reached at t = 0.06 s about. However, turning our attention to the outlet static pressure evolution in Fig. 11, it can be shown that a plateau is still not reached at that time. Following instant t = 0.0077 s, several tests are performed modifying T,,pdare and K to evaluate their impact on the final state, and show that several final states are likely to occur depending on these parameters.

O.lS

Outlet Massflow evolution -0.16 AvJVV,,,. . . ...... , ... _..-. ... ,-�, -----

� 0.14 •v :;,t O. UH•• •••• •••••••i••••••••••••• ••• i • •••••••••••••••••••• •• ••• •••••• • • ••••••••• • •••• •• ••••••••••••••••i•• •••• ••• • ••••••i••••••••••• •••� i 0.10 q:: o. oeH••••••••••••••'•••••••••••• ••• ' ' •••••••••••••••' '•••• •••••••••••' • •••••••• • ••••• ••'•••••••••••••••'•••••••••••••••'•••••••••••••• � :li 0 . 06 •••• • ••• •••••••'•••••••••••••• • '' • ••••••••••••••' '•••• •••••••••••' • •••••••• • ••••• ••'•••••••••••••••'••• •••• ••• ••••••'•••••••••••• ••• i 0.04 •••• • ••• •••••••••••••••••••••• •• •
•• • ••••••• ••• •••••••••••••••••• ••• ••• •• ••••••••••••••••• ••• o. 02 •••• • •••••••••••••••••••••••• ••• • ••••••••••••••••'••••• ••••••••••' • •••••••• • • ••••• • •••••••••••••• ••••• •••••••• •••••••••••••••••• •••• o.oo�-������-�-�-�-�
0.000 0. 002 0 . 004 0.006 0. 008 0.010 0.012 0.014 0.016 lime (s) 7sooo-----O� u =t l e a = t �Pre=s = ; For comparison, the instantaneous Mach number at t = 0.0077 s is shown in Fig. 12, where a contour of M = 1 is drawn in black. The inner vane is seen to be adapted as expected, induc ing a mass-flow rate independent of the outlet specified pressure. As expected for supersonic outlet conditions, tbe proposed relax ation method becomes inappropriate and the outlet static pressure must be known in advance.

s u = re�•� vo' "' l u = t io�n----� 10000 �-• •••••••••••�•••••••••••• •• •l• ••••••••••••••••••• ••• • •••••••••••••••••••••• •• ' ' ••••••••••••••• ' ' •••••••••••••••; .............. � - 1 i l � 6SOOO ! ••• ••• ••••••• j • .. •••••••••• •• j

Rotor-stator configuration

The proposed overall methodology is finally tested on a 3D rotor-stator configuration representative of an actual industrial problem where the inlet conditions, tbe rotational speed and the mass-flow rate are known, while the outlet pressure is unknown. The transient from initial state where NJan = 0 rpm to the oper ating point where NJan = 11614 rpm is considered in tbis sub section, which is a delicate part of the simulation, especially for compressors where a small modifi cation of tbe outlet pressure may significantly modify the operating point. The confi guration under study is the fan belonging to the DGEN-380 designed by Price Induction, that is a small turbofan of a small private jet engine. The fan bas a bypass ratio of 6, is composed of 14 ro tor blades and 40 stator vanes, and has a diameter of 352 mm. More details about the geometry and the operating points may be found in Garcfa Rosa et al. [START_REF] García Rosa | Experimental Analysis of the Global Performance and the Flow Through a High-Bypass Turbofan in Windmilling Conditions[END_REF]. To alleviate tbe computational cost, the numerical domain has been modified to simulate 42 sta tor vanes, thus ensuring a 360/14 ° periodicity. The stator chord has been scaled to preserve tbe original solidity (Rai and Mada van [START_REF] Rai | Multi-Airfoil Navier-Stokes Simulations of Turbine Rotor-Stator Interaction[END_REF]). The mesh considered in this study is a coarse mesh consisting of 50 millions cells. Wall-laws are used to mode! the boundary layers. More details about the numerical predictions of the average flow and turbulence characteristics on this confi gu-ration using the MISCOG methodology may be found in Odier et al. [30, 5 l]. At the inlet, the proposed characteristic formula tion is used imposing P, = 98803 Pa, I; = 28 l K, and a normal flow direction, white Q,ef = 0.8428 kg/s is targeted at the out let Note that for the entire engine, the total mass-flow rate is Q = 14 * 0.8428 = 11.8 kg/s. At the inlet, the experimental tur bulence intensity is less than 1 %. For this reason, no turbulence is injected in the simulation. The configuration under study is presented in Fig. 13, together with the inlet and outlet boundary condition locations. The outlet relaxation methodology is not considered to start the engine and a fixed static pressure of P s = 100000 Pa is arbitrarily imposed at the outlet for the fi rst rotation. Tests have shown that starting the simulation from NJan = 0 rpm to N 1 Qfl = 11614 rpm with the proposed relaxation methodology is not adequate, since significant transients occur and a shock is likely to occur in the stator vane resutlting in an inappropriate functionning of the proposition as discussed previously. The pro posed relaxation methodology is only considered after an entire rotation, when the flow is established to some arbitrary work ing point. Note that another methodology may be to initialize the simulation with a first RANS solution. Following this ini tial state, the temporal evolution of the outlet pressure and mass fl ow rate through the outlet boundary are respectively depicted in Figs. 14 and15, together with the relaxation coefficient K used along the simulation.

Large relaxation coefficients are first considered to rapidly modify the working point and reduce CPU cost. This results in a rapid mass-flow rate convergence to the target value, as shown in Fig. 15 between rotation n ° l and rotation °2. However, it should be noted that the mass-flow rate evolution should not be the only convergence indicator, and that paying attention to the outlet

•••••••••••••• t •••••••••••••••••••• ............. ••••• ( •••••••••••••••••••• ( ••••••••••••••••••• rtl .
. . Once the outlet mass-flow rate has reached a value close to the target, the relaxation coefficient K may be progressively decreased If K is not decreased, the prescribed outlet pressure is likely to experi ence oscillations, and a reduce in K allows a smooth convergence to the target At rotation n °4, Fig. 14 shows that the outlet static pressure bas nearly reached a converged value, and a constant static pressure P s = 98000 Pa is finally imposed. As shown in Fig. 15, this pressure allows the establishement of the mass-flow rate reaching the final target state. The total computational cost for these 4.5 rotations is 14000 hCPU using 360 MPI processes.

• • • • • • • • • • • • • • • • • r • • • • • • • • • • • • • • • • • • • t-• • • • • • • • • • • • • • • • • • • t-• • • • • • • • • • • • • • • • • • • 6 0000 0 1 2 3 4

Number of rotations

CONCLUSION

Inlet and outlet boundary conditions for compressible and unsteady turbomachinery simulations have been addressed in this paper for its use in an industrial context. A characteristic formulation to impose inlet total pressure, total temperature and flow direction is first presented. A relaxation methodology to automatically adjust the outlet static pressure is then presented and evalutated on a turbine blade configuration. A full rotor-stator configuration is finally studied, taking into account the initialisation phase. A best-practise guideline is proposed to users in order to deal with the relaxation coefficients both for the inlet and the outlet boundary conditions.

ACKNOWLEDGMENT

Dr Pierre Aillaud is acknowledged for sharing his turbine blade configuration and for preliminary code implementation. Dr Nicolás García Rosa and Adrien Thacker from ISAE-Supaero, as well as Price Induction are acknowledged for fruitfull discussions and support regarding the DGEN-380 configuration. The computational resources were provided by the GENCI network, and were part of the allocation no. A0032B10157 (CINES-OCCIGEN).

FIGURE 1 :

 1 FIGURE 1: VELOCITY VECTOR → U , AND CORRESPOND-ING ANGLES θ AND φ .

PeriodicityFIGURE 2 :

 2 FIGURE 2: TEST CASE UNDER STUDY

FIGURE 3 :

 3 FIGURE 3: EVOLUTION OF THE INLET TOTAL PRESSURE P, FOR DIFFERENT VALUES OF RELAXATION COEFFI CIENTS ax

  ............... , .....................

  . � ...................... � ..................... � ................. � ..................... ; ................... .

FIGURE 4 :

 4 FIGURE 4: EVOLUTION OF THE INLET TOTAL TEMPERA TURE T, FOR DIFFERENT VALUES OF RELAXATION CO EFFICIENTS ax

06 FIGURE 5 :

 065 FIGURE 5: EVOLUTION OF THE INLET MACH NUMBER M FOR DIFFERENT VALUES OF RELAXATION COEFFI CIENTS ax

FIGURE 6 :

 6 FIGURE 6: TURBINE BLADE CONFIGURATION.

  100000�-�--�=-=c.c..=..-'---�-� 9 5000 .............. . . ..............•

FIGURE 7 :

 7 FIGURE 7: SUBSONIC TURBINE BLADE: OUTLEf PRES SURE EVOLUTION.

FIGURE 8 :

 8 FIGURE 8: SUBSONIC TURBINE BLADE: OUTLET MASS FLOW EVOLUTION.

FIGURE 9 :

 9 FIGURE 9: SUBSONIC TURBINE BLADE: FINAL STATE.

� T �4 =12.Se ◄ 9 ,FIGURE 10 :

 910 FIGURE 10: CHOIOED TURBINE BLADE PASSAGE: MASS FLOW EVOLUTION.

  ............... , ...... . . ........ , ............... • , ............... , , ... ............ , ............ .. , I 60000 �

FIGURE 11 :

 11 FIGURE 11: CHOIOED TURBINE BLADE PASSAGE: OUT LET PRESSURE EVOLUTION. DASHED LINES: PRE SCRIPTION, FULL LINES: ACTUAL PRESSURR

FIGURE 12 :

 12 FIGURE 12: CHOIOED TURBINE BLADE PASSAGE: IN STANTANEOUS SOLUTION AT t = 0.0077 s.

FIGURE 13 :

 13 FIGURE 13: ROTOR-STATOR CONFIGURATION: DGEN-380. INLET AND OUTLET BOUNDARY CONDIDONS, AND RESPECTIVE PRESCRIPTIONS.

  e;_ 100000 ............ .................... , .................... : .................. ) .................. .. ............. + . . . . . . . . . . . . . . . . . . . . ( . . . . . . . . . . . . . . . . . . . ( . . . . . . . . . . . . . . . . .

FIGURE 14 :FIGURE 15 :

 1415 FIGURE 14: ROTOR-STATOR CONFIGURATION: OUTLET STATIC PRESSURE EVOLUTION

=2.Se-5 s, K =2.e -+s


Amadou Ciss
email: amadou.cisse@univ-lorraine.fr

Mohamed Boutayeb
email: mohamed.boutayeb@univ-lorraine.fr

Backstepping Based State Observer Design Vlasov-Poisson Equation

In this work, we present some new results on state observers design, in infinite dimension, for the 1D × 1D-Vlasov-Poisson transport equation. We show that the observers gains may be deduced from the backstepping technique. The obtained kernel equations are solved by the prime integrals of the ordinary differential system and, by the way, lead to assure convergence of the state observer. To show high performances, numerical simulations of the distribution function and its state observer are presented at the last section.

Introduction

We start by defining the Vlasov-Poisson equation 1D × 1D:

∂ t f (t, x, v) + v∂ x f (t, x, v) -∂ x ϕ(t, x)∂ v f (t, x, v) = 0, (1)
where

• (t, x, v) ∈ Ω T × Ω x × Ω v • f (t, x, v
) is the electron distribution function in the phase space Ω = Ω x × Ω v (with the mass normalized to one and the charge of an electron to plus one) • ϕ(t, x) Electrostatic potential of the E(t,x) field. In the following, we will take

x ∈ Ω x = [0, L], t ∈ Ω T = [0, T ], v ∈ Ω v = [-V c , V c ] with V c > 0.
In order to describe the motion of charged particles in plasmas, we must calculate the force field from the macroscopic density of the particles

ρ (t, x) = Ω v f (t, x, v) dv. (2)
We calculate the force field from the Poisson equation,

-∂ xx ϕ = ρ (t, x) -1, (t, x) ∈ Ω T × Ω x , (3)
with a positive initial data,

f (t 0 , x, v) = f 0 (x, v), (x, v) ∈ Ω x × Ω v . (4)
We impose periodic boundary conditions in x:

ϕ(t, 0) = ϕ(t, L), ∂ x ϕ(t, 0) = ∂ x ϕ(t, L) (5)
and

f (t, 0, v) = f (t, L, v) = q(t, v), (t, v) ∈ Ω T × Ω v . (6)
We assume that q(t, v) is C 2 -function and T -periodic.

A condition of compatibility signifying the "global neutrality" of the plasma is given by

Ω x R f (t, x, v)dvdx = 1. (7)
To uniquely determine the electric field E(t, x), we add a condition of zero mean we have

Ω x E(t, x)dx = 0, t ∈ Ω T , (8)
The output of the system is restricted to x = 0, this is due to the fact that the sensor is placed at the beginning of the membrane. we have :

y(t, v) = f (t, 0, v) = q(t, v), (t, v) ∈ Ω T × Ω v . (9)
T -periodic function q(t, v) ̸ = 0 guarantee the observability of the system (1)-(9).

Objetive

Consider the Luenberger state Observer Vlasov equation:

∂ t f (t, x, v) = -v∂ x f (t, x, v)+∂ x ϕ(t, x)∂ v f (t, x, v)+l(t, x, v)(y-y). (10) Periodic boundary conditions in x f (t, 0, v) = q(t, v) + l 0 (t, v)(y -y), (v, t) ∈ Ω v × Ω T , (11) f (t, L, v) = q(t, v), (v, t) ∈ Ω T × Ω v , (12)
Initial condition

f (t 0 , x, v) = f 0 (x, v), (x, v) ∈ Ω x × Ω v . (13)
While l(t, x, v) and l 0 (t, v) denote the observer gains to be determined.

Introducing the observer error

e(t, x, v) = f (t, x, v) -f (t, x, v)
. The observer error dynamics satisfies

∂ t e(t, x, v) = -v∂ x e(t, x, v) + ∂ x ϕ(t, x)∂ v f -∂ x ϕ(t, x)∂ v f -l(t, x, v)(y -y). ( 14
)
Boundary conditions

e(t, 0, v) = -l 0 (y -y), (15)
e(t, L, v) = 0 (16)
Initial condition

e(t 0 , x, v) = e 0 (x, v), (x, v) ∈ Ω x × Ω v (17)
where ϕ x is given by

-∂ xx ϕ = ρ (t, x) -1, (t, x) ∈ Ω T × Ω x , (18)
Boundary Conditions Poisson equation

ϕ(t, 0) = ϕ(t, L), ∂ x ϕ(t, 0) = ∂ x ϕ(t, L) (19)
The electrostatic field is calculated in an explicit way then one can consider without ambiguity that:

∂ x ϕ(t, x) = ∂ x ϕ(t, x)
The dynamic error (14) becomes

∂ t e(t, x, v) = -v∂ x e(t, x, v)+∂ x ϕ(t, x)∂ v e(x, v, t)-l(t, x, v)(y-y). ( 20
)

Methodology

The idea is to determine l 0 and l by applying the Backstepping method. For that we define the target system of the dynamic error. The choice of the system must be made in such a way that the boundary conditions for the transformation of the invertible volterra integral are verified. It must also guarantee the existence of the integral. Now we consider the following target linear system :

∂ t ẽ(t, x, v) = -v∂ x ẽ(t, x, v) + x∂ v ẽ(t, x, v). (21)
Boundary Conditions

ẽ(t, 0, v) ̸ = 0, (22)
ẽ(t, L, v) = 0. (23)
Initial condition

ẽ(t 0 , x, v) = ẽ0 (x, v), (x, v) ∈ Ω x × Ω v ( 24
)
The system is a transport equation whose resolution is simple with the method of characteristics. The general solution of this equation is:

ẽ = Θ vt -x, x 2 + v 2 , ( 25
)
where Θ is a continuously differentiable function.

With the following lemma We prove that the equations system (21)-( 24) is asymptotically stable in the norm L 2 . Lemma 1 The transport equation (21)-( 24) is exponentially stable in L 2 -norm, if

2V c - γL 2 2 - sup Ω α(t, x, v) 2γ > 0 with V c , γ > 0.

Results and discussion

The backstepping approch is utilized to connect the linear observer error system (20) and a suitable target system in the state ẽ(t, x, v) by the invertible Volterra integral transformation

e(t, x, v) = ẽ(t, x, v) - x 0 k(t, x, ξ, v)ẽ(t, ξ, v)dξ (26)
where k(t, x, ξ, v) is the Kernel transformation. By calling Leibniz's differentiation rule and differentiating the transformation (26) we obtain:

l(t, x, v) = (∂ x ϕ(t, x) + x -v)k(t, x, 0, v) ( 27
)
l 0 is given by (15)

l 0 (t, v) = -1 (28)
We also obtain the following Kernel equation:

∂ t k(t, x, ξ, v) + 2v∂ x k(t, x, ξ, v) -x∂ v k(t, x, ξ, v) = 0 (29)
Choose ξ ∈ (0, L) and x ∈ (ξ, L) we get the Boundary Conditions

(∂ x ϕ(t, x) + x -v)k(t, x, x, v) = 0 (30)
Choose x = L in (26) we have:

k(t, x, L, v) = 0 (31)
then we have:

k(t, L, L, v) = 0 (32) k(t, 0, 0, v) = 0 (33) Initial condition k(t 0 , x, ξ, v) = k 0 (x, ξ, v) (34)
In order to analyze the exponential stability of the observer error dynamicics it is necessary to consider the inverse back-stepping transformation from the e(t, x, v)-system to ẽ(t, x, v)-system in the terms of

ẽ(t, x, v) = e(t, x, v) + x 0 g(t, x, ξ, v)e(t, ξ, v)dξ (35)
where g(t, x, ξ, v) is the transformation Kernel. We proceed in the same way as ( 29)-(33) to determine the kernel-PDE for g(t, x, v) and is given by

∂ t g(t, x, ξ, v) + 2v∂ x g(t, x, ξ, v) + x∂ v g(t, x, ξ, v) = 0 (36)
for ξ ∈ (0, L) and x ∈ (ξ, L) we get the Boundary Conditions

(∂ x ϕ(t, x) -x -v)g(t, x, x, v) = 0 (37)
By taking x = L in (35) we have :

g(t, x, L, v) = 0 (38)
thus we have:

g(t, L, L, v) = 0 (39) g(t, 0, 0, v) = 0 (40)
Initial condition

g(t 0 , x, ξ, v) = g 0 (x, ξ, v) (41)
with the function g 0 (x, ξ, v) being consistent with (38)-(40) at t = t 0 . Based on the solution properties for the two kernels k(t, x, ξ, v) and g(t, x, , ξ, v), the exponential stability of the linearized observer error dynamics with output injection weights l(t, x, v) and l 0 (t, v) can be verified. Theorem 1 Let l(t, x, v) and l 0 (t, v) be as in (27), (28). Then the equilibrium e(t, x, v) ≡ 0 of (14)-( 20) is exponentially stable for all t ∈ R + t 0 in the L 2 -norm.

Numerical results

In this section we present the numerical results of the vlasov-Poison system (1)-(6) and observer (10)-( 13) by the Time-Splitting method (see [START_REF] Afeyan | Son-274 nendr ücker. Simulations of kinetic electrostatic electron nonlinear (keen) waves with275 variable velocity resolution grids and high-order time-splitting[END_REF]). The initial condition is chosen as being a perturbation of a homogeneous equilibrium

f 0 (v) f 0 (x, v) = f 0 (v)(1 + α cos(sx)), (42)
s being the wave number. The equilibrium f 0 (v) can be a Maxwelienne (Landau test case)

f 0 (v) = 1 √ 2π exp(- v 2 2 ). (43)
We suppose then that we have L = 2π/s, s = 0.5 and that the perturbation is equal to α = 0.5 (Landau nonlinear).

The numerical parameters are chosen as follows: the domain Ω v = [-6, 6] and the time step is chosen equal to ∆t = 0.05. The phase domain is meshed by (N x , N v ) with N v = 2N x where N x , N v are respectively the points of the grids of Ω x , Ω v . We will then present below two results of the systems (1) -( 6) and ( 10) -( 13) at times t = 40

Conclusion

In this work we have established the Luenberger observer of the Vlasov-Poisson system 1D × 1D of the phase space.

We have studied the stabilization of the dynamic error of the target linear system (transport equation). The observer output injection weights are determined by make use of a bacstepping-based design.

Figure 1 : 1 √

 11 Figure 1: Solutions f (t, x, v) and f (t, x, v) in phase domain (256, 512) at time t = 40 respectively left and right, l(t, x, v) = ( 1 √ π sin( π 2 x + 2πt) + x -v)x and l 0 (t, v) = -1 .


Zoé Krug
email: zoe.krug@isae.fr

Romain Guillaume
email: romain.guillaume@irit.fr

Olga Battaïa
email: olga.battaia@isae.fr

Lexicographic R * Criterion For Decision Making Under Uncertainty in Reverse Logistics

Keywords: Supply chain, Uncertainty, Reverse logistics

Reverse supply chains (RSC) deal with End-Of-Life products (EOL). They provide the benefits of reducing solid waste, creating new jobs, and generating income from the sales of re-manufactured products and savings in raw materials. At the same time, their design is more challenging in comparison with forward supply chains because of hardly predictable reverse flows of EOL products. Conventional approaches to address the uncertainty are usually risk-oriented and do not consider eventual opportunities. In order to overcome this drawback, we propose to use the R * criterion in the decision making process in a context of complete ignorance with a discrete set of equally possible scenarios. This criterion integrated in a new lexicographic algorithm is used in order to take the Decision-Maker's (DM) behavior into account in different ways regarding if the uncertainty is considered as a risk or as an opportunity. Its performances are evaluated on the case study of facility location problem in reverse logistics.

INTRODUCTION

When information about the distribution of uncertain parameters is not available, it is common to define a set of possible scenarios, each of them corresponding to a possible realization of the uncertain parameters. Robust approach is then applied in the case of a discrete set of equally possible scenarios. Its goal is to maximize the minimum profit over all scenarios. It is very conservative in the way that too much weight is given to the worst case scenario, preventing from finding satisfying solutions when good case scenarios happen [START_REF] Bertsimas | The price of robustness[END_REF]). In this paper, we propose to use a different criterion in order to consider the DM optimism while searching for opportunistic solutions.

Several criteria have been proposed in the literature in order to consider the DM's optimism in the case of discrete scenarios. The most common ones are the Hurwicz criterion [START_REF] Hurwicz | Optimality criteria for decision making under ignorance[END_REF]) and the Ordered Weighted Average (OWA) criterion [START_REF] Yager | Generalized OWA aggregation operators[END_REF]). With the OWA criterion, costs of the scenarios are ranked from the best one to the worst one, then a weighted average is calculated on the ranked costs (the weights depend on the ranking). Hurwicz criterion is a special case of OWA when the best case scenario has a weight equals to α and the worst case scenario has a weight equals to 1-α and all other scenarios weights equal 0. These two criteria are compensatory, i.e. a good case scenario will soften the effect of a bad case scenario and reciprocally.

In order to take into account the psychological evidence that a DM will not react to uncertainty in the same way if considered as a risk or as an opportunity [START_REF] Grabisch | Aggregation on bipolar scales[END_REF]), we are interested in considering the DM optimism with a non compensatory criterion, in order to distinguish areas of risks and opportunity.

In the case of the design of a RSC, the DM has to deal with a lot of uncertainty: on the demand for reprocessed products, the quantity of EOL products returned, the quality of the products, etc... Since the concept of RSC is relatively new, there is a lack of knowledge preventing from efficient probability distribution of possible scenarios. Many studies have already been conducted in this area [START_REF] Govindan | Supply chain network design under uncertainty: A comprehensive review and future research directions[END_REF]), but existing models are most of the time risk-oriented and neglect the opportunities that can exist. However, it is very important to take into account these opportunities in order to build more realistic models and show the benefits the RSC may provide. Recent studies show the positive impact of the reverse logistics on economy. For instance, the McKinsey Center for Business and Environment has demonstrated that the choice of a circular economy could reduce the cost of raw materials in Europe by 2030 by 1,800 billion Euro (25%). Economic benefits can be obtained by the reuse of materials and/or components [START_REF] Toffel | Strategic management of product recovery[END_REF]). Furthermore, the management of EOL products improve the image of the company and protect its market which is an indirect source of additional economic benefits [START_REF] De Brito | Modelling product returns in inventory controlexploring the validity of general assumptions[END_REF]).

In order to help the DM to explore the opportunities of RSC, we develop a new lexicographic R * criterion. This criterion is a "bipolar" operator able to take into account uncertainty in the form of risk but also in the form of opportunity. The classic R * criterion's formulation is given in Section 2. The Lexicographic algorithm with R * is detailed in Section 3. In Section 4, it is applied to a facility location problem for a RSC. The conclusions of the study are presented in Section 5.

BACKGROUND

As aforementioned, there is a psychological evidence that in many cases the DM does not consider uncertainty in the same manner if it is perceived as a risk or an opportunity. To take into account this fact, bipolar operators have been introduced by [START_REF] Grabisch | Aggregation on bipolar scales[END_REF]. They consider a neutral value (for instance 0 or e) below which a score is considered as "bad" and above which a score is considered as "good". Then, "good" and "bad" scores are not aggregated in the same way. [START_REF] Yager | Uninorm aggregation operators[END_REF] presented the family of Uninorm operators which are a subgroup of bipolar operators and are defined as follows:

A Uninorm R is a mapping R : R + × R + → R + having the following properties:

(1) R(a, b) = R(b, a) (Commutative) (2) R(a, b) ≤ R(c, d) if a ≤ c and b ≤ d (Monotonous) (3) R(a, R(b, c)) = R(R(a, b), c) (Associative) (4
) There exists some element e ∈ R + , called the identity element, such that for all

x ∈ R + R(x, e) = x
On [0, e] a Uninorm behaves like a t-norm, and on [e, +∞[ like a t-conorm. We can thus consider that [0, e] is an interval of hazards while [e, +∞[ is an interval of opportunities.

The threshold e therefore corresponds to the optimism level of the DM.

In the considered context of RSC, we are particularly interested in the Uninorm which takes the minimum for tnorm and takes the maximum for t-conorm. The objective is to minimize the risks and maximize the opportunities. This Uninorm is called R * and is defined as following:

Let F (x, s) be the evaluation of the objective function for solution x over a scenario s ∈ S,

R * ((F (x, s)) s∈S , e) = min s∈S F (x, s) if ∃ F (x, s) ≤ e max s∈S F (x, s) otherwise
(1) The solution space of the R * criterion in the case of the maximization of the profit and two scenarios is shown in Fig. 1. f (x, s 1 ) is the profit when scenario s 1 happens and f (x, s 2 ) when scenario s 2 happens. We consider x 1 , x 2 two solutions. The zone of hazards is represented in grey and the zone of opportunities in white. In the first figure, both solutions x 1 and x 2 are in the zone of hazards. In that case, to select the final solution we choose the one with the maximum minimum scenario (here x 2 ) called robust. In the second figure, one solution is in the zone of hazards and another is in the zone of opportunities. In this case, the solution in the zone of opportunities (here x 1 ) will always be preferred. If all the solutions are in the zone of opportunities, the final solution selected is the one with the maximum value of the objective function (or optimistic solution).

f (x, s 1 ) f (x, s 2 ) 0 e e x 1 x 2 • • f (x, s 1 ) f (x, s 2 ) 0 e e x 1 x 2 • • Fig. 1. R * criterion solution space
R * is an interesting non compensatory criterion, but it only takes into account the best case and worst case scenarios and neglects all scenarios in between. It means that two solutions with the same best and worst case scenario won't be differentiated by this criterion. We propose a lexicographic approach of R * in order to fix this drawback.

The lexicographic R * (lexiR * ) approach is the combination of a lexicographic maximum (leximax) approach and of a lexicographic minimum (leximin) approach. The former searches to maximize not only the largest profit but also to maximize the second largest profit, the third largest, and so on. The latter searches to maximize not only the smallest profit but also to maximize the second smallest profit, the third smallest and so on... (see [START_REF] Ogryczak | On the lexicographic minimax approach to location problems[END_REF]; [START_REF] Ogryczak | On solving linear programs with the ordered weighted averaging objective[END_REF]). In lexiR * approach, we use leximin in the areas of hazards and leximax in the areas of opportunities.

In a more formal way, it can be written as follows (we denote S -= {s ∈ S|F (x, s) ≤ e} and S + = S \ S -):

LexiR * ((F (x, s)) s∈S , e) = Leximin s∈S -F (x, s) Then Leximax s∈S + F (x, s) (2)
In order to make the process of resolution more clear for the reader, we present an example below.

Example 1. Let S be a discrete set of scenarios with S = {s1, s2, s3, s4}. Let F (x, s) = (2, 3, 8, 10) (resp F (y, s) = (2, 5, 7, 10)) be a vector of objective values for solution x (resp y) of a decision problem regarding the set of scenario S. In order to lighten the notation, we note F x s = F (x, s) and F y s = F (y, s) in the rest of the example. Our goal is the maximization of the objective. The two vectors cannot be differentiated with R * as they both have the same worst case scenario (i.e. 2) and the same best case scenario (i.e. 10). We will study the resolution of the problem with LexiR * criterion for 3 different cases:

(1) In the first case, we consider e = 12. In that case, the DM is pessimistic and both solutions have all their objective values in the area of hazards. Thus, a leximin is applied in order to choose the solution (see Fig. 2(a)). The evaluation of the objective in the worst case scenario being the same for the two solutions, they cannot be differentiated at the first iteration. Yet, at the second iteration, F y s2 is higher than F x s2 . In this case, solution y is chosen.

(2) In the second case, we consider e = 1. In that case, the DM is optimistic and both solutions have all their objective values in the opportunity area. Thus, a leximax is applied in order to choose the solution (see Fig. 2(b)). The evaluation of the objective in the best case scenario is the same for the two solutions, they are differentiated at the second iteration, where F x s3 is better than F y s3 .

(3) In the last case, we consider e = 6. The DM is neither optimistic nor pessimistic. Two objective values of each solution are in the area of hazards, and two are in the opportunity area (see Fig. 2(c)). In that case, LexiR * starts by considering the risky scenarios with leximin. The evaluation of the objective in the second worst case scenario is in the risky area for both solutions, however, it is better with y than with x, thus, solution y is chosen. If F y s and F x s had been equal for all risky scenarios then the differentiation of the solution would have been on the opportunity area with leximax. 1, the lexicographic approach can be described as follows. Firstly, the model is solved with R * to select the first set of solutions. Then, at the first iteration, the scenario on which the selection of the solution has been made (either the maximum best case scenario or the maximum worst case scenario) is recorded and removed. Then, R * criterion is reapplied to the selected set of solutions with one less scenario. A second set of solutions is thus defined, and the process is repeated until it is possible to differentiate the best solution.

P rof it 0 s 1 s 2 s 3 s 4 e F x s1 /F y s1 F x s2 F x s3 F x s4 /F y
Algorithm 1. Algorithm 1 operates in the following way.

• Step.0 k ← 0 N ← {s 1 , ..., s m } i = 0..k • Step.1 Solve Model (3) max R k + r k (3) S.t (a) R k ≤ F ((x, y), s) ∀s ∈ N, (b) R k ≤ e (c) F ((x, y), s) ≥ -M * Y s + e(1 -Y s ) ∀s ∈ N (d) F ((x, y), s) ≤ e * Y s + (1 -Y s )M ∀s ∈ N (e) r k ≤ (1 -Y s )M ∀s ∈ N (f ) S s=1 δ k s = 1 (g) r k ≤ F ((x, y), s) + (1 -δ k s )M ∀s ∈ N • Step.2 If R k + r k ≤ e then R k = R * k , k=k+1 and go to Step.3. Otherwise, define R k = e, r 1 = r * k , k = 2 and go to Step.5 • Step.3: Solve Model (4) max R k + r k (4) S.t (a) R i ≤ F ((x, y), s) + M (1 -γ i s ) ∀s ∈ N, i ∈ 1..k (b) S s=1 γ i s = S -i ∀i ∈ 1..k (c) R k ≤ e (d) F ((x, y), s) ≥ -M * Y s + e(1 -Y s ) ∀s ∈ N (e) F ((x, y), s) ≤ e * Y s + (1 -Y s )M ∀s ∈ N (f ) r k ≤ (1 -Y s )M ∀s ∈ N (g) S s=1 δ k s = 1 (h) r k ≤ F ((x, y), s) + (1 -δ k s ) * M ∀s ∈ N • Step.4: If S -k = 0 then go to Step.7 Otherwise: If R k ≤ e, then R k = R * k , k = k + 1 and return to Step.3 If R k + r k > e, then R k = e,
Step 0 is the initialization of the algorithm: we define an index k that will be the number of the iteration, a set of all scenarios considered N and an index i.

In

Step 1, the model is solved with the R * criterion (only the best and worst case are considered).

At

Step 2, the cases minimax and maximax are distinguished in order to proceed differently at the second iteration: if the found solution is lower than e, we want to optimize the second worst case scenario; in that case we go to Step 3. On the contrary, if the found solution is higher than e, we want to optimize the second best case scenario, thus we go to Step 5.

Step 3 chooses the second maximum minimum scenario, then the third and so on as long as the solution at the previous iteration stays lower than e.

Step 4 makes a transition between a maxmin case and a maxmax case if the solution of Model.4 at the previous iteration become higher than e.

Step 5 and 6 select the second maximum maximum scenario, then the third and so on and so forth as long as the solution at the previous iteration is higher than e

When there is no more scenario to treat, Algorithm stops in Step 7.

CASE STUDY

We apply our LexiR * criterion to a RSC design problem.

The results obtained are compared with the robust solutions.

Description of the case study

The application is based on a facility -location -transportation problem applied to RSC design for a third party logistic company.

Our goal is to locate the different facilities (collection center, re-manufacturing center and disposal) and optimize the transport between them. The objective is the maximization of the profit. The uncertainty is related

Mathematical model

Indexes, parameters and variables are presented here: In order to simplify the presentation of the model, we introduce the following expressions.

Indexes k = 1..K Index of customers m = 1..M Index of spare markets c = 1..C Index of collection centers r = 1..R Index of re-manufacturing centers d = 1..D Index of disposal sites s = 1..S Index
The total income: it includes all sales revenues. It is scenario dependent and can be formulated as:

Incomes = M m=1 ( R r=1 (RSPm * XRMr,m,s))) (6)
The total operational cost: it includes all production costs, assembling costs, buying costs, dismantling costs or distribution costs from/to all centers of the chain. It is scenario dependent and can be defined as follows:

OpCosts = K k=1 ( C c=1 (Cophc * XKC k,c,s )) + C c=1 ( R r=1
(Cdisr * XCRc,r,s))

+ F f =1 ( C c=1 (Ceco * (XRD r,d,s + XCD c,d,s )))) (7)
The total fixed cost: it is the sum of the set-up costs of facilities:

F ixedCost = C c=1 (CF Cc * Y Cc) + P p=1 (CF Rr * Y Rr) + D d=1 (CF D d * Y D d ) (8)
The total transportation costs: it is the sum of travel costs between connected points of the Supply Chain. It is scenario dependent and can be written as:

T rtCosts = K k=1 ( C c=1 (T C * DKC k,c * XKC k,c,s )) + R r=1 ( C c=1 (T C * DCRc,r * XCRc,r,s)) + D ( C c=1 (T C * DCD c,d * XCD c,d,s )) + R r=1 ( D d=1 (T C * DRD r,d * XRD r,d,s )) + R r=1 ( D d=1 (T C * DRSMr,m * XRMr,m,s)) (9)
Objective The objective is to maximize the total profit calculated as:

T otalP rof its = Incomes -OpCosts -F ixedCost -T rtCosts (10)

Constraints

(1) Constraints (1) to (3) verify that the capacities of all centers are respected. Constraint (4) verifies that the demand is never over-satisfied. However, the demand can remain unsatisfied and considered lost in this case. Constraint (5) calculates the quantity of collected EOL products. Constraints ( 6) and ( 7) calculate the quantity of dismantled, repaired and recycled products. Constraints (8) and ( 9) are the flow balance constraints.

K k=1 XKC k,c,s ≤ CapCc * Y Cc ∀c ∈ C, s ∈ S (2)

Results: Robust model versus LexiR * model

In order to show the advantages of our model, we compare it with a Robust model, and with a classic R * model as presented in Section 2. The results of the numerical experiment are reported in this section. The scale of the numerical experiment is as follows: the number of potential locations for establishing the collection centers, re-manufacturing centers and disposal is 10, the number of customers and spare market customers is 10. The data used for the resolution are presented in Table 2. The transportation costs are defined per product and per 1 kilometer.

We define four different scenarios from the worst case to the best case. They are presented in Table 3. The results obtained show that:

-If e is equal to or above the robust "MinMax" solution, the found solution is the same with the Robust model and with the R * model while the LexiR * model improves the two middle case scenarios compared to the robust solution. Thus, new opportunities are disclosed without any risks for the DM. -With the R * and the LexiR * models, the more e decreases compared to the "MinMax" solution, the more the two worst case scenarios are degraded but in return, opportunities are revealed in two best case scenarios. The gains are much more important than the losses in all considered cases and are slightly superior with the LexiR * model than with the R * model. By choosing the value for e the DM accepts to take a little more risks, but it is possible to keep the control over the robustness of the solution while the opportunity area is much better explored than with the Robust model. -The LexiR * model gives better results for two worst case scenarios compared to the R * model. Therefore, the risks taken by the DM are better controlled.

CONCLUSION

The implementation of RSC is essential in order to reduce waste. A successful RSC design relies on appropriate modeling of uncertainty not only in terms of risks, but also opportunities. In this paper, we suggest a new modeling approach using LexiR * to take DM optimism into account in both hazard and opportunity zones. The results obtained show that the use of the LexiR * criterion makes it possible to better explore the opportunity zone compared to a robust model without loosing control over robustness.

It also provides the possibility to better control the risks taken in comparison to the R * model without lexicographic approach.

This study reveals many research perspectives. For instance, we are planning to examine the case of a discrete set of scenarios with imprecise probabilities. The presented algorithm could also be improved by considering a fuzzy value for e. Its extension to a continuous set of scenarios should also be examined.

Fig. 2 .

 2 Fig. 2. Examples of resolution with LexiR * criterion

••

  r 1 = r * k define k = 2 and go to Step.5 StepR i ≤ F ((x, y), s) + M (1 -γ i s ) ∀s ∈ N, ∀i ∈ 1r i ≤ F ((x, y), s) + M (1 -δ i s ) ∀s ∈ N, ∀i ∈ 1Step.6: If S -k -k -1 = 0 then go to Step.7 Otherwise r k = r * k , k = k + 1, and go to Step.5• Step.7: Stop.

Fig. 3 .

 3 Fig. 3. Reverse Supply Chain to the quantity of returned EOL products, the quantity of demand for re-manufactured EOL products and the reprocessing time of EOL products.

  ,s * T remanus ≤ CapRr * Y Rr ∀r ∈ R, s ∈ S

Table 1 .

 1 NotationsWith the use of notations in Table

	Notation	Description	Notation	Description
	R k	Linearization of min	Ys	Binary variable
	r k k	Linearization of max Number of iteration	δ k s γ k s	Binary variable Binary variable
	e	Risk threshold	M	Big value

  The numerical experiment was conducted with IBM-ILOG CPLEX 12.6.3 on an Intel Core 2.60 gigahertz machine with 15 gigabyte RAM. The results from R * are presented in Table4and the results from LexiR * are presented in Table5. The two tables are built in the same way: Column 1 is the model used to solve the problem, Column 2 is the chosen value for e when the R * or the LexiR * model

	Parameter	Value
	CapC	Uniform(55000,200000)
	CapR	Uniform(28000,56000)
	Distances	Uniform(0,500)
	Coph	Uniform(1,6)
	Cdis	Uniform(2,8)
	Ceco	1
	T C	0.003
	Rd	80%
	Rr	80%
	RSP	Uniform(90,120)
	CF C	Uniform(10000,90000)
	CF R	Uniform(40000,60000)
	CF D	Uniform(40000,60000)

Table 2 .

 2 Deterministic Parameters

	Parameter	s 1	s 2	s 3	s 4
	D ( * 1000)	U(22,27.5) U(12,17.5) U(22,27.5) U(22,27.5)
	R ( * 1000)	U(12,17.5) U(12,17.5) U(22,27.5) U(22,27.5)
	T remanu	U(1,2)	U(5,6)	U(5,6)	U(1,2)

Table 3 .

 3 Uncertain Parameters are used (they are indicated in relative percentage of the "MaxMin" solution). Column 3 to 6 are the results for each scenario. Every scenario for which the R * or the LexiR * model brings an improvement compared to the Robust model are colored in green and every scenario for which the R * or the LexiR * model brings a degradation are colored in red.

	Model	e	s 1 (e)	s 2 (e)	s 3 (e)	s 4 (e)
	Robust -	1552422 1899237 1969349 1981001
	R *	-0% 1552422 1899237 1969349 1981001
	R *	-3% 1538872 1896560 2662113 2902355
	R *	-5% 1482802 1666845 3204561 3212546

Table 4 .

 4 Compared Profit between Robust approach and R * approach for the case of 4 scenarios

	Model	e	s 1 (e)	s 2 (e)	s 3 (e)	s 4 (e)
	Robust	-	1552422 1899237 1969349 1981001
	LexiR *	+3% 1552422 1899994 1970978 1981001
	LexiR *	-0%	1552422 1899994 1970978 1981001
	LexiR *	-3%	1538875 1897844 2663437 2902355
	LexiR *	-5%	1482862 1669413 3211150 3212546

Table 5 .

 5 Compared Profit between Robust and LexiR

* approach for the case of 4 scenarios