MECHANICAL BEHAVIOUR OF AN AUSTRALIAN MUDSTONE SUBJECTED TO HEATING TREATMENT AT HIGH-TEMPERATURES

Keywords: Mudstone, Heating Treatment, High Temperature, Mercury intrusion porosimetry, SEM

This study falls in the context of underground coal fires (UCF) where burning coal can elevate the temperature of a rock mass in excess of thousand degrees. The objective of the research was to experimentally investigate the change in mechanical behavior of a rock when subjected to temperatures up to 900 °C. Several blocks of an Australian mudstone were collected from a mining site in New South Wales, Australia. X-ray diffraction and thermal-gravimetric analysis were conducted on powdered rock in order to determine the critical temperatures at which the main mineralogical changes occur. Then, cylindrical specimens were prepared from intact rock blocks for unconfined compression tests. These rock specimens were then treated by heating in a temperature-controlled furnace at temperatures varying from 100 °C to 900 °C. The heating process consisted of two stages: raising temperature with very slow rate of 1 °C per hour followed by a holding period of 24 hours and then a slow cooling process to room temperature. The compression tests were carried out on treated rock specimens to characterize the change in compressive strength. The macroscopical study was complemented by microstructural investigations via optical microscope, scanning electron microscope (SEM) and mercury intrusion porosimetry (MIP). Results show that the unconfined compressive strength of mudstone decreases with increasing heating temperature up to 450 °C, beyond which it increases. This can be explained by the combined effects of micro cracks and mineral thermal reactions caused by heating process.

INTRODUCTION

Underground coal fires (UCFs) are a serious environmental problem of global dimension. These fires, which can result from human activity (in most cases through mining activities) or occur naturally, are problematic for many reasons including loss of valuable non-renewable resources, mine safety, damage to infrastructure and damage to the natural environment [START_REF] Stracher | Coal fires burning out of control around the world : thermodynamic recipe for environmental catastrophe[END_REF][START_REF] Sinha | Spontaneous Coal Seam Fires: A Global Phenomenon[END_REF][START_REF] Kuenzer | Uncontrolled coal fires and their environmental impacts: investigating two arid mining environments in northcentral China[END_REF]. This is why underground coal fires have been of great concern to scientists and engineers [START_REF] Voigt | Spontaneous Coal Seam Fires: The Sino-German Coal Fire Research Initiative[END_REF]. UCFs are characterized by very slow burning kinetics due to the limited availability of oxygen underground and very elevated temperature at the core of the fire because of surrounding ground preventing heat dissipation. Temperatures at combustion center of UCFs of up to 1300 °C have been reported in the literature [START_REF] Stracher | Coal fires burning out of control around the world : thermodynamic recipe for environmental catastrophe[END_REF]. The continuous combustion of the coal seam typically leads to a progressive collapse of the overburden soil/rock mass, that is left unsupported, and to the formation of large cracks, which provide new air circulation paths allowing the fire to propagate [START_REF] Ide | Fissure formation and subsurface subsidence in a coalbed fire[END_REF].

It has been recognized that the existing sets of discontinuities tend to influence the roof collapse and patterns of crack propagation. However, the effect of heat generated by the fire on the mechanical properties of the rock mass has been overlooked. The present experimental work was initiated to investigate the influence of high temperatures on the mechanical behaviour of some Australian rocks, that could possibly be associated with underground coal fires (Liu et al. 2012a, b).

Many studies have been carried out over the last several decades on the impact of high temperatures on the physical and mechanical properties of rocks and this was in different engineering fields such as underground coal gasification or nuclear waste disposal [START_REF] Dmitriev | Physical properties of rocks at high temperatures[END_REF][START_REF] Somerton | Thermal properties and temperature-related behavior of rock/fluid systems[END_REF][START_REF] Sygała | High temperature versus geomechanical parameters of selected rocks-the present state of research[END_REF], to name few). It has been found that apart from the external stress, two key factors govern the changes in rocks when exposed to high temperatures: mineral thermal expansion and mineral thermal reaction [START_REF] Tian | Physical properties of sandstones after high temperature treatment[END_REF]. The former may induce great internal thermal stresses as a consequence of the differences in thermal-expansion characteristics of minerals, resulting in thermal cracking and strength degradation [START_REF] David | Influence of stress-induced and thermal cracking on physical properties and microstructure of La Peyratte granite[END_REF][START_REF] Keshavarz | Damage and changes in mechanical properties of a gabbro thermally loaded up to 1,000 C[END_REF]. On the other hand, mineral thermal reaction may result in some new bonding agents, which can strengthen the rocks [START_REF] Luo | High-Temperature Mechanical Properties of Mudstone in the Process of Underground Coal Gasification[END_REF]. The predominance of one phenomenon over the other depends on the mineral composition of rocks and thermal loading conditions (temperature and its change rate, heating and cooling cycle) as well as external stress state, the mechanical response of the heated rocks would be different as summarized by [START_REF] Zhang | Experimental study on the mechanical properties of rocks at high temperature[END_REF][START_REF] Ranjith | Transformation plasticity and the effect of temperature on the mechanical behaviour of Hawkesbury sandstone at atmospheric pressure[END_REF][START_REF] Brotóns | Temperature influence on the physical and mechanical properties of a porous rock: San Julian's calcarenite[END_REF]. Most studies focused on the change at macroscopic level with only a few studies attempting to investigate the fundamental physic-chemical mechanisms from a microscopic perspective [START_REF] Hajpál | Changes in sandstones of historical monuments exposed to fire or high temperature[END_REF][START_REF] Keshavarz | Damage and changes in mechanical properties of a gabbro thermally loaded up to 1,000 C[END_REF][START_REF] Yavuz | Thermal effect on the physical properties of carbonate rocks[END_REF]).This study investigates the change in unconfined compressive strength of a sedimentary mudstone subjected to heating treatment up to 900 °C. Microstructural analysis techniques (MIP and SEM) were combined to thermal-gravimetric analysis (TGA) to enhance the understanding of the mechanisms of heating treatment on the mechanical response of the mudstone tested. In order to try to simulate the very slow heating process of rocks in UCFs, a relatively slow heating process was used: the heating period proceeded at a rate of 1 °C per hour and was followed by a holding period of 24 hours, before cooling was initiated. This was also deemed to reduce the differential thermal stresses inside the rock generated by rapid heating and thereby the development of microcracks.

Material

A sedimentary mudstone was collected from a mining site in Hunter Valley, New South Wales, Australia. It is referred to as HV mudstone. The basic physical properties of the intact (i.e. no heat treated) HV mudstone and its main mineralogical composition are presented in Table 1. The presence of dominant minerals such as Quartz and Muscovite was confirmed by microscopic imaging analysis on a thin section of HV mudstone (Fig. 1). The microscopic images also indicate that the maximum particle size is less than 100 µm and there are some evident bedding planes existing in the mudstone (Fig. 1b).

EXPERIMENTAL METHODS AND TESTING PROGRAM Experimental Methods

Thermal-Gravimetric Analysis

To highlight the mineralogy changes of HV mudstone in response to heating, thermal-gravimetric analysis was performed on the oven-dried mudstone powder with particle size below 75 µm. The tests were run under Argon environment whilst raising the temperature of 20 °C to 1100 °C at a rate of 5 °C/min. The mass loss and heat flow were recorded during the heating process, which enables the determination of the critical temperatures at which the dominant mineral thermal reactions take place.

Specimen preparation and mechanical testing

The rock samples were cored perpendicular to the bedding planes (as shown in Fig. 1b) with a diameter of 12 mm and a height of 24 mm. Note that these dimensions comply with the ASTM standard D-4543 (2008), which states that the diameter should be at least ten times the size of the biggest particles. With the particles smaller than 0.1 mm, 12 mm diameter satisfies this criterion. Such small dimensions helped limiting the variability of the material. The specimens were first dried in an oven at a temperature of 100 °C (reference state). Then, the dried specimens were subjected to heating treatment at different target temperatures (300, 450, 600, 750 and 900 °C) by using a temperature-controlled furnace. The heating treatment consisted of a slow heating period from 20 °C to target temperatures at a rate of 1 °C/hour and a holding period of 24 hours followed by a cooling process at a rate of 25 °C/hour back to 20 °C. Note that the rock specimens were periodically weighted during the heating process, in order to monitor the mass loss. Uniaxial compression tests were conducted on the heated specimens according to ASTM standard D7012-2014. The loading rate of these tests was 0.8 mm/hour.

Mercury intrusion porosimetry (MIP) and Scanning electron microscope (SEM) observations

The MIP analyses were conducted using a Micromeritics AutoPore (IV 9500) with the mercury pressure incrementally increased to the maximum value of the apparatus (228 MPa), which corresponds to an equivalent entrance pore diameter of approximately 0.0065 µm. As for the SEM observations using a Philips XL30 SEM, Oxford ISIS EDS, the specimens were coated by gold powder in order to obtain the better contrast of SEM image prior to SEM observations. Both MIP tests and SEM observations were carried out on the specimens after heating treatment

Testing Program

Table 2 summarizes all the tests performed in this study, including TGA, UCS tests and microstructural analysis.

RESULTS AND DISCUSSION

Variation of Mass Loss Ratio with Temperature

Fig. 2 presents the mass loss ratios (in percent of initial mass) of HV mudstone in response to an increase in temperature, which were respectively determined by thermal-gravimetric analysis (20 to 1100 °C) and from heating treatment using a temperature-controlled furnace (100 to 1200 °C). In order to more accurately determine the critical temperatures at which the dominant mass loss occurred, the differential thermal-gravimetric (DTG) evolution was also shown in Fig. 2, which was obtained by simply differentiating the mass loss ratio curve.

It can be noticed that the mass loss ratio from TGA slightly decreases over the temperature range below 400 °C with the sharpest drop located between 400 and 600 degrees. 3. Thermal reactions of the dominant minerals of HV mudstone [START_REF] Dmitriev | Physical properties of rocks at high temperatures[END_REF][START_REF] Somerton | Thermal properties and temperature-related behavior of rock/fluid systems[END_REF] Interestingly, a small proportion of mass gain can be noticed when the heating temperature exceeds 900 °C. A similar trend was observed for the specimens burnt in the furnace although there exists a clear offset of about 1%. This offset is believed to be due to the type of specimen: the TGA requires powdered rock with a higher surface are of material exposed to the heat, resulting in higher mass loss ratio.. To understand the mechanism related to phenomenological mass loss, possible temperature alterations on the dominant minerals of the mudstone tested is summarized in Table 3 as shown by [START_REF] Dmitriev | Physical properties of rocks at high temperatures[END_REF][START_REF] Somerton | Thermal properties and temperature-related behavior of rock/fluid systems[END_REF]. We can clearly see that the main mass loss is likely to be associated with the loss of structural water from the muscovite and kaolinite, which are characterized by the DTG peak of 502 °C (Fig. 2). These reactions can possibly produce some bonding agents enabling to strength the mudstone as shown by Hajpál & Török (1998) and [START_REF] Hajpál | Changes in sandstones of historical monuments exposed to fire or high temperature[END_REF].

Variation of Unconfined Compressive Strength (UCS) with Temperature

Figure 3 shows the evolution of UCS and corresponding axial strain at peak strength of the burnt HV mudstone with increase in temperature. Interestingly, the evolution is non-monotonic: the UCS first decreases from 100 °C to °C before increasing over the reference strength (i.e. at 100 degrees). The values of strain at failure are more scattered and no clear trend can be defined. The mechanical response of the burnt mudstones may be explained by the change in their mineralogical composition and their microstructure. The former is governed by the thermal reactions of two dominant minerals (muscovite and kaolinite), which may create some cementation effects during heating at temperature above 450 °C, as discussed above. For the temperatures below 450 °C, the significant reduction in UCS may be imputable to the creation of micro cracks by heating treatment as suggested in some studies (e.g. [START_REF] Wolf | The interaction between underground coal fires and their roof rocks[END_REF]). This will be further verified by the microstructural analysis discussed in the following sections. Figure 4 presents the MIP results obtained on the mudstone specimens after heating treatment at temperatures of 100, 450, 600 and 900 °C while the corresponding SEM images appear in Figure 5. Figure 4 clearly shows that the total void ratio intruded by mercury intrusion increases by a factor 2 (from 0.08 to 0.17) as the heating temperature rose from 100 to 900 °C. This corresponds to a shift in dominant pore size from 0.1 µm to 1 µm. These changes can be explained by two basic mechanisms: (1) the differential thermal expansion of the different minerals (Table 3), (2) thermal reactions of the dominant minerals (Muscovite and kaolinite), resulting in the loss of structural water. The first one may induce the creation of microcracks inside the mudstone that are visible on the SEM images (Fig. 5b,5d), leading to the strength degradation. As for the thermal reactions of the two minerals, it tends to create the sintered phases shown in Fig. 5c and5d, which results in a strengthening effect. The predominance of one mechanism over the other governs the final strength changes of the mudstones in response to heating treatment. For instance, the second mechanism seems to play a more important role in strength changes at temperatures over 450 °C, at which the onset of dominant thermal reactions took place as shown in Table 3 and Fig. 2. This explains the non-monotonic UCS evolution and the transition from the decrease to the increase in UCS.

CONCLUSIONS

This study aimed at investigating the thermal effects on UCS of a sedimentary mudstone collected from a Hunter Valley site close to an existing underground coal seam fire. The compression tests were combined to thermal-gravimetric analysis and microstructural analysis that were employed to understand the mineralogical and microstructural changes associated with increasing temperatures. The tests showed that the UCS changes non-monotonically with temperature and this can be explained by two dominant mechanisms: differential thermal expansion of the minerals (Quatz, Muscovite, kaolinite and Albite) and thermal reactions of muscovite and kaolinite over the temperature range of interest. The predominance of one mechanism over the other governs the final strength changes of the mudstones in response to heating treatment. The first mechanism tends to induce microcracking, and increase in porosity as evidenced respectively by SEM observations and mercury intrusion porosimetry (Fig. 5, Fig. 4).

The second mechanism consists of dehydration and dehyoxylation processes, in other words, the loss of strongly adsorbed water and the structural water. The dehyoxylation occurs at temperatures above 450 °C and can create some sintered phases, which results in strengthening of the mudstone. Indeed, the UCS doubled when the temperature went from 450 °C to 750 °C. Future work will include investigations on another sedimentary rock (sandstone) and on the thermal effects on the stress-strain relationships under triaxial conditions as well as tensile strength using Brazillian test. This will provide a more comprehensive dataset on failure criterion of sedimentary rocks under high temperatures.

Figure 1 -

 1 Figure 1-(a) Microscopic view of minerals in a thin section and (b) Macroscopic view of HV mudstone

Figure 2 .

 2 Figure 2. Comparison of mass loss ratio of HV mudstone determined from TGA and heating treatment using temperature-controlled furnace Table3. Thermal reactions of the dominant minerals of HV mudstone[START_REF] Dmitriev | Physical properties of rocks at high temperatures[END_REF][START_REF] Somerton | Thermal properties and temperature-related behavior of rock/fluid systems[END_REF] Dominant minerals Critical temperatures (°C) Thermal reactions

Figure 3 .

 3 Figure 3. Variation of UCS and axial strain at peak strength of burnt HV mudstone with temperature Change in Microstructure with Temperature

Figure 4 .

 4 Figure 4. MIP results of four mudstone specimens burnt at temperatures of 100, 450, 600 and 900 °C: (a) Variation of total mercury intrusion void ratio with entrance pore diameter (b) Pore size distribution

Table 1 .

 1 Some physical properties and dominant minerals of HV mudstone Specific gravity (G

s ) Bulk density (g/cm 3 ) Initial void ratio (e) Mineralogical composition from Semiquantitative XRD analysis 2.6567 ± 0.002 2.38 ± 0.03 0.11 ± 0.01 Quartz (66%), Muscovite (12%), Albite (7%), Kaolinite (9%), Anorthite (4%), Vermiculite (< 1%)

Table 2

 2 Testing program (TGA, UCS test, MIP and SEM)

	Specimen heating temperature (°C)	TGA	Unconfined compression test	MIP	SEM observation
	100	1	2	1	1
	300	-	1	-	-
	450	-	1	1	1
	600	-	2	1	1
	750	-	2	-	-
	900	-	1	1	1


Florent Vieux-Champagne

Stéphane Grange

Yannick Sieffert

Patrice Garcia

Carole Faye

Jean-Charles Duccini

Laurent Daudeville

NUMERICAL ANALYSIS OF TIMBER-FRAME STRUCTURES WITH INFILL UNDER SEISMIC LOADING

Keywords: Timber-frame structures with earth and stone infill, Multi-scale approach, Seismic behavior, Constitutive model

Timber frame with infill have proved to be efficient structures when submitted to seismic hazard. No building code and very few scientific studies are available to design this kind of structures. Therefore, the present work aims at predicting the behavior of timber-frame structure with earth infill under seismic loading. Four scales are defined: the joint, the elementary cell, the shear wall and the entire house. At each scale, both an experimental and a numerical analysis are carried out in order to bring extensive knowledge with regard to the behavior of these structures under seismic loading and to validate a numerical tool able to predict their response.

INTRODUCTION 123

Timber frame structures with infill represent an important share of the traditional built environment over the world and especially in earthquake-prone areas ([ DUT 12], [ OZT 10], [VIE 14]). Their construction is environmentally friendly due to the proximity of the building materials and coherent with the respect of local culture. Furthermore, during recent earthquakes (Turkey 1999, Greece 2003, India 2005 and Haiti 2010), timber-framed structures filled with natural materials (earth, stones) often showed a better resistance to earthquakes compared with imperfectly built reinforced concrete structures ([ AUD 11], [LAN 08]). Nevertheless, few scientific studies (e.g. [ AHM 12], [CEC 10], [ VAS 13]) and no building code are available to assess the seismic vulnerability of this kind of building. Therefore the present work aims at predicting the behavior of timber-frame structures with infill under seismic loading from the specific case of a Haitian timbered masonry structure, called "Kay peyi", which is reused in one of Haiti's ongoing reconstruction projects (Figure 1).

The study is based on the definition of different scales of analysis, from the connections to the global structure, as already proposed in [HUM 14], [RIC 02] and [XU 09] for timber-frame shear wall. Herein, the same approach is defined by means of four scales (cf. Figure 2): the joint, the elementary cell, the shear wall and the entire house. For each one, both an experimental and a numerical analysis are performed. The scale 1 allows calibrating Florent Vieux-Champagne, Email: florent.vieux.champagne@gmail.com 1 CRAterre, AE&CC Research Unit, National School of Architecture of Grenoble, France 2 Univ. Grenoble Alpes, 3SR, F-38000 Grenoble, France, CNRS, 3SR, F-38000 Grenoble, France 3 FCBA, French Technological Institute for Forestry, Cellulose, Timber Construction and Furniture numerical models that are used to model the structure at the scale 2 and predict its behavior regarding the experimental results. From this scale, a simplified model (macro-element) is developed and is used to model the structure at the scale 3 and the scale 4. At each one, the model is validated from the experimental results. One of the aims of this work is to develop a simplified numerical tool able to predict the dynamic behavior of a timber frame structure with infill. This tool has been implemented in a Finite Element Code developed in Matlab®. Built for a reconstruction project

DESCRIPTION OF THE STRUCTURE

The timber structure is braced by X-cross and the spaces are filled with a stone masonry linked by an earth mortar. The earth is composed of a 1/1 volume ratio of a claycalcareous mix (D < 125 mm), a 2/1 ratio of sand (D < 2 mm), a 1/2 ratio of water (depending on the water content of sand) and a 1/1 ratio (very coarse measure) of sisal fibers (easily available in Haiti). To prevent the mortar from cracking, the water content must be limited in the earth and moreover sisal fibers are to be added. The claycalcareous mix is chosen to obtain a mortar as close as possible to what can be produced from the local earth.

The fibers serve to limit the appearance of cracking and its propagation in mortar. Plain shank nails are thus placed inside wood triangles to improve the bond between both parts.

EXEPRIMENTAL AND NUMERICAL ANALYSIS 3.1 SCALE 1: JOINTS

More than 30 experimental tests on coiled strapping joints have been performed to provide input data for the numerical models of structural elements. These connections have been tested through monotonic and cyclic tests, both in tension and shear, according to [EN-12512 02].

The results of the tests are used to calibrate the constitutive model of Humbert ([HUM 14]) (cf. Figure 3). Two levels of calibration are performed. The first level is a direct calibration, which consists in reproducing one particular test Fig. 3a). The second level is an average calibration, which consists in calibrating the parameters to reproduce

SCALE 2: ELEMENTARY CELL

This structure is an elementary part of the shear wall. The goal of studying this scale is motivated by two facts:

• Firstly, because of the lack of knowledge regarding the structural behavior of filled timber frame structures, a parametric study has been performed to assess the impact, on the lateral load capacity of the structure, of the presence of infill, of the kind of infill (stones or adobes masonry), of the presence of bracing. The results show that the strength of the structure is not governed by these parameters but more by its lateral stiffness and its dissipation energy capacity. This is explained by the lattice behavior-like depicted by the elementary cell ([VIE 14]).

• Secondly, these results allow validating the development of the FE models.

A detailed FE model of elementary cell is developed using beam elements for the frame and two-node spring-like finite elements for the joints and for the experimental friction (Figure 4a). The experimental behavior of the elementary cell, represented in Figure 4b by the positive (+) and the negative (-) average (A) envelop curves, is then predicted from the detailed model. From the results of this analysis, a simplified FE model (macro-element) is developed and calibrated (Figure 5a). The macro-element is composed of a frame of bars ensuring the hypothesis of a parallelogram-like deformation of the wall, therefore modeling only the shearing behavior. The in-plane horizontal degree of freedom of the simplified element is modeled by the previously described constitutive hysteretic behavior law, whose parameters are identified thanks to the results of the detailed model (Figure 5b).

SCALE 3: SHEAR WALL

Experimental campaign has been realized at the CNR-Ivalsa (Trento, Italy) where three filled walls (two reversed cycle tests and one monotonic) and two empty walls (one reversed cycle test and one monotonic) were tested:

At this scale, the objectives of the experimental study are the following:

• Firstly, to perform a comparison of the behavior of a timber shear wall with and without filling. The tests highlighted the influence of the filling on the deformation shape of the shear wall and hence its influence on the resistance of the wall ([VIE 14]). • Secondly, to validate a simplified FE-model built by the combination of the macro-elements. A simplified representation of the finite element model of the wall is depicted in Figure 6a. The vertical posts are modeled by continuous beam which have their extremities hinged by static condensation. The horizontal members in the middle of the wall are modeled by bars. The prediction of the behavior of the shear wall is shown in Figure 6b.

The experimental behavior is well predicted by the modeling and hence validates this simplified approach.

SCALE 4: HOUSE

A dynamic test has been performed on a timber-frame entire house with walls filled with earth and stones. The test aimed at demonstrating the safe behavior of this structure under seismic loading.

Experimental set-up Shaking table

The dynamic tests were performed on the uniaxial shaking table of the Mechanical Laboratory of the French Institute FCBA. The table is composed of a 6m × 6m aluminium platform moved by a 250kN servo-hydraulic actuator.

Specimen

The house is a square of 4.65m × 4.65m wide, 3.20m high (at the ridge) and symmetrical according to the N-S axe (solicitation direction, see Fig. 8). The shear walls have six vertical posts with a window and a door respectively between the post 2-3 and 4-5 in the N-S direction. The transversal walls are identical with also six vertical posts and have one window between the posts 3-4. According to a preliminary numerical study, the specimen was placed on the shaking table in order to solicit its weaker direction. It was anchored in steel beams by means of 40 bolts (2 between 2 vertical posts) which are themselves bolted on the table.

Experimental procedure

As said above, the test aimed at demonstrating the safe behavior of this structure under seismic loading. For this purpose, the 100% seismic signal was designed according to the empirical green's function method to be representative of the 2010 Haiti earthquake (see Fig. 7).

Then, in order to analyze the non-linear behavior of the structure, the house was submitted to the 200% and 300% signal. To assess the fundamental frequency evolution, the structure was submitted to low-level white noise (0.03g and 0.5mm RMS (Root Mean Square)) between each seismic loading. The speed camera (Phantom v641 camera with 28mm lens, 4 megapixels) was used to stare the East shear wall, at a speed of 150 fps (frames per second), in order to obtain the absolute displacement field of the panel by DIC. The acquisition rate was determined so that it is higher than the frequencies involved in the seismic signal (< 30 Hz that is to say minimum 30 fps). The speckle pattern was made by the projection of black, gray and white paint. The average speckle dot diameter was about 1cm. The lighting (4000W, ). The accuracy of the measurement is about 0.1 pixel or 0.22mm across the building.

Experimental results

Main results obtained from the analysis of the seismic tests are presented below.

Mode shape and evolution of the fundamental frequency of the building

Before the seismic tests, a modal analysis has been performed by measuring three axis acceleration response at each node of the structure to a white noise excitation. It allowed obtaining the three main mode shapes, natural frequencies and damping ratios (see Tab. 1). The two first modes are similar and are the transversal walls racking mode. The small frequency difference is due to the different positions of the walls (one near the doors (more flexible), one near the windows) and due to the intrinsic differences of the two walls workmanship (weight, geometrical dimensions, etc.). The second main mode is the vibration mode of the shear walls. After each seismic test, the structure was submitted to a low-level white noise. An analysis of the structure modification resulting from the successive seismic tests was performed by computing the Frequency Response Function (FRF) between the output acceleration at node 7 (ACC7) and the input table acceleration (ACC10) obtained during each low-level white-noise input test (Fig. 8). As a reminder, the FRF is a representation in the frequency domain of the ratio between an output and input acceleration signals. The parameters observed are the variation of the natural frequencies values, shape and amplitude. The successive FRF obtained during each lowlevel white-noise input test are presented in Figure 9. On this graph, we can notice that progressively after each seismic test:

• The natural frequency, representative of the structure stiffness, decreases. • The width curve around the natural frequency peak, representative of the structure damping, increases. • The peak amplitude decreases.

DIC results versus DWDS measures

As explained above, it was therefore necessary to validate the results obtained by DIC. The Figure 10 depicts the comparison of the displacement obtained on the shaking table by the LVDT and by the DIC (for the 300%, Fig. 11). The correlation between the results from both measurement methods is excellent and thus validates the use of the CIN.

Numerical modeling

Using a simplified FE-model based on macro-elements defined at scale 2, the behavior of the house is predicted and compared to the experimental results for the three earthquakes in a row. The Figure 8a shows this comparison for the top displacement of the shear wall (Fig. 7) for the full duration of the 200% signal and for the time interval of this signal during which the acceleration is maximal (Figure 8b). These results validate the multi-scale approach and validate the numerical tool used to predict the behavior of the timber frame structure with infill under seismic loading.

CONCLUSIONS

This study presents both experimental and numerical results at different scales of analysis to assess the behavior of timber-frame structures with stones and earth infill under seismic loading. The first scale provides extensive information on the nonlinear dissipative phenomena in the structure. The second scale provides a numerical non-linear model of the elementary cell that consists in a braced square filled with earth and stones. The model of a shear wall or an entire house by means of macro-elements defined at scale 2 and whose constitutive model is built with the knowledge of joints behavior, provides a computationally efficient and yet relatively accurate tool for the analysis of the seismic vulnerability of timberframed houses with earth infill.

Figure 1 :

 1 Figure 1: Rural Haitian house: (a) Traditional "Kay peyi" (b) Built for a reconstruction project

  Coiled strapping joints consisting of a punched steel strip surrounding the wood parts (in a T-shape) and fastened by few nails that govern behavior of the connection. The timber strength class is C18 with a density = 380 kg/m according to European standard EN 338 (see [EN-338 03]). The strip is a FP30/1,5/50 from Simpson Strong Tie®, 30 mm wide and 1,5 mm thick. Eight 3 mm × 70 mm plain shank nails are used (see [EN-10230-1 00]) to fasten the strip. The bracings are made from Saint Andrew's crosses with one continuous diagonal and another diagonal divided into two parts. They are connected to the vertical post using two 70-mm long wooden nails.

Figure 2 :

 2 Figure 2: Multi-scale approach

Figure 3 :Figure 4 :Figure 5 :

 345 Figure 3: Reversed cyclic tests performed on a connection: (a) Fitting of the Humbert model's parameters from the results of a shear reversed cyclic test (b) Fitting of the Humbert model's parameters from the results of a tension monotonic test

Figure 6 :

 6 Figure 6: Scale 3: simplified modeling of the shear wall by macro-element combination and results: Experimental and numerical results

Figure 8 :

 8 Figure 8: Measurement devices positionMeasurement devicesMeasurement devices used were placed as depicted on the Figure8. The accelerometers allowed performing an analysis of the natural frequencies evolution, the draw wire displacement sensors (DWDS) and the LVDTs measurement were compared to the results of the Digital Image Correlation (DIC) made through high speed camera pictures and compared to the results of the numerical modeling. To the knowledge of the authors, it is the first time that a DIC is used to analyze this kind of structure at this scale. For this reason, a focus is made on the set-up used to perform the DIC.

Figure 7 :

 7 Figure 7: Seismic signal used: (a) Time history acceleration (b) Pseudo-acceleration elastic response spectra (5% damping)

Figure 9 :

 9 Figure 9: Frequency Response Function

Figure 10 :Figure 11 :

 1011 Figure 10: Points whose displacement is analyzed by the DIC

Table 1 :

 1 Natural frequencies and the equivalent viscous damping ratio

	Poles	Element	Natural	Damping
			frequency (Hz)	ratio (%)
	1a	Transversal wall	5.2	5.2
	1b	Transversal wall	5.9	5.7
	2	Shear wall	10.9	5.3

The authors also thank the technical staffs of the FCBA, of the Trento/CNR Ivalsa and of Grenoble/UJF for the tests they conducted.

Gratitude is extended to Simon Pla, Patrice Clavel and Laurent Munier for their valuable assistance during the experimental program at the scale 2 and 4.

ACKNOWLEDGEMENT

The authors wish to thank and acknowledge the support of the French National Research Agency (ANR) under reference code ANR-10-HAIT-003, the Haitian association GADRU, the participating associations of the PADED platform and all local partners for their involvement and participation, contributing to this research project.


F Marinelli

Y Sieffert

R Chambon

Hydromechanical modelling of an initial boundary value problem: Studies of nonuniqueness with a second gradient continuum

A non-uniqueness study for a hydromechanical boundary value problem is performed. A fully saturated porous medium is considered using an elastoplastic constitutive equations to describe the mechanical behavior of the skeleton. A real hydromechanical experiment which consists in a hollow cylinder test on a Boom Clay sample is modelled. It is shown that the time step discretisation of the numerical problem has an effect on the initialisation of the Newton-Raphson algorithm on a given time step. Different solutions for the same initial boundary value problem can consequently be found.

Introduction

In the simulation of initial boundary value problems using constitutive equations for geomaterial behavior, it is well known that some difficulties can arise, particularly if degradation of the materials occurs. These problems have been studied for the case of single phase materials and some theoretical results have been established. Within the small strain assumption and for a rate problem, the uniqueness can be proved using the so-called Hill exclusion functional [START_REF] Hill | Aspects of invariance in solid mechanics[END_REF]. From a numerical standpoint, in case of non-uniqueness, the solution can be influenced by different numerical inputs such as the spatial discretization of the problem, the size of the time step, the mesh spacing, some tolerance values and so on (see [START_REF] Chambon | An algorithm and a method to search bifurcation points in non-linear problems[END_REF]).

In this study, a real hydromechanical experiment is modelled using a finite element code Lagamine (University of Liège). A local hydromechanical second gradient model (see [START_REF] Collin | A finite element method for poro mechanical modelling of geotechnical problems using local second gradient models[END_REF]) is considered to obtain regularised solutions of this initial boundary value problem. No other regularization is necessary in terms of the fluid phase or the solid-fluid coupling.

We show that, changing the time discretisation of the numerical problem, different solutions can be found. In function of the shear bands path during the time, a classification of all computed solutions is proposed.

Local hydromechanical second gradient model

It is well known that classical continuum medium cannot be considered to model the localization of plastic strains. An internal length has to be introduced into the model in order to obtain mesh independent solutions of an initial boundary value problem. In this respect, the problem can be considered regularised. The regularization (that means the introduction of an internal length into the model) cannot restore the uniquess of the solution as it is proved in [START_REF] Chambon | One-dimensional localisation studied with asecond grade model[END_REF] for a local second gradient continuum model.

In the framework of microstructure continuum [START_REF] Germain | The method of the virtual power in continuum mechanics. Part 2: Microstructure[END_REF], a local hydromechanical second gradient model can be considered [START_REF] Collin | A finite element method for poro mechanical modelling of geotechnical problems using local second gradient models[END_REF]. The problem is solved by enforcing (in a weak form) the balance equations for the mixture (equation 1), and the mass conservation for the fluid (equation 2). Both equations are written in the current solid configuration denoted (updated Lagrangian configuration).

(

(

In equation ( 1), Terzaghi's principle of effective stress is assumed to compute the Cauchy stress of the mixture . The mass density of the mixture at the time t is defined as follows:

(3)

where is the solid phase density (assumed to be incompressible, i.e., is a constant); is the fluid density and is the porosity defined as . is the current volume of a given mass of the solid skeleton and the corresponding volume of fluid. The fluid is assumed to be compressible. All variables appearing in equations ( 1) and ( 2 To solve the equations ( 1) and ( 2) at time with a Newton-Rapshon algorithm a consistent linearisation of the system has to be performed. All the developments of the linearization for a local hydromechanical second gradient model are detailed in [START_REF] Collin | A finite element method for poro mechanical modelling of geotechnical problems using local second gradient models[END_REF]. The classical part of the constitutive equation, giving the effective stress of the skeleton, is described by a non-associated elasto-plastic model based on the Drucker-Prager yield field (Plasol model, [START_REF] Barnichon | Finite element modelling in structural and petroleum geology[END_REF]). Cohesion and friction angle for a compressive triaxial path , represent the hardening/softening variables for this elasto-plastic constitutive equation. The only internal variable that characterizes the hardening/softening behavior of the material is the Von Mises equivalent plastic strain , which describes degradation due to the shear. Dilatancy is considered as a constant during plastic flow. A hyperbolic variation of the friction angle and cohesion is proposed as hardening law. The elastic behavior is modelled by isotropic linear elasticity. This constitutive equations is computed with a return mapping algorithm [START_REF] Ortiz | An analysis of a new class of integration algorithms for elastoplasticconstitutive relations[END_REF]. The other part (giving the double stress) is a particular case of a more general isotropic linear relation derived by [START_REF] Mindlin | Micro-structure in linear elasticity[END_REF] that involves five independent parameters, in the interests of simplicity only one parameter (namely D, see [START_REF] Bésuelle | Switching deformation modes in postlocalizationsolutions with quasi-brittle material[END_REF]) is used.

Non-uniqueness studies

In order to perform our non-uniqueness studies, we use realistic material properties related to the behavior of Boom Clay [START_REF] Horseman | Geotechnical characterization of boomclay in relation to the disposal of radioactive waste (page 87)[END_REF]. In our study we choose to model a real hydromechanical hollow cylinder experiment performed in Lausanne [START_REF] Labiouse | Hollow cylinder simulation experiments of galleries inboom clay formation[END_REF]. This experiment is divided in two phases, that are described in figure 1. A non-uniqueness study is carried out for Plasol constitutive equations by considering phase A of boundary conditions that represents a hydromechanical unloading. In phase B, the evolution of the shear bands with time is modelled. All computations presented are performed under plane strain hypothesis to model the behavior of a horizontal section of the hollow cylinder.

From a numerical point of view, non-uniqueness can be studied by changing the initialization of the Newton-Raphson algorithm for a given time step. In our case, the initialization of the numerical algorithm is obtained from the nodal velocity field (that one of the previous step) at the first iteration of the current step. [START_REF] Chambon | One-dimensional localisation studied with asecond grade model[END_REF] where represents the step size, is the nodal coordinate at the previous step and is the nodal coordinate for the first iteration of the current step. Different sequences of time steps can cause a change in the initialisation of the numerical problem for each time step. In this way, non-uniqueness can be studied for the same initial boundary value. In order to visualize the non-uniqueness of this initial boundary value problem, the evolution of the shear strains are shown in terms of the loading index. This quantity is defined as follows: if a Gauss point undergoes plastic loading, a small blue square is plotted, otherwise (in the case of elastic unloading or reloading of the Gauss point), a small yellow square is plotted. Basically, all the localized deformation patterns present two kinds of shear band configurations: the first is characterized by some crossing shear bands (the socalled conjugated shear bands), while the second is characterized by spiral shear bands (the so-called periodic shear bands).

Among the great number of computations performed five solutions, which are representative of the non-uniqueness problem, are shown in figure 2. All the solutions are presented in Figure 3 and can be classified qualitatively according to shear bands evolution as follows:

 Computations with radial symmetry ISO-0.

 Computations with shear bands that start from a thin plastified region around the inner hollow cylinder, which then grow towards the external part of the sample (solutions ISO-1 and ISO-4).

 Computations with shear bands that appear by an elastic unloading in a thick plastified zone of the hollow cylinder at the end of phase A (solutions ISO-2 and ISO-3).

The development of the shear bands during the phase B of the boundary conditions is shown in figure 3.

Conclusions

In this paper, the different solutions of the same initial boundary value problem (hollow cylinder experiment) for a fully saturated porous medium, modelled with an elastoplastic model, have been studied. All solutions have been obtained changing the time discretization of the problem. The solutions obtained with Plasol models can be classified according to the evolution of shear bands configurations in time.

It can be concluded from the present work that, as soon as degradation is incorporated in a constitutive equation, uniqueness of the solution of initial boundary value problems is questionable even in coupled problems. This work extends to the case of coupled problems the finding that the use of an enhanced second gradient model is not sufficient to restore the uniqueness of the solutions (see [START_REF] Chambon | One-dimensional localisation studied with asecond grade model[END_REF]).

Fig. 1 .

 1 Fig. 1.The initial boundary conditions on the inner and external side of the hollow cylinder: a continous line represents boundary conditions on the inner hollow cylinder, a dot line on the external hollow cylinder.

Fig. 2 .Fig. 3 .

 23 Fig. 2.Loading index of five solutions during phase A of boundary conditions


Thanh-Phong Pham
email: thanh-phong.pham2@gipsa-lab.grenoble-inp.fr

Olivier Sename
email: olivier.sename@gipsa-lab.grenoble-inp.fr

Luc Dugard
email: luc.dugard@gipsa-lab.grenoble-inp.fr

Design and Experimental Validation of an H ∞ Observer for Vehicle Damper Force Estimation

Keywords: Semi-active suspension, H ∞ observer, damping force estimation, Lipchitz condition

The real-time estimation of damper force is crucial for control and diagnosis of suspension systems in road vehicles. In this study, we consider a semi-active electrorheological (ER) suspension system. First, a nonlinear quarter-car model is proposed that takes the nonlinear and dynamical characteristics of the semi-active damper into account. The estimation of the damper force is developed through a H ∞ observer whose objectives are to minimize the effects of bounded unknown road profile disturbances and measurement noises on the estimation errors of the state variables and nonlinearity through a Lipschitz assumption. The considered measured variables, used as inputs for the observer design, are the two accelerometers data from the sprung mass and the unsprung mass of the quarter-car system, respectively. Finally, the observer performances are assessed experimentally using the INOVE platform from GIPSAlab (1/5-scaled real vehicle). Both simulation and experimental results emphasize the robustness of the estimation method against measurement noises and road disturbances, showing the effectiveness in the ability of estimating the damper force in real-time.

INTRODUCTION

Nowadays, semi-active suspensions are widely used in vehicle applications due to their advantages compared to active and passive suspensions [START_REF] Savaresi | Semi-active suspension control design for vehicles[END_REF] and references therein). Central issues of these applications include dynamic modeling and control designs based on a reduced number of sensors to improve vehicle comfort and road holding. Depending on the characteristics of the adjustable shock absorber, models have been derived using several methodologies with different complexity and accuracy. Main models may be classified in terms of static and dynamic characteristics. Static models include Bingham model with Coulomb friction (see [START_REF] Stanway | Nonlinear modelling of an electro-rheological vibration damper[END_REF]), hysteresis based model (see [START_REF] Guo | Dynamic modeling of magnetorheological damper behaviors[END_REF], [START_REF] De J Lozoya-Santos | Magnetorheological damperan experimental study[END_REF]). Dynamic models group considers the Bouc-Wen model in [START_REF] Wen | Method for random vibration of hysteretic systems[END_REF], [START_REF] Ahmadian | No-jerk skyhook control methods for semiactive suspensions[END_REF] [START_REF] Spencer | Phenomenological model for magnetorheological dampers[END_REF]). Based on several approximators such as neural network [START_REF] Chang | Neural network modeling of a magnetorheological damper[END_REF], [START_REF] Chang | Neural network emulation of inverse dynamics for a magnetorheological damper[END_REF]), fuzzy [START_REF] Schurter | Fuzzy modeling of a magnetorheological damper using anfis[END_REF]), polynomial [START_REF] Du | Semi-active h control of vehicle suspension with magneto-rheological dampers[END_REF]) and among others [START_REF] Savaresi | Identification of semi-physical and black-box non-linear models: the case of mr-dampers for vehicles control[END_REF], proposed black box model can be also divided into static or dynamic groups, depending on the typical model. Generally, the common drawback of the above-mentioned dynamic models is the high complexity for controller and observer synthesis. [START_REF] Guo | Dynamic modeling of magnetorheological damper behaviors[END_REF] introduced a semi-phenomenological model with high accuracy and control-oriented design. Noticeably, Guo's model only characterizes the static nonlinear behavior of semi-active suspensions. Therefore, how to extend the model taking into account both nonlinear and dynamic behaviors of the ER suspension system is needed.

Based on the above models, many control designs were proposed in the literature (see [START_REF] Poussot-Vassal | A new semi-active suspension control strategy through lpv technique[END_REF], [START_REF] Priyandoko | Vehicle active suspension system using skyhook adaptive neuro active force control[END_REF] and a review in [START_REF] Poussot-Vassal | Survey and performance evaluation on some automotive semi-active suspension control methods: A comparative study on a single-corner model[END_REF]). Some control design methodologies considered the damper force as the control input of the suspension system, then using an inverse model for implementation (see for instance [START_REF] Do | An lpv control approach for semi-active suspension control with actuator constraints[END_REF], [START_REF] Nguyen | Semi-active suspension control problem: Some new results using an lpv/h state feedback input constrained control[END_REF]). Others use the force tracking control schemes in order to attain control objectives [START_REF] Priyandoko | Vehicle active suspension system using skyhook adaptive neuro active force control[END_REF]). Indeed, the damper force signal plays an important role in control synthesis. For this reason, several damper force estimation methodologies were presented (see [START_REF] Estrada-Vela | H observer for damper force in a semi-active suspension[END_REF], [START_REF] Reichhartinger | Robust estimation of forces for suspension system control[END_REF][START_REF] Tudon-Martinez | Parameter-dependent h filter for lpv semi-active suspension systems[END_REF], [START_REF] Koch | Nonlinear and filter based estimation for vehicle suspension control[END_REF], [START_REF] Rajamani | Adaptive observers for active automotive suspensions: theory and experiment[END_REF]), since in practice the damper force measurement is difficult and expensive. In [START_REF] Koch | Nonlinear and filter based estimation for vehicle suspension control[END_REF]), parallel Kalman filters were developed to estimate the damper force without considering the dynamic behavior of the semi-active damper. [START_REF] Estrada-Vela | H observer for damper force in a semi-active suspension[END_REF] introduced the H ∞ damping force observer using a dynamic nonlinear model of the ER damper, while requiring three sensors as inputs of observer. In order to reduce the number of the sensors and maintain a consideration of dynamic and nonlinear characteristics of MR damper, [START_REF] Tudon-Martinez | Parameter-dependent h filter for lpv semi-active suspension systems[END_REF] proposed a LPV-H ∞ filter to estimate the damper force using deflection and deflection velocity signals, which are difficult and expensive to measure in practice. Moreover, based on accelerometers, a full-car observer using the linearized model of the damper is proposed in [START_REF] Dugard | Full vertical car observer design methodology for suspension control applications[END_REF]) and gives interesting results both in simulations and experiments. Despite these achievements, the damping force observer based on the dynamic nonlinear model of semi-active ER suspension system and using low-cost sensors is still in demand.

To deal with the above problem, the H ∞ observer for Lipschitz nonlinear systems is a potential candidate since the nonlinear term in the ER model satisfies the Lipschitz condition. Over the years, many LMI-based observer designs for Lipschitz nonlinear systems were widely investigated in [START_REF] Rajamani | Observers for lipschitz nonlinear systems[END_REF], [START_REF] Zemouche | On lmi conditions to design observers for lipschitz nonlinear systems[END_REF], [START_REF] Pertew | H/sub/spl infin//observer design for lipschitz nonlinear systems[END_REF], [START_REF] Abbaszadeh | Robust h observer design for a class of nonlinear uncertain systems via convex optimization[END_REF], [START_REF] Darouach | H observers design for a class of nonlinear singular systems[END_REF], [START_REF] Koenig | Observers design for unknown input nonlinear descriptor systems via convex optimization[END_REF]). An interesting solution of designing observer for Lipschitz nonlinear system in the absence of unknown inputs was presented by [START_REF] Phanomchoeng | Observer design for lipschitz nonlinear systems using riccati equations[END_REF]), based on the Sprocedure lemma. Here, we aim to extend the approach to the Lipschitz system in the presence of sensor noises and unknown input disturbances and apply it to estimate the damping force of a semi-active suspension system.

In this paper, an H ∞ observer using two accelerometers is proposed in order to estimate the damper force in the presence of unknown road input and measurement noises. The design of the observer is based on a nonlinear suspension model consisting of a quarter-car vehicle model, augmented with a first order dynamical nonlinear damper model. Such a model captures the main behaviour of the ER dampers in an automotive applications. The major contribution of this paper are as follows:

• An H ∞ approach for Lipchitz nonlinear system is developed to design an observer minimizing, in an L 2induced gain objective, the effect of unknown inputs (road profile and measurement noises). • The proposed observer has been implemented on a real scaled-vehicle test bench, through the Matlab/Simulink real-time workshop. The observer performances are then assessed with experimental tests

The remainder of this paper is organized as follows:

• Section 2: Semi-active suspension modeling.

• Section 3: Observer design.

• Section 4: Analysis of the observer design: frequency and time domains • Section 5: Experimental validation • Section 6: Conclusion

SEMI-ACTIVE SUSPENSION MODELING AND QUARTER-CAR SYSTEM DESCRIPTION

Semi-active suspension modeling

In the sequel, the dynamic nonlinear model for semi-active ER Damper is described. According to [START_REF] Guo | Dynamic modeling of magnetorheological damper behaviors[END_REF], a phenomenological model of semi-active suspension can be represented by the following nonlinear equation:

F d = k 0 x d + c 0 ẋd + f c .u.tanh(k 1 x d + c 1 ẋd ) (1)
where F d is the damper force; c 0 , c 1 , k 0 , k 1 , f c are constant parameters. x d , ẋd are deflection and deflection velocity of the damper, respectively. u is control input (in the available test bench, this is the duty cycle of the PWM signal that controls the application).

Fig. 1. 1/4 car model with semi-active suspension

From (1), the damper force F d is separated into the passive and controlled parts as follows:

F d = k 0 x d + c 0 ẋd + F nl (x, u) F nl (x, u) = f c .u.tanh(k 1 x d + c 1 ẋd ) (2)
Now in order to take into account the dynamical behavior of the ER fluid, it is important to complete the above model by including a first-order dynamical equation in the controlled part F nl (x, u):

τ Ḟer + F er = F nl (x, u) (3)
Therefore, the complete nonlinear damper dynamical model is given as

   F d = k 0 x d + c 0 ẋd + F er Ḟer = - 1 τ F er + f c τ Φ(x, u) (4)
where Φ(x, u) = u.tanh(k

1 x d + c 1 ẋd )
It is noted that linear and nonlinear indentification methodologies are used to determine all the parameters of the above model (shown in table 1). They are not described here since it is out of the scope of this paper.

Quarter-car system description

This section introduces the quarter-car model with the semi-active ER suspension system depicted in Fig. 1. The well-known model consists of the sprung mass (m s ), the unsprung mass (m us ), the suspension components located between (m s ) and (m us ) and the tire which is modelled as a spring with stiffness k t . From second law of Newton for motion, the system dynamics around the equilibrium are given as:

m s zs = -F s -F d m us zus = F s + F d -F t (5)
where

F s = k s (z s -z us ) is the spring force, F t = k t (z us - z r
) is the tire force, and the damper force F d is given as in (4) with deflection

x d = z def = z s -z us .
Fig. 2. Schematic diagram of testbed of semi-active suspension system Subtituting (4) into (5), one easily obtains

               zs = - 1 m s [(k s + k 0 )(z s -z us ) + c 0 ( żs -żus ) + F er ] zus = 1 m us [(k s + k 0 )(z s -z us ) + c 0 ( żs -żus ) + F er -k t (z us -z r )] Ḟer = - 1 τ F er + f c τ Φ(x, u) (6
) where z s and z us are the displacements of the sprung and unsprung masses, respectively; z r is the road displacement input. By selecting the system states as x = [x 1 , x 2 , x 3 , x 4 , x 5 ] T = [z s -z us , żs , z us -z r , żus , F er ] T ∈ R 5 and the measured variables y = [z s , zus ] T ∈ R 2 , the system dynamics in the state-space representation can be written as follows

ẋ = Ax + BΦ(x, u) + D 1 ω y = Cx + D 2 ω (7)
where

A =           0 1 0 -1 0 - (k s + k 0 ) m s - c 0 m s 0 c 0 m s - 1 m s 0 0 0 1 0 (k s + k 0 ) m us c 0 m us - k t m us - c 0 m us 1 m us 0 0 0 0 - 1 τ           C =    - (k s + k 0 ) m s - c 0 m s 0 c 0 m s - 1 m s (k s + k 0 ) m us c 0 m us - k t m us - c 0 m us 1 m us    B =       0 0 0 0 f c τ       , D 1 =      0 0 0 0 -1 0 0 0 0 0     
, D 2 = 0 0.01 0 0.01 ω = żr n , in which, żr is the road profile derivative and n is the sensor noises.

The control input function Φ(x, u) of the system (7) can be rewritten under the following form

Φ(x, u) = u.tanh      [k 1 c 1 0 -c 1 0]      x 1 x 2 x 3 x 4 x 5           = u.tanh(Γx) (8)
where

Γ = [k 1 , c 1 , 0, -c 1 , 0] Therefore, Φ(x, u) satisfies the Lipschitz condition in x Φ(x, u) -Φ(x, u) Γ(x -x) , ∀x, x (9)
Note that the measured outputs y = [z s , zus ] T can be obtained easily from on board sensors (accelerometers). The H ∞ observer for the quarter-car system ( 7) is defined as follows ẋ = Ax + L(y -C x) + BΦ(x, u) (10) where x is the estimated states. The observer gain L will be determined in the next steps

The estimation error is given as

e(t) = x(t) -x(t) (11)
Differentiating e(t) with respect to time and using ( 7) and ( 10), one obtains

ė = ẋ - ẋ = Ax + BΦ(x, u) + D 1 ω -Ax -L(y -C x) -BΦ(x, u) = (A -LC)e + B(Φ(x, u) -Φ(x, u)) + (D 1 -LD 2 )ω (12)
Assuming the Lipschitz condition (9) for Φ(x, u), the H ∞ observer design objective is stated below

• The system (12) is stable for ω(t) = 0 • e(t) L2 < γ ω(t) L2 for ω(t) = 0
The following theorem solves the above problem into an LMI framwork. Theorem 1. Consider the system model ( 7) and the observer (10) . Given positive scalars γ and l . The system (12) is asymptotically stable for ω = 0 and

e(t) L 2 ω(t) L 2
< γ for ω(t) = 0 if there exist a symmetric positive definite matrix P and a matrix Y satisfying

  Ω P B P D 1 + Y D 2 * -l I d 0 n,d * * -γ 2 I   < 0 ( 13
)
where

Ω = A T P + P A + Y C + C T Y T + l Γ T Γ + I n
The observer matrix will be then L = -P -1 Y Proof. Consider the following Lyapunov function

V (t) = e(t) T P e(t) (14)
Differentiating V (t) along the solution of (12) yields V (t) = ė(t) T P e(t) + e(t) T P ė(t) = [(A -LC)e + B(Φ(x, u) -Φ(x, u))

+ (D 1 -LD 2 )ω] T P e + e T P [(A -LC)e + B(Φ(x, u) -Φ(x, u)) + (D 1 -LD 2 )ω] =   e T (Φ(x, u) -Φ(x, u)) T ω T   T ×   (A -LC) T P + P (A -LC) P B P (D 1 -LD 2 ) B T P 0 0 (D 1 -LD 2 ) T P 0 0   × e Φ(x, u) -Φ(x, u) ω (15)
For brevity, define η = e Φ(x, u) -Φ(x, u) ω

, then one obtains

V (t) = η T M η (16) where M =   Ω 1 P B P (D 1 -LD 2 ) B T P 0 0 (D 1 -LD 2 ) T P 0 0   with Ω 1 = (A -LC) T P + P (A -LC)
From ( 9), the following condition is obtained

(Φ(x, u) -Φ(x, u)) T (Φ(x, u) -Φ(x, u)) e T Γ T Γe ⇔ e T (Φ(x, u) -Φ(x, u)) T T × -Γ T Γ 0 0 I e (Φ(x, u) -Φ(x, u)) 0 ⇔η T Qη 0 ( 17
)
where

Q =   -Γ T Γ 0 0 0 I 0 0 0 0  
In order to satisfy the objective design w.r.t. the L 2 gain disturbance attenuation, the H ∞ performance index is defined as:

J = e T e -γ 2 ω T ω = e T ω T T I 0 0 -γ 2 I e ω = η T Rη (18) where R =   I 0 0 0 0 0 0 0 -γ 2 I  
By applying the S-procedure [START_REF] Boyd | Linear matrix inequalities in system and control theory[END_REF] to the two contraints and J 0, V (t) < 0 if exists scalar l 0 such that V (t)l (η T Qη) + J < 0

⇔η T M -l (η T Qη) + η Rη < 0 ⇔η T (M -l Q + R)η < 0 (19)
The condition ( 19) is equivalent to

M -l Q + R < 0 ⇔   Ω 1 + l Γ T Γ + I P B P (D 1 -LD 2 ) B T P -l I 0 (D 1 -LD 2 ) T P 0 -γ 2 I   < 0 (20)
Let define Y = -P L and substitute into (20), the LMI ( 13) is obtained.

-140 -120

-100 In this section, the synthesis result of the H ∞ observer is presented and some simulation results are given.

Synthesis results and frequency domain analysis

Solving theorem 1 with l = 2, we obtain the L 2 gain γ = 1.1938 and the observer gain

L =      -0.
    
The resulting attenuation of the sensor noises and road profile disturbance on the estimation error, subject to the minimization problems, is shown in Figure 4. These results emphasize the attenuation level of measurement noises and unknown road profile effect on the 5 estimation errors, since the largest sensor noise and road profile disturbance amplification of the 5 errors, over the whole frequency range, are -30dB and -77dB, respectively.

Simulation

To demonstrate the effectiveness of the proposed design, the simulations are made with the nonlinear quarter-car model ( 7). The block-scheme given in Figure 3 illustrates how the simulations of the H ∞ observer are done. The following initial conditions of the proposed observer are considered:

x 0 = [0 0 0 0 0] T x0 = [0.01, -0.4, 0.001, -0.15, 2]

T Two simulation scenarios are used to evaluate the performance of the observer as follows:

Scenario 1:

• The road profile is a sequence of sinusoidal bumps z r = 15sin(4πt)(mm). • The duty cycle u = 0.1 is chosen Scenario 2:

• An ISO 8608 road profile signal (Type C) is used.

• The duty cycle u = 0.1

The simulation results of two tests are shown in the Fig. 5 and Fig. 6. It can be clearly observed in Fig. 5 and Fig. 6 that the damping force is estimated with a satisfactory accuracy.

Robustness analysis

In this section, the µ-tool is used to analyse the robustness of the H ∞ observer. The estimation error dynamics (12), including the system model ( 7) and the observer (10), is analysed here. The system ( 23) is written into the LFT representations (see Figure 7 and Figure 8). In which, Σ(s) is the system (23) with nominal values of parameters, shown in the Table 2. ∆ r (s) represents the parametric uncertainties shown in the Table 2. ∆ L (s) is the uncertain Lipschitz condition (22). The performance objectives are represented by the fictitious uncertainties ∆ f (s). Both robust stability and performance analysis are based on the LFT representations of the uncertain systems (shown in Figure 7 and Figure 8, respectively). These results are presented in the Fig. 9, where the upper and lower bounds of µ are always less than 1. Therefore, the Fig. 8. Uncertain system for robust performance analysis robust stability and robust performance of the observer system are guaranteed. Note that the observer is stable even for larger uncertainties. 2). The proposed observer system is implemented on the host PC using Matlab/simulink with the sampling time T s = 0.05s. Note that the experimetal platform is fully equipped sensors to measure its vertical motion. Each corner of the system has a DC motor to generate the road profile.

In this study, the damping force estimation algorithm is applied for the rear-left corner whose available sensors are the unsprung mass accelerometer zus , the sprung mass accelerometer zs , the damping force sensor F d , and the position sensors to measure suspension deflection z def , road profile z r and unsprung mass position z us . As previously mentioned, only both unsprung mass acceleration zus and sprung mass acceleration zs are used as inputs of the proposed observer. The following block-scheme illustrates the experiment scenario of the observer (shown in Fig. 11)

In this experiment senario, the duty cycle of PWM signal is constant u = 0.1 and the real road profiles are sequence of sinusoidal bumps and ISO 8608 road, shown in Fig. 12( This paper presented an LMI-based H ∞ observer to estimate the damper force, based on the dynamic nonlinear model of the ER damper. For this purpose, the quartercar system is represented by considering a phenomenological model of damper. Using two accelerometers, an H ∞ observer is designed, providing a good estimation result of the damping force. The estimation error is minimized accounting for the effect of unknown inputs (road profile disturbance and measurement noises) and the nonlinearity term bounded by a Lipchitz condition. The robust statbiliy and robust performance properties are ensured by using the µ tool. Both simulation and experiment results assess the ability and the accuracy of the proposed models to estimate the damping force of the ER semi-active damper.

Fig. 4 .

 4 Fig. 4. Transfer e/ω -Bode diagrams

Fig. 6 .

 6 Fig. 5. Simulation scenario 1: Estimated force v.s real force

Fig. 7 .

 7 Fig. 7. Uncertain system for robust stability analysisThe two main types of uncertainties are the Lipschitz nonlinear matrix Γ and the time-varying parametric uncertainty, concerning the unsprung mass, the sprung mass, the spring stiffness and the damper model coefficients. In particular, the parametric uncertainties are shown in the Table2and a Lipschitz nonlinear uncertainty is assumed as followsΦ ∆ (x, u) = Φ(x, u) + ∆Φ(x, u) (21) where ∆Φ(x, u) -∆Φ(x, u) ∆Γ(x -x) (22) with ∆Γ = Γ ∆ -ΓThe estimation error dynamics (12) in the uncertain domain is rewritten as followsė∆ = (A ∆ -LC ∆ )e + B ∆ (Φ ∆ (x, u) -Φ ∆ (x, u)) + (D 1∆ -LD 2∆ )ω(23) where A ∆ , B ∆ , C ∆ , D 1∆ , D 2∆ are the system matrices considering the uncertainties.

Fig. 9 .

 9 Fig. 9. (a) Robust stability analysis, (b) Robust performance analysis

Fig. 10 .

 10 Fig. 10. The experimental testbed INOVE at GIPSA-lab (see www.gipsa-lab.fr/projet/inove)

  a) and Fig. 13(a), respectively. The experiment results of the observer are presented in Fig. 12 (b) and Fig. 13 (b). The result illustrates the accuracy and efficiency of the proposed observer. To further describe this accuracy, Table 3 presents the normalized root-mean-square errors, considering the difference between the estimated and measured forces, for the simulation and experimental results presented in the Fig. 5, Fig. 12 (b) and Fig. 13 (b).

Fig. 12 .

 12 Fig. 12. Experiment test 1: (a) Road profile, (b) Real force vs. Estimated force

Table 1 .

 1 Parameter values of the quarter-car model equipped with an ER damper

	Parameter	Description	value Unit
	ms	Sprung mass	2.27	kg
	mus	unsprung mass	0.25	kg
	ks	Spring stiffness	1396 N/m
	kt	Tire stiffness	12270 N/m
	k 0	Passive damper stiffness coefficient	170.4 N/m
	c 0	Viscous damping coefficient	68.83 N.s/m
	k 1	Hysteresis coefficient due to displacement 218.16 N.s/m
	c 1	Hysteresis coefficient due to velocity	21 N.s/m
	fc	Dynamic yield force of ER fluid	28.07	N
	τ	Time constant	43	ms
		3. OBSERVER DESIGN
	In this section, an H

∞ observer is developed to estimate the damping force accurately. The unknown input ω (road profile disturbance and measurement noise) is considered as an unknown disturbance. Therefore, an H ∞ observer is proposed to minimize the effect of the accounting for unknown disturbance ω on the state estimation errors and to bound the nonlinearity by Lipschitz constant.

Table 2 .

 2 Parameter values of the quarter-car model equipped with an ER damper

	Uncertain parameters	Variation
	ms	2.27 ∓ 50%(kg)
	mus	0.25∓ 20%(kg)
	ks	1396∓ 20%(N/m)
	kt	12270 ∓ 20%(N/m)
	k 0	170.4∓ 20%(N/m)
	c 0	68.83∓ 20%(N.s/m)
	k 1	218.16∓ 20%(N.s/m)
	c 1	21∓ 20%(N.s/m)
	fc	28.07∓ 20%
	τ	43∓ 20%(ms)

Table 3 .

 3 Normalized Root-Mean-Square Errors

	Road Profile Sequence	Simulation Experiment
	Sinusoidal bumps	0.0234	0.0896
	ISO 8608 road	0.0366	0.1015
	This test-bench which involves 4 semi-active ER suspen-
	sions is controlled in real-time using xPC target and a
	host computer. The target PC is connected to the host
	computer via Ethernet communication standard (see Fig-
	ure

This work has been partially supported by the 911 scholarship from Vietnamese government. The authors also thank for financing project EMPHYSIS.