Florence Podwojewski
email: florence.podwojewski@ifsttar.fr

Mélanie Ottenio

Philippe Beillas

Gaëtan Guerin

Frédéric Turquier

D Mitton
email: david.mitton@ifsttar.fr

M Otténio

G Guérin

Mechanical response of animal abdominal walls in vitro: Evaluation of the influence of a hernia defect and a repair with a mesh implanted intraperitoneally

Keywords: Biomechanics, mesh repair, digital image correlation, pressure, contact loading

come

Mechanical response of animal abdominal walls in vitro:

Evaluation of the influence of a hernia defect and a repair with a mesh implanted intraperitoneally

Introduction

Treatment of incisional hernias is a common surgical procedure [START_REF] Cobb | Incisional herniorrhaphy with intraperitoneal composite mesh: a report of 95 cases[END_REF]. The laparoscopic repair with mesh implantation improves the treatment of incisional hernia by reducing the rate of recurrence and risk of wound complication (Cobb et al., 2005a). The use of composite mesh, that limits adhesions with abdominal organs, makes the laparoscopic repair easier by putting the mesh in the intraperitoneal location [START_REF] Cobb | Incisional herniorrhaphy with intraperitoneal composite mesh: a report of 95 cases[END_REF]. However, problems of recurrence can still happen. Problems of restriction of the abdominal wall mobility or pains may also occur to patients [START_REF] Mclanahan | Retrorectus prosthetic mesh repair of midline abdominal hernia[END_REF][START_REF] Müller | Abdominal wall compliance after Marlex mesh implantation for incisional hernia repair[END_REF].

Studies on the interaction of mesh with the abdominal wall have been performed. Factors identified as possible causes of recurrence include the size of the implanted mesh compared to the size of the defect (the overlap), and the method of mesh fixation [START_REF] Binnebösel | Biomechanical analyses of overlap and mesh dislocation in an incisional hernia model in vitro[END_REF][START_REF] Schwab | Biomechanical analyses of mesh fixation in TAPP and TEP hernia repair[END_REF]. It has also been suggested that the mesh structure should be improved in order to have similar biomechanical behaviour to the abdominal wall to reduce discomforts (Konerding et al., 2011a;[START_REF] Hernández | Mechanical and histological characterization of the abdominal muscle. A previous step to modelling hernia surgery[END_REF].

Numerical models could be useful to assess the influence of a mesh on the behaviour of the abdominal wall. In particular, such models could be used to evaluate the possible effects of mesh design changes or implantation procedures very early in the design process. However, experimental data on the mechanical response of the abdominal wall are a prerequisite to develop models. [START_REF] Hernández | Mechanical and histological characterization of the abdominal muscle. A previous step to modelling hernia surgery[END_REF] conducted uniaxial tensile tests on the flat abdominal muscles of rabbits in order to characterize the passive behaviour of the abdominal wall and develop a numerical model of a healthy and a partially herniated repaired abdominal wall of rabbit [START_REF] Hernández-Gascón | Mechanical behaviour of synthetic surgical meshes: Finite element simulation of the herniated abdominal wall[END_REF]. [START_REF] Förstemann | Forces and deformations of the abdominal wall -a mechanical and geometrical approach to the linea alba[END_REF] performed uniaxial tensile tests on the linea alba to develop an analytical model of the abdominal wall. While these studies can provide useful data to help in building a model, they cannot be used to evaluate or validate the performance of a whole abdominal wall model subjected to mechanical load. Thus, Konerding et al. (2011a) performed a human cadaver study to validate the model proposed by [START_REF] Förstemann | Forces and deformations of the abdominal wall -a mechanical and geometrical approach to the linea alba[END_REF]. [START_REF] Song | Elasticity of the living abdominal wall in laparoscopic surgery[END_REF] and [START_REF] Szymczak | Investigation of abdomen surface deformation due to life excitation: Implications for implant selection and orientation in laparoscopic ventral hernia repair[END_REF] studied the complete human abdominal wall response respectively to pressure during a laparoscopic surgery and during different physical exercises. While such studies provide a useful starting point for the validation of model, it only included data on the response of the external surface of the intact human abdominal wall . This is an important limitation for the validation of a model, whose purpose is to study the interaction of the mesh with the wall when implanted on the internal surface. No study comparing the behaviour of the external wall (skin side) and the internal wall (peritoneal side) could be found in the literature.

Furthermore, many studies provide experimental results on the mechanical properties of the complete abdominal wall [START_REF] Song | Elasticity of the living abdominal wall in laparoscopic surgery[END_REF]Konerding et al., 2011a;[START_REF] Junge | Elasticity of the anterior abdominal wall and impact for reparation of incisional hernias using mesh implants[END_REF] and its components [START_REF] Hernández | Mechanical and histological characterization of the abdominal muscle. A previous step to modelling hernia surgery[END_REF][START_REF] Förstemann | Forces and deformations of the abdominal wall -a mechanical and geometrical approach to the linea alba[END_REF][START_REF] Rath | The abdominal linea alba: an anatomo-radiologic and biomechanical study[END_REF][START_REF] Hollinsky | Measurement of the tensile strength of the ventral abdominal wall in comparison with scar tissue[END_REF], or on the characterization of meshes only [START_REF] Junge | Elasticity of the anterior abdominal wall and impact for reparation of incisional hernias using mesh implants[END_REF][START_REF] Hernández-Gascón | Mechanical behaviour of synthetic surgical meshes: Finite element simulation of the herniated abdominal wall[END_REF]. However, only a few studies were found regarding the mechanical response of the abdominal wall with an implanted mesh, and how it may differ from the response of a healthy wall. [START_REF] Müller | Abdominal wall compliance after Marlex mesh implantation for incisional hernia repair[END_REF] studied the difference of abdominal wall mobility between a healthy control group and a group of patients with incisional hernia, who was treated by mesh implantation. [START_REF] Hernández-Gascón | Mechanical behaviour of synthetic surgical meshes: Finite element simulation of the herniated abdominal wall[END_REF] observed the mechanical response of intact abdominal wall and abdominal wall that was repaired using different meshes on rabbits. This study was performed on small animals and only the lateral part of the abdominal wall was considered. Also, samples were characterized by uniaxial test that do not reproduce the physiological loading conditions of the abdominal wall. To the authors' knowledge, no studies compared the behaviour of a same intact abdominal wall, and then implanted with a mesh.

Therefore, the current study focused on the anterolateral abdominal wall, where incisional hernias may occur. The objective of this work was to study the response of a same abdominal wall under three different states: intact, with a hernia defect, and repaired. Particular attention was paid to study simultaneously the internal and external surfaces. Boundary and loading conditions were defined to ensure reproducibility and control.

Materials and methods

Specimen & preparation

Anterolateral abdominal walls of six female pigs, aged 4 to 5 months and weighing about 45 kg, were used for the current study. The abdominal walls were removed from the animals less than 30 minutes after euthanasia at the VetAgro Sup, Veterinary Campus of Lyon (Marcy l'Etoile, France), and then kept frozen at -20°C until testing. The abdominal walls were cut along the xiphoid process and the costal margins and along the pubic bones and the iliac crests. The lateral incisions were done between the iliac spines and the lower part of the rib cage. Thus the abdominal wall had a triangular shape (Fig. 1). All the layers were preserved: muscles, aponeuroses, adipose tissue, skin and peritoneum. The abdominal walls were thawed at room temperature 16 hours before the test. Just before testing, the external surface of the abdominal wall was shaved. The thickness of each abdominal wall was also measured (Fig. 1) and their average thickness is reported in Table 1. The specimens were sprayed with saline solution to enable their hydration.

Experimental setup

First, the abdominal wall was placed on a hemispherical support (diameter of 9 cm) in order to induce a curvature. Then, it was put on an aluminum plate with a triangular hole exposing the anterolateral abdominal wall. A rubber sheet with the same hole was added to cover the section of the abdominal wall outside the hole. This sheet was clamped using custom designed clamps positioned all around the hole (Fig. 2a). As the abdominal wall thickness was not constant over its circumference, the clamps were adjusted to provide adequate tightening in over the whole circumference and prevent local sliding during loading. The abdominal wall was not removed from this fixture until the very end of the experiment.

Then, the abdominal wall was positioned on a custom designed aluminium table mounted on a testing machine (INSTRON 8802, High Wycombe, England) (Fig. 2b). The external face of the wall was directed downwards, leading to a natural curvature due to gravity.

Physiologically, the abdominal wall can be loaded by contact of the hollow organs or by pressure in case of laparoscopic surgery. So, during the experiment, the abdominal wall was submitted to two mechanical loading cases: contact and air pressure.

First, the air pressure was applied on the internal surface (Fig. 3a). A Plexiglas plate was mounted on the top of the aluminium plate in order to create a closed cavity. The Plexiglas plate was transparent, allowing to observe the internal abdominal wall. Compressed air and a manual valve were used to control the pressure in the cavity. At a laparoscopic pressure of 12 mmHg, a smaller displacement of the abdominal wall than the one measured by [START_REF] Song | Elasticity of the living abdominal wall in laparoscopic surgery[END_REF] was observed due to strong boundary conditions. So, the pressure was increased until it reaches 50 mmHg, which is in the physiological range of intra-abdominal pressure (Cobb et al., 2005b). Then, it was let to return to atmospheric pressure by opening the valve. The pressure loading cycle was repeated six times. The first 5 cycles enabled to precondition the wall in order to reach a steady state, limiting variability to assess various conditions. Only the sixth cycle of pressurization was used for the data analysis.

Then, a contact loading was applied to the internal side using a rigid sphere, to better control the loading conditions for further numerical model (Fig. 3b). After the pressurization test, the Plexiglas plate was removed and the abdominal wall was loaded directly by the rigid Plexiglas sphere (12 cm in diameter) fixed to the actuator of the testing machine. First, the sphere was moved to get in contact with the abdominal wall. Then, the specimen was preconditioned with 5 cycles of 20 mm amplitude at a frequency of 0.5 Hz, close to the respiratory frequency. A displacement of 35 mm, which is in the range of physiological displacement of the abdominal wall [START_REF] Müller | Abdominal wall compliance after Marlex mesh implantation for incisional hernia repair[END_REF], was finally applied at a velocity of 80 mm per minute. The velocity was slow enough to characterize the quasi-static response and also led to a short duration of test that limit the quite long duration of all the protocol.

An incision was made in the middle of the linea alba over a length of 5 cm, which is the size of a medium incisional hernia. The skin was kept intact. The incision was filled with Vaseline and covered with a latex film to avoid air leaks below the skin during pressure loading. The incised abdominal wall was re-loaded. Then a surgical repair was performed using a 10*15 cm² Parietex ® Composite mesh, centred on the defect and fixed on the wall with 20 tackers (AbsorbaTack ® ). The tackers were located at one centimetre of the border of the implant and spaced apart 2 cm. A plastic film (15*20 cm², thickness 10μm) was put on the mesh to avoid air infiltrations between the implant and the peritoneum. Finally the repaired abdominal wall was loaded in the same conditions.

In summary, the abdominal walls were loaded consecutively by pressure and contact in the three following states:

-intact abdominal wall, -incised abdominal wall (after creating a defect simulating an incisional hernia), -repaired abdominal wall (finally after surgical repair by mesh implantation).

Testing duration influence on an intact abdominal wall behaviour

Beforehand, a healthy abdominal wall was tested to assess the effect of time and of the loading sequence on the abdominal wall response. An abdominal wall was subjected to the protocol twice within 4 hours, by switching the different loading cases (pressure, contact and pressure loading again). Between the two series of tests, the specimen was covered with moist gauzes.

Measurements

The pressure and the force applied by the sphere to the abdominal wall were measured by a 7 bar ENTRAN EPX-N02 pressure sensor and a 1000 N INSTRON force sensor (accuracy 0.5%), respectively.

Four synchronized SA3 PHOTRON black and white video cameras (Tokyo, Japan) were used to record videos of the abdominal wall during the deformation. Two cameras equipped with 35 mm Zeiss lenses (Oberkochen, Germany) were set on the internal surface, and the two others (equipped with 24-70mm Sigma lenses (Tokyo, Japan) in 24 mm position) were set on the external surface. The resolution of the cameras was 1024 by 1024 pixels, which led to approximately 3 pixels per mm in the area of interest. The acquisition frequency was 10 frames per second.

Before testing, the two faces of the abdominal wall were covered with white make-up to make the background uniform. Then a random speckle pattern was applied with a black paint spray .

The speckle enabled to determine 3D local displacement and 3D strain fields by digital image correlation using the VIC3D ® stereo-correlation software (Correlated Solution, South Carolina, USA).

Data analysis

For the contact loading case, the force versus displacement curves were plotted for the last (6 th ) loading cycle, and a stiffness in N/mm was worked on. This stiffness was calculated as the slope of the curve determined by linear regression between 26 and 30 mm of displacement (Fig. 4). This interval of displacement defined a common quasi linear area for all the tested abdominal walls. The relative differences of stiffness were calculated between two states, 100 S S S state state state s s

     1 1 2 1 2 tate tate
, where S is the stiffness.

The displacement and strain fields were obtained by stereo correlation for all configurations and for both loading cases on the external wall surface. The strain fields could only be computed on the intact internal wall surface for pressure loading because the Plexiglas sphere hid the internal wall surface during the test. An average strain value was calculated in the central area of the abdominal wall (~70mm*105mm between the nipples) from the Lagrange strain values computed by Vic 3D in first principal direction (E1) (Fig. 5).Average strains along longitudinal (linea alba) and transverse lines of the external abdominal wall surface were also calculated for the pressure loading case.

Statistical analysis

A statistical analysis was performed in order to assess the influence of the state (intact, incised, and repaired) of the abdominal wall on its mechanical response. A Wilcoxon nonparametric test (Wilcoxon, 1945) for paired samples was used. A value of p < 0.05 was selected to indicate statistical significance. Parameters studied are the stiffness and the average strain. Statistical analysis was performed using Unistat® software (London, England).

Results

Testing duration influence on an intact abdominal wall behaviour

For both mechanical loading, the mechanical response of the abdominal wall was found to be non-linear, with a toe zone followed by a more linear zone . The unloading curve did not overlap with the loading curve. As an example, the force-displacement curves for contact loading are given Fig. 4. A similar mechanical response was obtained whatever the time of test for the pressure and the contact loading, with respectively difference of stiffness of about 5% and 10%.

Relation between internal and external surface strains in the intact abdominal wall

The strain fields of the internal and the external abdominal wall surfaces exhibited different patterns as illustrated on Fig. 5. The location of the maximum strain on the internal surface did not match the location of the maximum strain on the external surface. On the internal surface, a region of greater strain near the long edge was observed. On the external surface, there was less strain near the edges and especially at the three corners of the triangle.

Table 1 displays the mean value and the standard deviation of the local strains in the central area for the internal and external surfaces for the pressure loading case at 50 mmHg. The average strain in the central area of the external surface (13.7 (2.1) %) was almost 2.6 higher than the mean strain of the internal surface (5.3 (0.7) %), with variation of 23.3%.

Influence of a defect and of mesh repair on the behaviour of the abdominal wall

The curves of the contact loading case for the three studied states are displayed on Fig. 6.

Mean stiffness calculated for each state are presented in Fig. 7. The difference of stiffness between intact walls and incised walls was statistically significant (p = 0.03). However, the stiffness obtained for intact and repaired states were similar (p = 0.43). Relative differences between the abdominal walls states are summarized Table 2. The stiffness of the incised cases was lower by about 25% than the stiffness of intact cases. There was a smaller relative difference of stiffness between intact and repaired walls; the stiffness of intact cases was 6% greater than the stiffness of repaired cases.

The mean strains values computed in the central area for the external surfaces of the abdominal wall for each configuration and loading case are displayed in Table 3. For both pressure and contact loading cases, the defect had a significant effect on the strain of the abdominal wall (p = 0.03). The defect increased the average strain of about 74% and 35% respectively for the pressure and contact loading. For both loading cases , significant differences of strains were found between incised and repaired states (p = 0.03). For the pressure loading case, significant differences of strains were found betweenintact and repaired states (p = 0.03) whereas no significant differences were found for the contact loading case (p = 0.81). More elongation was observed in the longitudinal direction than in the transverse direction for each state (e.g. for the repaired state: 15.9 (3.0) % and 15.3 (2.8) % respectively)

The defect doubled the average strains along the two directions. No differences were found between incised and repaired states.

Discussion

In this study, an experimental protocol was developed to characterize the mechanical response of animal passive abdominal wall in three states (intact, incised and repaired) and for two loading cases (increased air pressure similar to a laparoscopic procedure and contact loading).

The test of an intact abdominal wall on a day showed that the evolution of the mechanical response of a wall is low (relative difference of 10% for the stiffness). Thus, the duration of the test on a same abdominal wall has limited influence on its behaviour. So, it is expected that differences observed for each configuration can be attributed to the change of state and not to the duration of the protocol.

This study was conducted on porcine abdominal walls. An animal model was used to develop the protocol since the anatomical samples were easier to access than those of human. Pig was chosen, because the dimensions of its abdominal wall are close to those of human. Swine are also often used in biomedical research since they have large anatomical and physiological similarities with humans [START_REF] Swindle | Comparative anatomy and physiology of the pig[END_REF]. In addition, the experiments on porcine abdominal walls allowed us to validate the repeatability of our protocol. For six specimens with small inter-individual variability (pigs of same weight, about 45kg, and of same age, between 4-5 months) and with very constrained testing boundary conditions, relatively small variability of the stiffness data was observed, with variations in the order of 12.3% for the intact state, 16.5% for the incised state, and 6 % for the repaired state.

This study also shows the influence of a defect and a mesh repair on the behaviour of the abdominal wall. The defect decreased the stiffness of the intact abdominal wall by 25% on average and increased the average strain of the abdominal wall by about 74% and 35% respectively for the pressure and contact loading cases. The presence of the mesh in immediate postoperative configuration on an incised abdominal wall enabled to restore the stiffness to values that were not significantly different from those of an intact wall. Other studies interested in the repaired abdominal wall, but when the mesh is integrated. Konerding et al. (2011b) compared the response of different types of meshes that are integrated to the abdominal wall without referring to intact abdominal wall behaviour. This is not the case for the study of [START_REF] Hernández-Gascón | Long-term anisotropic mechanical response of surgical meshes used to repair abdominal wall defects[END_REF] that assessed the response of the repaired abdominal wall and compared it to the response of intact abdominal wall. However, the response was obtained by uniaxial test.

For the first time, this study provides data on the strain fields of the abdominal wall for the internal and external surfaces. [START_REF] Szymczak | Investigation of abdomen surface deformation due to life excitation: Implications for implant selection and orientation in laparoscopic ventral hernia repair[END_REF] highlight the interest of studying the internal abdominal wall surface since it is where incisional hernias occur. However, as they performed in vivo experiments, they only considered the external abdominal surface. The mean strains computed on the external surface of the abdominal wall were found to be 2.6 times greater than those computed on the internal surface of the abdominal wall. For in vivo study, only the external surface of the abdominal wall can be measured. It is unclear if the mean ratio found between internal and external strains can be applied because of the difference of boundary conditions and of the muscular activity. Also, specimen to specimen variations on this ratio were relatively high (23.3 %). This suggests that relatively large uncertainties would be associated with the estimation of internal strains solely based on external strains.

. Regarding the present study, some limitations can be mentioned. The effect of freezing was not considered however the literature is often contradictory on the effect of freezing of soft tissues [START_REF] Van Ee | Quantifying skeletal muscle properties in cadaveric test specimens: effects of mechanical loading, postmortem time, and freezer storage[END_REF][START_REF] Clavert | Effects of freezing/thawing on the biomechanical properties of human tendons[END_REF][START_REF] Rubod | Biomechanical properties of vaginal tissue. Part 1: new experimental protocol[END_REF]. During harvesting, there was no opportunity to measure the curvature of the abdominal wall and the shortening of tissues. Thus, the initial curvature and tension could not be reproduced in vitro. The initial tension and boundary conditions were not representative of physiological conditions.

Perspectives of this study are to approach the physiological conditions with in vivo experiments. The anisotropy of the strain was assessed in this study but the analysis is limited by the shape of the current device which is not symmetric. However, these findings correspond with the studies analysing the surface deformation of the in vivo human abdominal wall [START_REF] Song | Elasticity of the living abdominal wall in laparoscopic surgery[END_REF][START_REF] Szymczak | Investigation of abdomen surface deformation due to life excitation: Implications for implant selection and orientation in laparoscopic ventral hernia repair[END_REF]. Boundary conditions were well controlled in order to develop a numerical model. Due to the strong boundary conditions, the results of stereo-correlation showed artefact near the edges for the pressure loading case. So, the strains patterns analysis was restricted to the central part of the abdominal wall limited by the nipples.

The main results of this study regarding the stiffness variations between states are currently used to test the validity of numerical models of a porcine abdominal wall. In the future, the same methodology will be used to characterize the passive response of human abdominal walls. Table 1. Mean Lagrange first principal strains E1 in % calculated for the internal and external surfaces of each abdominal wall for the pressure loading case at 50 mmHg.

Table 2. Relative difference of stiffness between the abdominal wall states for the contact loading case.

Table 3. Mean Lagrange first principal strains E1 in % calculated for the external surface of each abdominal wall for the pressure at 50mmHg and for the contact at 165N.

Figure legends : Fig. 1 .

 :1 Figure legends:

Fig. 2 .

 2 Fig. 2. Experimental device: (a) Tightening of the abdominal wall between an aluminium

Fig. 3 .

 3 Fig. 3. Loading mechanisms: (a) Pressure loading, (b) Contact loading

Fig. 5 .

 5 Fig. 5. Lagrange first principal strain fields of the internal (a) and external (b) surfaces of an

Fig. 6 .Fig. 7 .

 67 Fig. 6. Displacement-Force curves for the contact loading for the three states. Error bars:

Fig

  Fig. 1.

Table 2 Δ

 2 Intact-Incised (%) Δ Intact-Repaired (%) Δ Incised-Repaired (%)

	Table 1
		Thickness	E1 (%)	E1 (%)	Ratio
		(mm)	Internal surface	External surface	E1 Ext. / E1 Int.
	Wall	28	4.4	13.7	3.1
	Wall	25	4.5	17.3	3.8
	Wall	29	5.6	13.4	2.4
	Wall	27	5.0	10.8	2.1
	Wall	30	6.6	14.9	2.3
	Wall	30	5.5	12.0	2.2
	Mean	28	5.3	13.7	2.6
	SD	1.7	0.7	2.1	0.6

Acknowledgements

This study was partly funded by Covidien Company.

Conflict of interest statement

This study was partly funded by Covidien Company. G. Guérin and F. Turquier are employees of Covidien France.


Gou-Jen Wang
email: :giwang@dragon.nchu.edu.tw

Chi-Sheng Peng

CONTROL POROUS PATTERN OF ANODIC ALUMINUM OXIDE BY FOILS LAMINATE APPROACH

In this research, a novel, much simpler, and low cost method to fabricate porous pattern of the anodic aluminum oxide (AAO) based on the aluminum foils laminate approach were carried out. During our experiments, it was found that the pores of the AAO grew only on the upper foil, bi-directionally from both the top and the bottom surfaces. Experimental results further indicate that the upward porous pattern of the upper foil is determined by the surface structure of the bottom surface of the upper foil. The porous pattern of AAO can be controlled by a pre-made pattern on the bottom surface. Furthermore, no barrier removing process is required in this novel laminate method. the bottom of the being anodized aluminum sheet acts as a barrier layer during anodization, no barrier removing process is required in this novel laminate approach. The developed work in this study can be further applied to the fabrication of nanofunction devices.

INTRODUCTION

Ordered porous membrane is a desirable material for wide applications such as sensors, optical devices, catalysts, and microfabricated fluidic devices. Anodic aluminum oxide (AAO) membrane [1][2][3] is one of the most attractive ordered nanopore materials and has been demonstrated to be a relatively inexpensive and high throughput template material in fabrication nanodevices such as patterned carbon nano tube (CNT) arrays, photonic crystals, quantum-dot arrays, and nanosensors [4][5][6][7][8].

Key technology that leads to the successful applications of AAO aforementioned is the capability of growing desired porous pattern on the AAO membranes. Li et al. [9] originally integrated the two-anodization process with the lithographic technique to fabricate hexagonally ordered two-dimensional nanopore arrays in anodic alumina. The pore distance can be arranged by changing the anodic electrolyte and the applied voltage. Yan et al. [10] used silica as the anodization barrier to grow highly ordered, patterned through-hole nanopore arrays. Aligned and patterned CNTs based on AAO nano-template with feasible applications in field emitters have recently been reported [7,[11][12][13]. Masuda et al. [14][15] fabricated ordered hole configuration in anodic porous alumina by pretexturing the aluminum. For these patterned nanopore arrays based on AAO template, additionally lithographic technique was always required. Although the incorporation of the lithographic technique into anodic porous alumina can be easily conducted, supplementary equipments and procedures are needed. In this article, we propose a novel, much more convenient and cheaper method for the fabrication of porous pattern of the AAO by foils laminate approach. In this new foils laminate method, two tightly clamped aluminum foils, rather than the conventional single piece foil, are employed as the raw material for anodization. In general, an aluminum barrier layer exists at the bottom of the ordered AAO membrane after anodization. This barrier layer needs further removal to achieve a freestanding and through-pore AAO membrane. If an additionally aluminum sheet is attached to the bottom of the being anodized aluminum sheet and acts as a barrier layer during anodization, no additional barrier removing process is required. Furthermore, the porous pattern may be controlled by the structure of the two contacting surfaces.

CONTROL POROUS PATTERNS OF AN ANODIC ALUMINUM OXIDE

As shown in Figure 1, porous distribution of an anodic aluminum oxide depends on the surface roughness of the aluminum foil. An unpolished aluminum foil results in a relatively disordered porous pattern (Figure 1a). On the contrary, a polished aluminum foil can produce more orderly and complete porous pattern (Figure 1b). Form the atomic force microscopy (AFM) images of the aluminum foil surfaces as illustrated in Figure 2, it is found that the polished aluminum foil possesses flatter surface (Figure 2b), while the unpolished aluminum foil has parallel notches on its surface (Figure 2a). In contrast to Figure 1, it can be inferred that surface structure of an aluminum foil determines the porous pattern of the resulting anodic aluminum oxide. Based on this property, a novel and simple method to control porous patterns of AAO is presented.

Experimental Procedures

The experimental procedures include aluminum foil preparation, electropolishing, aluminum foils clamping, anodization, and aluminum foils separation.

(1) Aluminum foils preparation Aluminum foils (99.9995%, 175 µm thick) were used during the experiments. The aluminum was annealed at 400 °C for 3 h, vibrated by a super sonic vibrator for 1 m, then was cleansed with ethanol to degrease the surfaces.

(2) Electrolytic polishing Electropolished aluminum will produce a bright finish as well as level surface. The aluminum foil is dipped into a bath solution in which the aluminum metal is electrically anodic. This bath consists of about 96% sulfuric acid, 85% phosphoric acid, and DI water with 1:1:1 ratio. The bath runs at 40° C for 14 min, being agitated by a magnet to ensure a high quality polishing.

(3) Aluminum foils clamping The polished aluminum foil is vibrated with a super sonic vibrator for 1 m, and then is cleansed with ethanol to degrease the surfaces. Clamp two aluminum foils tightly together with a Teflon clamper as schematically illustrated in Figure 3.

Al

Experimental results and discussions

Five different experiments were conducted to investigate the dependency of the porous pattern on the surface structure of the laminate foils.

(1) Unpolished aluminum laminate Tightly clamp two unpolished aluminum foils together as shown in Figure 6, followed by anodization. SEM images of the resulting AAO are shown in Figure 7. It can be found that the porous pattern on the top view of the upper aluminum oxide foil (Figure 7a) is similar to the pattern depicted in Figure 1(b), where the pattern resulted from an unpolished aluminum metal processed by conventional approach. The bottom view of the upper aluminum oxide foil is shown in Figure 7(b). The porous pattern spreads along a certain direction and is different from that of the unpolished aluminum sheet under conventional processing. It implies that the porous pattern on the bottom surface of the upper foil is determined by the structure of the bottom surface. The cross-sectional image of the aluminum oxide is shown in Figure 8. It can be seen that the pores grew from both the top and the bottom surfaces. Those pores grew from the top surface are much longer than those from the bottom surface. Ideally, the bottom surface was tightly clamped together with the top surface of the lower metal sheet; therefore, there should not be pore at the bottom surface of the upper metal sheet. We deduce that leakage between the foils may a feasible cause to have the upward pores grow in the notches of the unpolished surface. (2) Lower aluminum foil polished laminate Tightly clamp an unpolished (upper sheet) aluminum foil and a polished (lower sheet) aluminum foil together as shown in Figure 9 and then process anodization. Figure 10 is the bottom view of the upper aluminum oxide foil. The porous pattern is similar to that of the unpolished aluminum laminate shown in Figure 7. It infers that the upward porous pattern of the upper foil has not much strong relationship with the surface structure of the lower metal sheet. (5) Upper aluminum foil pre-patterned laminate Based on the above experimental results, additional experiment with the bottom surface of the upper foil being pre-patterned (Figure 16) were carried out. Figure 17 is the bottom view of the upper aluminum oxide foil. The porous pattern is found to be influenced by the shape of the pre-made pattern on the bottom surface.

CONCLUSIONS

In this research, a series of anodic aluminum oxide fabrication experiments based on the aluminum foils laminate approach were carried out. During the experiments, we found that the pores of the AAO grew only on the upper foil, bi-directionally from both the top and the bottom surfaces. Experimental results further indicate that the upward porous pattern of the upper foil is determined by the surface structure of the bottom surface of the upper foil. The porous pattern of AAO can be controlled by a pre-made pattern on the bottom surface. In addition, since the lower aluminum sheet attached to

Figure 1 .Figure 2 .

 12 Figure 1. SEM image of the AAO made of (a) unpolished and (b) polished aluminum foils (a) (b)

Figure 3 .Figure 4 .

 34 Figure 3. Schematic illustration of the aluminum foils clamping (4) AnodizationAnodization is carried out under conditions of constant voltage 60 V in a 0.3 M oxalic acid solution at 0 °C for 7 h and being stirred by a magnet. After anodization (Figure4), the sample is rinsed again with DI water, and then is dried with ethanol. A milky thin film can be found covering the top surface of the anodized sample and could be removed by DI water.

Figure 5 .

 5 Figure 5. Aluminum foils separation

Figure 6 .Figure 8 .

 68 Figure 6. Schematic diagram of the unpolished aluminum laminate

Figure 9 .Figure 10 .

 910 Figure 9. Schematic diagram of the lower aluminum foil polished laminate

Figure 11 .Figure 12 .Figure 13 .Figure 14 .Figure 15 .

 1112131415 Figure 11. Schematic diagram of the upper aluminum foil polished laminate

Figure 16 .

 16 Figure 16. Schematic diagram of the Upper aluminum foil pre-patterned laminate


H S Ashour
email: hashour@alazhar-gaza.edu

A I Ass'ad

M M Shabat
email: shabat@mail.iugaza.edu.ps

THE DEFECT EFFECT ON THE ELECTRONIC CONDUCTANCE IN BINOMIALLY TAILORED QUANTUM WIRE

The paper considers the effect of the defects on the electronic transmission properties in binomially tailored waveguide quantum wires, in which each Dirac delta function potential strength have been weight on the binomial distribution law. We have assumed that a single free-electron channel is incident on the structure and the scattering of electrons is solely from the geometric nature of the problem. We have used the transfer matrix method to study the electron transmission. We found this novel structure has a good defect tolerance. We found the structure tolerate up to 20% ± in strength defect and 5% ± in position defect for the central Dirac delta function in the binomial distribution. Also, we found this structure can tolerate both defect up to 20% ± in strength and 5% ± in position dislocation.

INTRODUCTION

In the last decade, there was a growing interest in the electron conductance through one-dimensional scattering problems, especially in those cases where the potential is periodic structure with finite number of identical cells [1,2]. Because of the remarkable advances in nanotechnology and micro fabrication, it is possible to confine electrons in a conductor with a lateral extent of 100nm or less, resulting narrow quantum wire [3]. In these mesoscopic devices, the electron transport is best described by quantum mechanics. Miniature size of these devices eliminates the defect of scattering. At a low enough temperature, the motion of electrons through these devices is ballistic or quasiballistic and the electronphonon interaction can be neglected. So that, the phase coherence length enlarges enough when compared with the device dimension. Mesoscopic devices can be considered as a coherent elastic scatterer [3]. Therefore, the electron transport properties, solely depends upon the geometrical structure of the quantum waveguide. In recent years, there has been a growing interest in the electron transport through a sequence of Dirac delta function potential [4][5][START_REF] Fayad | Electronic Conductance through some quantum wire structure[END_REF][START_REF] Griffith | [END_REF][8][9]. The researcher used different methods to study the electron transport in a waveguide quantum wire [10][11][12]. Recently, Ashour et al [13] has proposed a novel structure which is the binomially tailored waveguide quantum wires, in which each Dirac Delta function potential strength has been weight on the binomial distribution law. In this paper, we study the defects effect on the electronic conductance on the novel structure proposed by [13].

In section 2, we outline the transfer matrix method which connects the solutions at the ends of the waveguide quantum wires. We introduce the novel structure of the binomially tailored waveguide quantum wires. In section 3, we explore defect effects on the electronic conductance spectrum through the binomially tailored quantum wire. In this section, we have studied strength defect and dislocation defects on the central Dirac delta function in the distribution. Section 4 has been devoted to the conclusions of this study.

TRANSMISSION MATRIX THOUGH

Periodic Structure

In this context, we consider a finite periodic structure of Dirac delta function potential (Dirac Comb). Also, we assumed that the structure is narrow enough so that just single channel of electrons can be considered. In this treatment, we neglected electron-electron interaction, and we assume the temperature low enough so that electronphonon interaction can be neglected as well. We assumed the scattering of electrons mainly form the geometrical structure of the potential. The potential can be written as follows:

( ) ( )

1 N j j j V x U x x δ = = - ∑ (1)
Thus, j U and j x represent the strength and the position of the j th delta function respectively, and N is the number of the Dirac delta functions in Dirac Comb. The distance between the adjacent barriers are given by

1 j j j d x x + = -.
The Schrödinger wave equation of one dimension can be written as follows:

( ) ( ) ( ) ( ) 2 2 * 2 2 d x V x x E x m dx ψ ψ ψ - + = h (2)
Thus, ( )

V x is the periodic potential given by equation ( 1), * m is the electron effective mass, which is considered approximately constant over the interaction range. The solution of Schrödinger wave equation for single Delta function potential can be found in the literature and also in the transfer matrix formulism [14][15][16][START_REF] Merzbacher | Quantum Mechanics[END_REF]. The transfer matrix for periodic structure has been used also to study the transmission of electron through Comb structure [4-6, 14-16, 17]. The transfer matrix, which is related to the input electron wave and the output, is given by [4-6]

1 1 0 0 1 0 0 j j j j ikx ikx j j ikx ikx j j i i C e e A i i D B e e β β β β - - - ⎛ ⎞ ⎛ ⎞ - - ⎛ ⎞ ⎛ ⎞ ⎛ ⎞ =⎜ ⎟ ⎜ ⎟ ⎜ ⎟ ⎜ ⎟ ⎜ ⎟ ⎜ ⎟ ⎜ ⎟ + + ⎝ ⎠ ⎝ ⎠ ⎝ ⎠ ⎝ ⎠ ⎝ ⎠ (3)
Thus, j andk is the wave number given by * 2 2 / m E h . So that the transfer matrix at a given single barrier can be written as, ( ) ( ) ( )

β is / 2 j k γ , where j γ is * 2 2 / j m U h ,
1 j j j j M S ikx S ikx β - = ⋅Γ ⋅ (4)
The total transfer matrix which represents the electron propagating through the entire device is just the repetitive product of the transfer matrix of a single barrier. We find

( ) ( ) ( ) ( ) ( ) 1 1 1 0 t N N j j j N M S ikd S ikd S β β - = - ⎡ ⎤ = ⋅Γ ⋅ ⋅Γ ⋅ ⎢ ⎥ ⎣ ⎦ ∏ (5)
Thus, j d is the periodic spacing between two adjacent Dirac Delta functions.

Then the transmission amplitude is given by [START_REF] Sprung | [END_REF],

( ) . A perfect transmission in this case is in general impossible as predicted by [4,22]. According to reference [START_REF] Takagaki | [END_REF]

The Binomially Tailored Quantum Wire (BTQW)

In this subsection, we reintroduce BTQW structure as shown in figure 2. The Dirac delta function has been equally spaced but their strength, j Ω , has been weighted according to the binomial distribution law, which is ( ) Thus, ( )

/ 2 , 0,1,...., N j j j N N N N N ⎛ ⎞ Ω = = ⎜ ⎟ ⎝ ⎠ (8) ( ) j N Ω j N
j N Ω
represents the strength of the Dirac delta potential, 1 N + represents the total number of Dirac delta function potentials in the quantum wire, and j N is the order of the Dirac delta potential. This novel structure of quantum wires can be released by putting metallic gates on top of a one dimensional electron gas and then by applying voltages, according to the binomial distribution law, to deplete the electron gas underneath. In this case, equation ( 9) is no longer valid for our new structure. So that the total transmission matrix can be written as follows:

( ) ( ) ( ) ( ) ( ) ( ) ( ) 0 . .......... 1 1 1 1 1 S ikd S ikd S ikd S M N N N N t ⋅ Γ ⋅ × × Γ ⋅ ⋅ Γ ⋅ = + + - β β β (9)
Notice that the potential strength is weighted according to equation (8). In figure 3, we show the conductance spectrum through a sequence of a binomially tailored Dirac delta function potentials. It is quite interesting to notice that we have reached a transmission through this structure approaches to unity in the allowed band region without any ripples after some k value. Here, we have a resonant tunneling due to coherent interference effects due to elastic scattering of electrons, which leads the transmission to reach unity and also to have constant value over the allowed band or conduction band. Also, we see that there is a forbidden band or conduction gap where the transmission is small. G Figure 3: Conductance spectrum G in the units of 2 2 / e has a function of / kd π for a BTQW, with N=5.

/

DEFECTS EFFECT ON THE ELECTRONIC CONDUCTANCE OF BTQW

Strength Defect

In this subsection, we study effect of strength defect of the central element of the binomial tailored quantum wire, and keeping the other elements and the spacing between the Dirac delta function potentials constant, on the electronic conductance through the BTQW. First, we consider the strength defect does not exceed 5% ± of the Dirac delta function potential strength. That is, when the central Dirac delta functions potentials strength is ( ) ( )

/ 2 1 0.05 / 2 1 j j N N Ω + ± Ω + .
In figure 4-a, we plot the electronic conductance spectrum for both strengths with 35 j N is and scaling factor of three. As can noticed there is slight difference between the two conductance spectrum curves, and defect free curves. In figure 4-b, we have increased the strength defect up to 20% ± , we have noticed some measurable differences between the two the conductance spectrum curves and defect free curves, but still the conduction band and the forbidden bands well defined, which is a very good feature for BTQW.

Dislocation Effect

In this subsection, we study dislocation defect effect on the position of the central element in the BTQW, and keeping the other elements and the spacing between the Dirac delta function potentials constant. First, we consider the position defect does not exceed 5% ± of the Dirac delta function potential spacing constant. That is, when the central Dirac delta function potentials spacing is 0.05 d d ±

. In figure 5-a, we plot the electronic conductance spectrum for both dislocations with 35 j N is and scaling factor of three. Compared to defect curves, as can noticed there is a difference between the two curves. The conduction band starts lose its flatness and the forbidden band shaper for increased spacing between the central Dirac delta function and the adjacent one. In figure 5-b, we increase the dislocation defect up to 20% ± , we have noticed measurable differences between the two curves and that of no defect case, but still the conduction band is well defined but the forbidden bands have a split compared to forbidden band in no defect curves. This splitting is due to resonant state in the forbidden energy band which leads to a bound state in the structure [3]. This is because the particle mode cannot propagate and hence get trapped.

CONCLUSION

We found the novel structure introduced by [13] has a good tolerance for strength reaches up to 20% ± and dislocation tolerance reaches up to 5% ± without losing the fascinating electronic transmission characteristics.

Figure 2 :

 2 Figure 2: Shows binomially tailored Dirac delta function potential.

jN

  values weighted by equation(8).

Figure 4

 4 Figure 4-b: The electronic conductance, in the units of 2 2 / e has

/

  Figure5-a: The electronic in the units of2 2 / e has

Figure 5 -/

 5 Figure 5-b: The electronic conductance, in the units of 2 2 / e has


Laurent Fribourg
email: fribourg@lsv.ens-cachan.fr

Euler's method applied to the control of switched systems

published or not. The documents may come

Euler's method applied to the control of switched systems

Laurent Fribourg

Introduction

In this paper, we present some recent results obtained for the control synthesis of nonlinear switched systems using the one-sided Lipschitz conditions of their dynamics. The main idea is to use "one-sided Lipschitz conditions" on the system vector fields in order to generate a sequence of balls enclosing the sets of trajectories. The method can be easily extended to take into account uncertainty and compositional synthesis. These results mainly originate from collaboration with A. Le Coënt, F. De Vuyst, L. Chamoin, J. Alexandre dit Sandretto and A. Chapoutot (see [START_REF] Le Coënt | Control synthesis of nonlinear sampled switched systems using Euler's method[END_REF][START_REF] Le Coënt | Distributed control synthesis using Euler's method[END_REF]).

The plan of the paper is as follows: in Section 2, we present the notions of switched systems and (R, S)-stability; in Section 3, we introduce a new error analysis for Euler's method, and explain how to use it for ensuring (R, S)stability in control synthesis of switched systems; we extend this control syn-thesis method to uncertain switched systems, and to compositional synthesis (Section 4); we conclude in Section 5.

2 Switched systems and (R,S)-stability

Switched systems

A hybrid system is a system where the state evolves continuously according to several possible modes, and where the change of modes (switching) is done instantaneously. We consider here the special case of hybrid systems called "sampled switched systems" where the change of modes occurs periodically with a period of τ seconds. We will suppose furthermore that the state keeps its value when the mode is changed (no jump). More formally, we denote the state of the system at time t by x(t) ∈ R n . The set of modes U = {1, . . . , N } is finite. With each mode j ∈ U is associated a vector field f j that governs the state x(t); we have: ẋ(t) = f j (x(t))

We make the following hypothesis:

(H0) For all j ∈ U, f j is a locally Lipschitz continuous map.

We will denote by φ j (t; x 0 ) the solution at time t of the system ẋ(t) = f j (x(t)),

x(0) = x 0 . (1)
The existence of φ j is guaranteed by assumption (H0). Let us consider S ⊂ R n be a compact and convex set, typically a "box" or "rectangular set", that is a cartesian product of n closed intervals. We know by (H0) that there exists a constant L j > 0 such that:

f j (y) -f j (x) ≤ L j y -x ∀x, y ∈ S. (2)
We also define, for all j ∈ U :

C j = sup x∈S L j f j (x) . (3)
Example 1. One consider the example (adapted from [START_REF] Girard | Low-complexity switching controllers for safety using symbolic models[END_REF]) of a two rooms apartment, with one heater per room. See Figure 1. There is heat exchange between the two rooms and with the environment. The objective is to control the temperature of the two rooms. The continuous dynamics of the system is given by the equation: A pattern π is a finite sequence of modes; e.g., the expression

˙ T1 T2 = -α21 -αe1 -α f j1 α21 α12 -α12 -αe2 -α f j2 T1 T2 + αe1Te + α f T f j1 αe2Te + α f T f j2 .
0 1 • 0 0 • 1 1 is a pattern in Example 1.
The (state-dependent) control synthesis problem consists in finding at each sampling time τ , 2τ , . . . , the appropriate mode u ∈ U (in function of the current value of x) to be selected for satisfying some objective, for example a safety property. More generally, the control synthesis problem (with a "time-horizon" bounded by a positive integer K) consists first in selecting at time 0 a pattern π1 of length, say 1 ≤ k1 ≤ K, according to the value of state x(0); then after k1τ seconds, selecting a new pattern π2, according to the value of x(k1τ ), and so on repeatedly. This induces a control (or switching) rule σ which is a piecewise constant function of time, with discontinuities occurring at sampling times. By convention, the control law σ is right-continuous.

(R, S)-stability

Among the classical objectives that one is generally aiming for, there are the reachability objective: given an initial region Rinit and a target region R, find a pattern which drives x(t) to R, for any initial state x 0 = x(0) ∈ Rinit; the stability objective: for any initial point

x 0 = x(0) ∈ R, find a pattern π ∈ U k (with 1 ≤ k ≤ K) which makes the trajectory return in R (i.e.: x(kτ ) ∈ R) while always maintaining x(t) in a neighborhood S = R + ε of R, (i.e.: x(t) ∈ S for 0 ≤ t ≤ kτ ).
The effect of such control rules is depicted on Figure 2.

For the sake of simplicity, we focus here on a property that we call "(R, S)-stability": given two rectangular sets (i.e., cartesian products of intervals) R and S with R ⊆ S ⊂ Fig. 2. Illustration of reachability (left) followed by stability (right) R n , called respectively "recurrence set" and "safety set", the (R, S)-stability control problem consists in finding a control σ ensuring, for all x(0) ∈ R 1. recurrence: the state of the system x(t) belongs to R for an infinite number of values of t; 2. safety: the state of the system x(t) belongs to S for all t ≥ 0.

The property of (R, S)-stability is illustrated in Figure 3 We now give the general scheme of control synthesis that has been proposed in MINIMATOR [START_REF] Fribourg | Finite controlled invariants for sampled switched systems[END_REF] for ensuring (R, S)-stability. This scheme consists in two steps:

1. cover R via a finite number m of subsets B 0 1 , B 0 2 , ..., B 0 m of S (with R ⊂ m i=1 B 0 i ⊆ S); 2. for each B 0 i (1 ≤ i ≤ m
), find a pattern πi of length ki ≤ K such that, starting at t = 0 from any point of B 0 i , the trajectory x(t) controlled by πi satisfies:

x(t) ∈ S for all t ∈ [0, kiτ ] ∧ x(t) ∈ R for t = kiτ.
Note that, when the system returns to R (after application of some pattern) at time, say t = t1, the state x(t1) belongs to B 0 i 1 for some 1 ≤ i1 ≤ m; the pattern πi 1 is then applied, which makes the system return to R at time t2 = t1 + ki 1 τ , and so on iteratively.

Remark 1. Let us give a rough estimation of the complexity of MINIMATOR scheme. Let N be the number of modes, n the state dimension, K the time-horizon (or maximum length of patterns), m = 2 nd the number of modes (assuming a uniform covering obtained by bisection of depth d); the MINIMATOR scheme consists essentially in enumerating all the possible patterns of length ≤ K until finding, for each B 0 i (1 ≤ i ≤ m) a safe recurrent candidate; a simple calculation shows that there are 2 nd N K candidate patterns; the complexity of the MINIMATOR scheme is thus exponential in n, d, K (note that the number of modes N may be itself exponential in the dimension n: for example, in a classical n-room heating example with one heater per room and two modes by heater, there are N = 2 n modes).

Remark 2. Note that the set of trajectories starting at points of R form a (positive) invariant set included into S. There are classical methods for generating (maximal) invariant sets included into S ( [START_REF] Asarin | Effective synthesis of switching controllers for linear systems[END_REF][START_REF] Blanchini | Set invariance in control: a survey[END_REF]). Unfortunately, these general methods are based on a backward reachability constructs, which, as noticed by I.M. Mitchell [START_REF] Mitchell | Comparing forward and backward reachability as tools for safety analysis[END_REF], "are more likely to suffer from numerical stability issues, especially in systems with significant contraction -the very systems where forward simulation and reachability are most effective". The forward analysis used by the MINIMATOR scheme (application of patterns) avoids such a difficulty.

Guaranteed integration

The MINIMATOR paradigm described in Section 2.2 relies implicitly on the existence of a process for overapproximating the set of trajectories originating from a subset B 0 i during a multiple of sampling periods. Such a process is called "guaranteed integration" (or "set-integration"). As said in [START_REF] Nedialko | A New Perspective on the Wrapping Effect in Interval Methods for Initial Value Problems for Ordinary Differential Equations[END_REF]: "Methods of guaranteed integration are methods capable to compute bounds that are guaranteed to contain the solution of a given ODE at points tj, j = 1, 2, . . . , m in the interval (t0, tm] for some tm > t0. These methods are usually based on Taylor series or extension of Hermite-Obreschkoff schemes to interval methods. They usually consist of two phases. On an integration step from tj-1 to tj, the first phase validates existence and uniqueness of the solution of (1) for all [tj-1, tj] and computes a priori bounds for this solution for al t ∈ [tj-1, tj], [19,20]; and the second phase compute tight bounds for the solution of (1) at tj. Note that a major problem in the second phase is the wrapping effect [16]. It occurs when a solution set that is not a box in R n , n ≥ 2, is enclosed, or wrapped, by a box on each integration step. (...) As a result of such a wrapping, an overestimation is often introduced on each integration step. Those overestimations accumulate as the integration proceeds, and the computed bounds may soon become unacceptably large. Many methods have been proposed to reduce the wrapping effect in the context of interval methods."

In order to avoid such a wrapping effect, we proposed an alternate method which, instead of using interval arithmetic [START_REF] Moore | Interval Analysis[END_REF] and higher order Taylor series, has simply recourse to the basic (forward) Euler method [START_REF] Le Coënt | Control synthesis of nonlinear sampled switched systems using Euler's method[END_REF]. This is made possible through a new error analysis of the Euler method via the notion of "one-sided Lipschitz constant".

3 Euler's method and error estimation

One-sided Lipschitz constant

As remarked in [START_REF] Abbaszadeh | Nonlinear observer design for onesided lipschitz systems[END_REF]:

"The Lipschitz constant of [many] functions is usually region-based and often dramatically increases as the operating region is enlarged. On the other hand, even if the nonlinear system is Lipschitz in the region of interest, it is generally the case that the available observer design techniques can only stabilize the error dynamics for dynamical systems with small Lipschitz constants but fails to provide a solution when the Lipschitz constant becomes large. The problem becomes worse when dealing with stiff systems. Stiffness means that the ordinary differential equation (ODE) admits a smooth solution with moderate derivatives, together with nonsmooth ("transient") solutions rapidly converging towards the smooth ones (...) This problem has been recognized in the mathematical literature and specially in the field of numerical analysis for some time and a powerful tool has developed to overcome this problem. This tool is a generalization of the Lipschitz continuity to a less restrictive condition known as one-sided Lipschitz (OSL) continuity."

Unlike Lipschitz constants, OSL constants can be negative, which express a form of contractivity of the system dynamics. Even if the OSL constant is positive, it is in practice much lower than the Lipschitz constant [START_REF] Dahlquist | Error analysis for a class of methods for stiff non-linear initial value problems[END_REF]. The use of OSL thus allows us to obtain an upper bound for the error associated with Euler's method that is more precise than by using Lipschitz constants [START_REF] Le Coënt | Control synthesis of nonlinear sampled switched systems using Euler's method[END_REF].

Let us denote by T a compact overapproximation of the image by φj of box S for 0 ≤ t ≤ τ and j ∈ U , i.e. T is such that

T ⊇ {φj(t; x 0 ) | j ∈ U, 0 ≤ t ≤ τ, x 0 ∈ S}.
The existence of T is guaranteed by assumption (H0). We now make the additional hypothesis that the vector fields fj of the system are one-sided Lipschitz (OSL) [START_REF] Donchev | Stability and Euler approximation of onesided lipschitz differential inclusions[END_REF]. Formally:

(HU ) For all j ∈ U , there exists a constant λj ∈ R such that fj(y) -fj(x), y -x ≤ λj y -x 2 ∀x, y ∈ T, where •, • denotes the scalar product of two vectors of R n . Remark 3. Constants λj as well as Lj and Cj (j ∈ U ) can be computed using constrained optimization algorithms. See Section 3.5 for details.

Euler approximate solutions

Given an initial point x0 ∈ S and a mode j ∈ U , we define the following "linear approximate solution" φj(t; x0 ) for t ∈ [0, τ ] by:

φj(t; x0 ) = tfj(x 0 ) + x0 . (4)
Formula ( 4) is nothing else but the explicit forward Euler scheme with"time step" t. It is thus a consistent approximation of order 1 in t of the exact solution of (2.1) under the hypothesis x0 = x 0 . More generally, given an initial point x0 ∈ S and pattern π of U k , we can define a "(piecewise linear) approximate solution" φπ(t; x0 ) of φπ at time t ∈ [0; kτ ] as follows:

-φπ(t; x0 ) = tfj(x 0 ) + x0 if π = j ∈ U , k = 1 and t ∈ [0, τ ],
and

-φπ(kτ + t; x0 ) = tfj(z) + z with z = φπ ((k -1)τ ; x0 ), if k ≥ 2, t ∈ [0, τ ], π = j • π for some j ∈ U and π ∈ U k-1 .
Fig. 4. Illustration of Euler's method (from Wikipedia)

We wish to synthesize a guaranteed control σ using approximate functions of the form φπ. We define the closed ball of center x ∈ R n and radius r > 0, denoted B(x, r), as the set {x ∈ R n | x -x ≤ r}. Given a positive real δ 0 , we now define the expression δj(t) which, as we will see in Theorem 1, represents (an upper bound on) the error associated to φj(t; x0 ) (i.e. φj(t; x0 ) -φj(t; x 0 ) ).

Definition 1. Let δ 0 be a positive constant. Let us define, for all 0 ≤ t ≤ τ , δj(t) as follows:

if λj < 0:

δj(t) = (δ 0 ) 2 e λ j t + C 2 j λ 2 j t 2 + 2t λj + 2 λ 2 j 1 -e λ j t 1 2
if λj = 0 :

δj(t) = (δ 0 ) 2 e t + C 2 j (-t 2 -2t + 2(e t -1)) 1 2
if λj > 0 :

δj(t) = (δ 0 ) 2 e 3λ j t + C 2 j 3λ 2 j -t 2 - 2t 3λj + 2 9λ 2 j e 3λ j t -1 1 2
Note that δj(t) = δ 0 for t = 0. The function δj(•) depends implicitly on parameter: δ 0 ∈ R>0. In Section 3.3, we will use the notation δ j (•) where the value of δ j (t) for t = 0 is implicitly a parameter denoted by (δ ) 0 . Theorem 1. Given an ODE system satisfying (H0 -HU ), consider a point x0 and a positive real δ 0 . We have, for all x 0 ∈ B(x 0 , δ 0 ), t ∈ [0, τ ]:

φj(t; x 0 ) ∈ B( φj(t, x0 ), δj(t)).
The proof of this theorem is given in [START_REF] Le Coënt | Control synthesis of nonlinear sampled switched systems using Euler's method[END_REF].

Remark 4. In Theorem 1, we have supposed that the step size h used in Euler's method was equal to the sampling period τ of the switching system. Actually, in order to have better approximations, it is often convenient to take a fraction of τ as for h (e.g., h = τ 10 ). Such a splitting is called "sub-sampling" in numerical methods.

Corollary 1. (one-step invariance) Given an ODE system satisfying (H0 -HU ), consider a point x0 ∈ S and a real δ 0 > 0 such that:

1. B(x 0 , δ 0 ) ⊆ S, 2. B( φj(τ ; x0 ), δj(τ )) ⊆ S, and 3.

d 2 (δ j (t)) dt 2 > 0 for all t ∈ [0, τ ].
Then we have, for all x 0 ∈ B(x 0 , δ 0 ) and t ∈ [0, τ ]: φj(t; x 0 ) ∈ S. ) with j = π(k ).

Let us now define the expression δ k π as follows: for k = 0: δ k π = δ 0 , and for 1 ≤ k ≤ k: δ k π = δ j (τ ) where (δ ) 0 denotes δ k -1 π , and j denotes π(k ). Likewise, the expression δπ(t) is defined, for 0 ≤ t ≤ kτ , by:

-for t = 0: δπ(t) = δ 0 , -for 0 < t ≤ kτ : δπ(t) = δ j (t ) with (δ ) 0 = δ -1 π , j = π( ), t = t -( -1)τ and = t τ . Note that, for 0 ≤ k ≤ k, we have: δπ(k τ ) = δ k π .
Following the MINIMATOR paradigm (see Section 2.2), we are now ready to synthesize a control σ ensuring (R, S)stability, using the approximate functions φπ.

Theorem 2. Given a sampled switched system satisfying (H0-HU ), consider a point x0 ∈ S, a positive real δ 0 and a pattern π of length k such that, for all 1 ≤ k ≤ k:

1. B(x k π , δ k π ) ⊆ S and 2. d 2 (δ j (t)) dt 2 > 0 for all t ∈ [0, τ ], with j = π(k ) and (δ ) 0 = δ k -1 π .
Then we have, for all x 0 ∈ B(x 0 , δ 0 ) and t ∈ [0, kτ ]: φπ(t; x 0 ) ∈ S.

Corollary 2. Given a switched system satisfying (H0 -HU ), consider a positive real δ 0 and a finite set of points x1, . . . xm of S such that all the balls B(xi, δ 0 ) cover R and are included into S (i.e. R ⊆ m i=1 B(xi, δ 0 ) ⊆ S). Suppose furthermore that, for all 1 ≤ i ≤ m, there exists a pattern πi of length ki such that:

1. B((xi) k π i , δ k π i ) ⊆ S, for all k = 1, . . . , ki -1 2. B((xi) k i π i , δ k i π i ) ⊆ R. 3. d 2 (δ j (t)) dt 2 > 0 with j = πi(k ) and (δ ) 0 = δ k -1 π i
, for all k ∈ {1, ..., ki} and t ∈ [0, τ ].

These properties induce a control σ 1 which guarantees -(safety): if x 0 ∈ R, then φσ(t; x 0 ) ∈ S for all t ≥ 0, and -(recurrence): if x 0 ∈ R then φσ(kτ ; x 0 ) ∈ R for some k ∈ {k1, . . . , km}.

A covering of R with balls as stated in Corollary 2 is illustrated in Figure 6 (left) with a pattern satisfying safety and recurrence in Figure 6 (right). Corollary 2 thus leads to the following method (inspired by the the MINIMATOR scheme described in Section 2.2), aiming for (R, S)-stability:

we (pre-)compute λj, Lj, Cj for all j ∈ U ; we find m points x1, . . . , xm of S and δ 0 > 0 such that R ⊆ m i=1 B(xi, δ 0 )) ⊆ S; we find m patterns πi (i = 1, . . . , m) such that conditions 1-2-3 of Corollary 2 are satisfied.

R S x1 x2 x3 x4 x5 x6 x7 x8 x9 R S δ B( x 3 1 , δ π3 1 ) x3 B( x 3 2 , δ π3 2 )
Fig. 6. Set of balls covering R (left) and safe recurrent pattern associated with one of these balls (right).

Avoiding wrapping effect with Euler's method

The problem of "wrapping effect" inherent to the method of interval analysis has been noticed from the outset: R. Moore [START_REF] Moore | Interval Analysis[END_REF] illustrates it on the simple rotation ẋ = 0 1 -1 0 x; x0 ∈ A for an initial set A which is rectangular. At each step, the rectangle is rotated and has to be wrapped by another one. At t = 2π, the blow up factor is by a factor e 2π ≈ 535, as the step size tends to zero (Figure 7: left). In contrast, the application of the Eulerbased method starting from a ball of radius 0.1 with step size 0.005, does not blow up on this example (Figure 7: right).

Numerical results

Our Euler-based synthesis method has been implemented by Adrien Le Coënt in the interpreted language Octave, and the experiments performed on a 2.80 GHz Intel Core i7-4810MQ CPU with 8 GB of memory. The computation of constants Lj, Cj, λj (j ∈ U ) are realized with a constrained optimization algorithm. They are performed using the "sqp" function of Octave, applied on the following optimization problems:

-Constant Lj:

Lj = max x,y∈S, x =y fj(y) -fj(x) y -x -Constant Cj: Cj = max x∈S Lj fj(x)
-Constant λj: λj = max

x,y∈T, x =y fj(y) -fj(x), y -x y -x 2 Fig. 7. left: guaranteed integration with interval method (from [START_REF] Moore | Interval Analysis[END_REF]); right: with Eulerbased method.

The convexity test

d 2 (δ j (t)) dt 2
> 0 can be performed similarly. Note that in some cases, it is advantageous to use a time sub-sampling to compute the image of a ball. Indeed, because of the exponential growth of the radius δj(t) within time, computing a sequence of balls can lead to smaller ball images. It is particularly advantageous when a constant λj is negative. We illustrate this with the example of the DC-DC converter [START_REF] Beccuti | Optimal control of the boost dc-dc converter[END_REF]. It has two switched modes, for which we have λ1 ≈ -0.014 and λ2 ≈ 0.14. In the case λj < 0, the associated formula δj(t) has the behavior of Figure 8 (a). In the case λj > 0, the associated formula δj(t) has the behavior of Figure 8 (b). In the case λj < 0, if the time sub-sampling is small enough, one can compute a sequence of balls with reducing radius, which makes the synthesis easier. Example 2. (Four-room apartment) We describe a first application on a 4-room 16switch building ventilation case study adapted from [START_REF] Meyer | Invariance and symbolic control of cooperative systems for temperature regulation in intelligent buildings[END_REF]. The model has been simplified in order to get constant parameters. The system is a four room apartment subject to heat transfer between the rooms, with the external environment, the underfloor, and human beings. The dynamics of the system is given by the following equation:

dTi dt = j∈N * \{i} aij(Tj-Ti)+δs i bi(T 4 s i -T 4 i )+ci max 0, Vi -V * i Vi -V * i (Tu-Ti), for i = 1, ..., 4.
The state of the system is given by the temperatures in the rooms Ti, for i ∈ N = {1, . . . , 4}. Room i is subject to heat exchange with different entities stated by the indices N * = {1, 2, 3, 4, u, o, c}. We have T0 = 30, Tc = 30, Tu = 17, δs i = 1 for i ∈ N . The (constant) parameters Ts i , V * i , Vi, aij, bi, ci are given in [START_REF] Meyer | Invariance and symbolic control of cooperative systems for temperature regulation in intelligent buildings[END_REF]. The control input is Vi (i ∈ N ). In the experiment, V1 and V4 can take the values 0V or 3.5V, and V2 and V3 can take the values 0V or 3V. This leads to 16 switching modes corresponding to the different possible combinations of voltages Vi. The sampling period is τ = 30s. Compared simulations are given in Figure 9. On this example, the Euler-based method works better than DynIBEX in terms of CPU time (see Table 1 1. Numerical results for the four-room example.

Bounded uncertainty

Let us now consider the case where the mode j is governed by the uncertain ODE:

ẋ(t) = fj(x(t), w(t)) with w(t) ∈ W
where W is a bounded set of diameter Let us suppose that the uncertain ODE satisfies the assumption:

(HU,W ) For all j ∈ U , there exist λj ∈ R and γj ∈ R ≥0 such that, for all x, x ∈ T , and all w, w ∈ W :

fj(x, w) -fj(x , w ), x -x ≤ λj x -x 2 + γj x -x w -w .
Definition 2. Let δ 0 be a positive real, and W a rectangular set of diameter |W |. We define, for all j ∈ U and 0 ≤ t ≤ τ , the expression δj,W (t) as follows:

-if λj < 0, δj,W (t) = C 2 j -λ 4 j -λ 2 j t 2 -2λjt + 2e λ j t -2 + 1 λ 2 j Cjγj|W | -λj -λjt + e λ j t -1 + λj γ 2 j (|W |/2) 2 -λj (e λ j t -1) + λj(δ 0 ) 2 e λ j t 1/2 (5)
-if λj = 0, δj,W (t) = C 2 j -t 2 -2t + 2e t -2 + Cjγj|W | -t + e t -1 +(γ 2 j (|W |/2) 2 (e t -1)+(δ 0 ) 2 e t ))) 1/2 (6) -if λj > 0, δj,W (t) = 1 (3λj) 3/2 C 2 j λj -9λ 2 j t 2 -6λjt + 2e 3λ j t -2 + 3λj Cjγj|W | λj -3λjt + e 3λ j t -1 + 3λj γ 2 j (|W |/2) 2 λj (e 3λ j t -1) + 3λj(δ 0 ) 2 e 3λ j t 1/2 (7)
Under assumption (HU,W ) instead of (HU ), one can naturally extend Theorem 1 and Corollary 1 to take the uncertainty set W into account, using δj,W (•) in place of δj(•). These extended results are useful to control systems with uncertainty, for example when the coefficients in the vector field definitions are known with a limited precision. Such extended forms of Theorem 1 and Corollary 1 can also be applied to control interconnected subsystems, each component regarding the input from the other one as a form of bounded uncertainty (see Section 4.2).

Application to distributed control synthesis

We now consider the distributed (or "compositional") approach which consists in splitting the original system into two sub-systems, in order to synthesize a controller σi (i = 1, 2) for each sub-system independently, then apply the control σ = (σ1|σ2) (by concurrent application of σ1 and σ2) to the global system. The interest of the approach is to break the exponential complexity of the original method w.r.t. the dimension of the system and the number of modes (see Section 2.2). We consider an ODE of the form ẋ = fj(x) with x ∈ R n , j ∈ U , which is of the form

ẋ1 = f 1 j 1 (x1, x2) (8) ẋ2 = f 2 j 2 (x1, x2) (9)
where the state x is of the form (x1, x2) , and a mode j = (j1, j2) ∈ U = U1 × U2, the solution of the ODE is now denoted by φ (j 1 ,j 2 ) (t; x 0 ), for all t ∈ [0, τ ]. The system (8-9) can be seen as the interconnection of a 1st sub-system [START_REF] Fribourg | Finite controlled invariants for sampled switched systems[END_REF] where x2 plays the role of an "input" given by ( 9), with a 2nd sub-system [START_REF] Girard | Low-complexity switching controllers for safety using symbolic models[END_REF] where x1 is an "input" given by [START_REF] Fribourg | Finite controlled invariants for sampled switched systems[END_REF]. Accordingly, the sets R, S and T are seen under their compositional form R = R1 × R2, S = S1 × S2, T = T1 × T2. We will denote by x m 1 (resp. x m 2 ) an arbitrary point of R1 (resp. R2), typically its central point. We denote by L 1 j 1 the Lipschitz constant for sub-system 1 under mode j1:

with x1 ∈ R n 1 , x2 ∈ R n 2 , n1 + n2 = n,
f 1 j 1 (x1, x2) -f 1 j 1 (y1, y2) ≤ L 1 j 1 x1 x2 - y1 y2
We introduce also the constant:

C 1 j 1 = sup x 1 ∈S 1 L 1 j 1 f 1 j 1 (x1, x m 2 )
Similarly, we define the constants for sub-system 2:

f 2 j 2 (x1, x2) -f 2 j 2 (y1, y2) ≤ L 2 j 2 x1 x2 - y1 y2 and C 2 j 2 = sup x 2 ∈S 2 L 2 j 2 f 2 j 2 (x m 1 , x2)
In the following, we assume that, for all j1 ∈ U1, there exist a real λj 1 and a nonnegative real γj 1 which make the 1st sub-system satisfy assumption (HU 1 ,W 2 ) for some overapproximation W2 of T2. Symmetrically, we assume that, for all j2 ∈ U2, there exist a real λj 2 and a non-negative real γj 2 which make the 2nd sub-system satisfy (HU 2 ,W 1 ) for some overapproximation W1 of T1.

Given two modes j1 ∈ U1, j2 ∈ U2, and two initial conditions x0 1 , x0 2 , we define the "decompositional" Euler approximate solutions φ1 j 1 and φ2 j 2 , for t ∈ [0, τ ], as follows:

φ1 j 1 (t; x0 1 ) = x0 1 + tf 1 j 1 (x 0 1 , x m 2 ) (10) φ2 j 2 (t; x0 2 ) = x0 2 + tf 2 j 2 (x m 1 , x0 2 ) (11)
We can now give the distributed version of Theorem 1.

Theorem 3. Given a distributed sampled switched system satisfying, suppose that the 1st and 2nd sub-systems satisfy, for all j1 ∈ U1 and j2 ∈ U2, the assumptions (HU 1 ,W 2 ) and (HU 2 ,W 1 ) respectively. Consider a point x0 1 and a positive real δ 0 . We have, for all

x 0 1 ∈ B(x 0 1 , δ 0 ), t ∈ [0, τ ], j1 ∈ U1: φ (j 1 ,j 2 ) (t; x 0 ) |1 ∈ B( φ1 j 1 (t, x0 1 ), δj 1 ,W 2 (t)) ∀j2 ∈ U2, ∀x 0 2 ∈ S2, x 0 = x 0 1 x 0 2 .
Likewise, we have, for all

x 0 2 ∈ B(x 0 2 , δ 0 ), t ∈ [0, τ ], j2 ∈ U2: φ (j 1 ,j 2 ) (t; x 0 ) |2 ∈ B( φ2 j 2 (t, x0 2 ), δj 2 ,W 1 (t)) ∀j1 ∈ U1, ∀x 0 1 ∈ S1, x 0 = x 0 1 x 0 2 .
The proof of this theorem is in [START_REF] Le Coënt | Distributed control synthesis using Euler's method[END_REF]. We can now state the distributed version of Corollary 2.

Corollary 3. Given a positive real δ 0 , consider two sets of points x1 1 , . . . , x1 m 1 and x2 1 , . . . , x2 m 2 such that all the balls B(x 1 i 1 , δ 0 ) and B(x 2 i 2 , δ 0 ), for 1 ≤ i1 ≤ m1 and 1 ≤ i2 ≤ m2, cover R1 and R2. Suppose that there exist patterns π 1 i 1 of length ki 1 for the 1st sub-system such that :

1. B((x 1 i 1 ) k π 1 i 1 , δ k π 1 i 1
) ⊆ S1, for all k = 1, . . . , ki 1 -1;

2. B((x 1 i 1 )

k i 1 π 1 i 1 , δ k i 1 π 1 i 1 ) ⊆ R1; 3. d 2 (δ j 1 (t)) dt 2 > 0 with j1 = π 1 i 1 (k ) and (δ ) 0 = δ k -1 π 1 i 1
, for all k ∈ {1, ..., ki 1 } and

t ∈ [0, τ ].
and symmetrically for the 2nd sub-system. These properties induce a control σ1 for the 1st sub-system, and σ2 for the 2nd sub-system such that the composed control σ = (σ1|σ2) ensures recurrence in R and safety in S, i.e.:

if x 0 ∈ R, then φσ(t; x 0 ) ∈ S for all t ≥ 0; if x 0 ∈ R, then φσ(ki 1 τ ; x 0 ) |1 ∈ R1 for some i1 ∈ {1, . . . , m1}, and symmetrically φσ(ki 2 τ ; x 0 ) |2 ∈ R2 for some i2 ∈ {1, . . . , m2}.

Example 3. We demonstrate the interest of the distributed approach by comparing it with respect to the (centralized) approach performed in Example 2. The main difficulty of this example is the large number of modes in the switching system, which induces a combinatorial issue. The centralized controller in Example 2 was obtained with 256 balls in 48 seconds, the distributed controller was obtained with 16+16 balls in less than a second. In both cases, patterns of length 2 are used. A sub-sampling of h = τ /20 is required to obtain a controller with the centralized approach (see Table 2). For the distributed approach, no sub-sampling is required for the first sub-system, while the second one requires a sub-sampling of h = τ /10 (see Table 3). Simulations of the centralized and distributed controllers are given in Figure 10, where the control objective is to stabilize the temperature in [20,22] 3. Numerical results for the distributed four-room example.

Final remarks

We have presented a simple method of control synthesis for switched systems using a new scheme of guaranteed integration based on Euler's method. Preliminary experiments show that, on some examples, the method avoids the wrapping effect occurring with interval-based integration methods. On-going work is done for adapting this Eulerbased method to the treatment of stochastic differential equations.

Fig. 1 .

 1 Fig. 1. 2-rooms example

  in the case of Example 1, with R = [18, 22] × [18, 22].

Fig. 3 .

 3 Fig. 3. (R, S)-stability

Corollary 1 is illustrated in Figure 5 .Fig. 5 .

 55 Fig. 5. Illustration of one-step invariance in S

3. 3 1 π+

 31 Application to control synthesis for (R, S)-stability Consider a point x0 ∈ S, a positive real δ 0 and a pattern π of length k. Let π(k ) denote the k -th element (mode) of π for 1 ≤ k ≤ k. Let us abbreviate the k -th approximate point φπ(k τ ; x0 ) as xk π for k = 1, ..., k, and let xk π = x0 for k = 0. It is easy to show that xk π can be defined recursively for k = 1, ..., k, by: xk π = xkτ fj(x k -1 π

Fig. 8 .

 8 Fig. 8. Behavior of δj(t) for the DC-DC converter with δj(0) = 0.045. (a) Evolution of δ1(t) (with λ1 < 0); (b) Evolution of δ2(t) (with λ2 > 0).

Fig. 9 .

 9 Fig. 9. Simulation of the four-room case study with Euler-based synthesis method (left) and with the synthesis method of [2] (right).

  the mode j is of the form (j1, j2), with j1 ∈ U1, j2 ∈ U2, U = U1 × U2.

Fig. 10 .

 10 Fig. 10. Simulation of the centralized (left) and distributed (right) Euler-based controllers from the initial condition (22, 22, 22, 22).

  ).

		Euler	DynIBEX
	R [20, 22] Number of balls/tiles 4096	252
	Pattern length	1	1
	CPU time	63 seconds 249 seconds
	Table

2 × [22, 24] 2 S [19, 23] 2 × [21, 25] 2 τ 30 Time subsampling No Complete control Yes Yes maxj=1,...,16 λj -6.30 × 10 -3 maxj=1,...,16 Cj 4.18 × 10 -6

  2 denoted by |W |. see: [Girard: Reachability of uncertain linear systems using zonotopes] [R. Alur, T. Dang, F. Ivancic, Reachability analysis of hybrid systems via predicate abstraction, Hybrid Systems : Computation and Control, C.J. Tomlin, M.R. Greenstreet (Eds), no . 2289 in LNCS, pp 35-48, 2002.] [E. Asarin, T. Dang, A. Girard, Reachability of non-linear systems using conservative approximations, Hybrid Systems : Computation and Control, O. Maler, A. Pnueli (Eds), no. 2623 in LNCS, pp 22-35, Spinger, 2003]

Table 2 .

 2 4 while never going out of[19, 23] 4 . Numerical results for centralized four-room example.

			Centralized
		R		[20, 22] 4
		S		[19, 23] 4
		τ		30
	Time subsampling		τ /20
	Complete control		Yes
	Error parameters	max j=1,...,16	λj = -6.30 × 10 -3
			max j=1,...,16	Cj = 4.18 × 10 -6
	Number of balls/tiles		256
	Pattern length		2
	CPU time		48 seconds
		Sub-system 1	Sub-system 2
	R		[20, 22] 2 × [20, 22] 2
	S		[19, 23] 2 × [19, 23] 2
	τ			30
	Time subsampling	No	τ /10
	Complete control	Yes	Yes
	Error parameters	max j 1 =1,...,4 max j 1 =1,...,4 λ 1 j 1 = -1.39 × 10 -3 max j 2 =1,...,4 γ 1 j 1 = 1.79 × 10 -4 max j 2 =1,...,4 λ 2 j 2 = -1.42 × 10 -3 γ 2 j 2 = 2.47 × 10 -4 max j 1 =1,...,4 C 1 j 1 = 4.15 × 10 -4 max j 2 =1,...,4 C 2 j 2 = 5.75 × 10 -4
	Number of balls/tiles	16	16
	Pattern length	2		2
	CPU time	< 1 second	< 1 second
	Table

Given an initial point x ∈ R, the induced control σ corresponds to a sequence of patterns πi 1 , πi

, . . . defined as follows: Since x ∈ R, there exists a a point xi 1 with 1 ≤ i1 ≤ m such that x ∈ B(xi 1 , δ 0 ); then using pattern πi 1 , one has: φπ i 1 (ki 1 τ ; x) ∈ R. Let x = φπ i 1 (ki 1 τ ; x); there exists a point xi 2 with 1 ≤ i2 ≤ m such that x ∈ B(xi 2 , δ 0 ), etc.

Acknowledgement. The results presented in this paper have been obtained through collaborations with Adrien Le Coënt, Florian De Vuyst, Ludovic Chamoin, Julien Alexandre dit Sandretto and Alexandre Chapoutot. I have also benefited from numerous insightful comments by Antoine Girard.