Hong Van Truong
email: hongvantruongiph@gmail.com

Elise Beck
email: elise.beck@ujf-grenoble.fr

Julie Dugdale
email: julie.dugdale@imag.fr

Carole Adam
email: carole.adam@imag.fr

Developing a model of evacuation after an earthquake in Lebanon

Keywords: Agent based simulation, pedestrian mobility, human behavior modelling, GAMA platform

This article describes the development of an agent-based model (AMEL, Agent-based Model for Earthquake evacuation in Lebanon) that aims at simulating the movement of pedestrians shortly after an earthquake. The GAMA platform was chosen to implement the model. AMEL is applied to a real case study, a district of the city of Beirut, Lebanon, which potentially could be stricken by a M7 earthquake. The objective of the model is to reproduce real life mobility behaviours that have been gathered through a survey in Beirut and to test different future scenarios, which may help the local authorities to target information campaigns.

INTRODUCTION

Every year earthquakes cause many victims. Although the survival of people is largely related to the resilience of buildings, the way that people behave after an earthquake also influences the number of fatalities (Rojo, Beck, Lutoff and Schoeneich, 2011), especially for low magnitude earthquakes. The rareness of strong earthquakes and the impossibility of forecasting them prevent us from assessing the effect that behaviours may have on survival rates. In this context, computer simulation is an appropriate and powerful tool that helps us to assess different scenarios, thus helping local authorities to develop their risk management policies and information campaigns.

The objective of the AMEL model (Agent-based Model for Earthquake evacuation in Lebanon) is to simulate pedestrians' mobility shortly after an earthquake. Although many victims in buildings are generally killed during the main shock, the mobility of any survivors can greatly affect their chances of surviving any aftershocks. We want to understand whether a better knowledge of the safety procedures and the location of safe places in a city could decrease people's vulnerability. Firstly of all, the model aims at reproducing the behaviours adopted by individuals. Secondly, the simulation tests some optimistic or pessimistic fictive scenarios, in order to evaluate the impact of individual behaviours on the number of fatalities or on the number of people in exposed and dangerous situations. This work is applied to the real case of a district of Beirut, Lebanon. This article presents a work in progress and focuses on explaining the model and how it was developed. Section 2 presents the overall approach and gives the state of the art, whilst section 3 describes the methodology used for developing AMEL. The model is described in section 4 and our first attempts at validation are detailed in section 5. Section 6 concludes the paper and presents future work.

Developing a model of evacuation after an earthquake in Lebanon

Proceedings of the ISCRAM Vietnam 2013 Conference -Institut de la Francophonie pour l'informatique (IFI), Ha Noi, Vietnam, October 30 th -November 1 st 2013 Serge Stinckwich, Ho Tuong Vinh, Julie Dugdale and Nguyen Hong Phuong, eds.

APPROACH AND STATE OF THE ART

Agent-based social simulation

Agent-based social simulation (ABSS) is a branch of distributed artificial intelligence and multi-agent systems whose focus is on understanding, modelling and simulating social behaviours. Rather than purely focusing on cognition, this approach recognises the social complexity of a society and looks at how social phenomena, such as cooperation, emerge through human behaviours.

The use of agents in modelling human systems has several advantages over other approaches [START_REF] Bonabeau | Adaptive Agents, Intelligence, and Emergent Human Organization: Capturing Complexity through Agent-Based Modeling[END_REF]. Firstly, agent based systems are able to capture emergent phenomena that are so representative of complex adaptive systems. Secondly, they provide a natural description of a system, which as Bonabeau notes makes the agent based approach much closer to reality. Finally, they are flexible, allowing us to study social systems at different levels of abstraction by varying the complexity of our agents or by aggregating agents into subgroups.

Hundreds of agent-based social simulators have been developed. These have been used for predicting future situations, as training tools, for developing and formalizing theories, or for testing new technological designs or new ways of organizing work [START_REF] Gilbert | Simulation for the Social Scientist[END_REF]). Here we focus on the use of a simulator for understanding and predicting human behaviour when faced with an earthquake, depending on their level of vulnerability.

Using agent based simulation as a way to assess emergency situations has become increasingly popular in recent years. The RoboCup Rescue Agents Simulation project has attracted many researchers who are interested in using the specially developed platform to develop intelligent agents that undertake the role of Police Forces, Fire Brigades, and Ambulance Teams in a disaster response scenario [START_REF] Skinner | The RoboCup Rescue simulation platform[END_REF]. Likewise, other works, such as the REScUE research project (Hawe, Coates, Wilson and Crouch, 2011) and EQ-Rescue [START_REF] Fiedrich | An HLA-based multiagent system for optimized resource allocation after strong earthquakes[END_REF] have been developed to evaluate different rescue plans, or to optimize resource allocation during an emergency. These works, like many of this kind of simulation, focus specifically on the response and rescue aspects, with the main accent being on modelling rescue agents, on coordination between different teams, or on the allocation of resources. Our work differs in this respect since we do not look at response activities, but just the behavior of the potential victims.

Pelechano and Badler present a model of building evacuation to study the influence of inter-agent communication and the effect of training [START_REF] Pelechano | Modeling Crowd and Trained Leader Behavior during Building Evacuation[END_REF]. The movements of agents in rooms are modelled by acceleration equations. The authors show that the evacuation may be made more efficient by introducing communication between agents to share their knowledge of blocked routes. Similarly the evacuation is also more efficient if there are a small number of leaders (trained staff) in the crowd that the other dependent agents group around and follow. Although we have designed AMEL to take into account leader/follower behaviours, our work is still in progress and therefore we have not yet experimented with this aspect. Concerning inter-agent communication, we hope to build on this work by considering neighbourhood evacuation.

Nguyen et al. present a hybrid model of the pedestrian flow on road networks, applied to the evacuation of Nhatrang (Vietnam) in case of a tsunami (Nguyen, Zucker, Nguyen, 2011). This model combines micro (agentbased, fine-grained but slow) and macro (equation-based, fast) models of pedestrians' movement, in order to improve the efficiency (speed) of simulations involving a large number of agents. This hybrid model was shown to be more efficient than a micro-model, and of better quality than a macro-model.

METHODOLOGY

A controversy in ABSS concerns the fundamental question of how to develop useful models of real-life social situations. Broadly there are two schools of thought. One follows the KISS (Keep It Simple, Stupid) philosophy where the aim is to develop simplistic models and where much of the real world detail has been abstracted away. Although there are obvious benefits, e.g. in terms of ease of constructing the models, the approach has been widely criticised. The arguments can be reduced to the idea that models that are too simple only address simplistic problems that are not representative of the real world. The other extreme is a KIDS (Keep It Descriptive Stupid) approach [START_REF] Edmonds | From KISS to KIDS: An "anti-simplistic" modelling approach[END_REF] where the model is constructed by taking into account the widest possible range of evidence, including anecdotal accounts and expert opinion. Although we obtain a much truer representation of reality it may be very difficult to obtain the data to build the model, implementation is more complicated, and validation of the model and simulator are problematic. The approach adopted in this work falls in between these two extremes and follows that proposed by Rosaria Conte: "Keep it Simple as Suitable" [START_REF] Conte | The necessity of intelligent agents in social simulation[END_REF]. Here models are abstract enough to achieve an adequate level of generality, but no less complex than what is required by the purpose of the simulation. Whilst the original KISS, KIDS, and reformulated KISS approaches provide advice on designing models it is very general and somewhat vague, and the approaches lack a complete modelling method. In response, several methods and modelling techniques have been proposed, for example GAIA [START_REF] Wooldridge | The GAIA Methodology for Agent-Oriented Analysis and Design[END_REF], VOWELS [START_REF] Demazeau | From interactions to collective behaviour in agent based-systems[END_REF], CoMoMAS [START_REF] Glaser | Contribution to Knowledge Acquisition and Modelling in a Multi-Agent Framework -The CoMoMAS Approach[END_REF], MMTS (Kinny, Geogeff and Rao, 1996), and Unified Approach (Sabas, Delisle and Badr, 1996). These all provide the standard framework for modelling the agent dimension, some taking into account the deliberative behaviour of agents. However, they are largely intended for developing general multiagent systems and are not specifically focused on modelling the social elements that are required in ABSS. Furthermore, they fail to provide a structure for analysing human agents in the design phase and for validating the model with respect to the observed human behaviour. Some methods have been developed that focus particular on agent based simulation. Notably the ODD protocol by Grimm and his colleagues aimed to standardising descriptions of models to aid understanding and ensure reproducibility (Grimm, Berger, Bastiansen, Eliassen, Ginot, Giske, Goss-Custard, Grand, Heinz, Huse, Huth, Jepsen, Jørgensen, Mooij, Müller, Pe'er, Piou, Railsback, Robbins, Robbins, Rossmanith, Rüger, Strand, Souissi, Stillman, Vabø, Visser and DeAngelis DL. 2006; Grimm, Berger, DeAngelis, Polhill, Giske and Railsback, 2010). Although very useful, the protocol only concerns model description and does not suggest how the model itself may be developed.

The methodology adopted in this work is shown in figure 1 and is adapted from that of Edmonds (Edmonds, 2000). However it has several important differences. Firstly, it focuses on analysing human behaviour in the real world situation through the use of extensive field studies. This provides a solid corpus of empirical data through video recordings and observations, etc. Secondly, it puts validation at the heart of the process ensuring that the results of the simulator can be more readily trusted (Dugdale, Bellamine-Ben Saoud, Pavard, Pallamin, 2001). Finally, it reinforces iteration; this allows us to revisit previous steps such as undertaking additional targeted field studies and refining the model and code as necessary.

FIELD STUDIES AND DATA COLLECTION

The studied district is characterized by a high density of buildings that have a differing number of storeys, date and construction material. Several types of data were collected using a multidisciplinary approach, involving the disciplines of geography, psychology and geotechnics. Concerning the geographic aspect, the buildings, streets, and green spaces correspond to spatial data represented by polygons and polylines. These objects correspond to geographic layers and were processed through geographic information systems. The data first comes from a national reference database (such as the French Geographic Institute database) and was corrected and updated using high spatial resolution satellite imaging. Concerning demographic data, the last Lebanon population census was carried out in 1932. Therefore, the population was estimated at 3.8 persons per apartment.

Because of the high heterogeneity of the building types, which leads to high variations of physical vulnerability to earthquakes, each of the 357 buildings was documented through a field survey that aimed at characterizing its vulnerability. For example, the survey form included information such as number of storeys, year of construction, construction material. Some geotechnical data was collected with material samples for different types of construction. This data was then computed with the FEMA's method called HAZUS (Kircher, Whitman and Holmes, 2006) in order to estimate the damage rate for each building in case of an M7 earthquake on the Yammouneh fault.

Finally, the model requires social data concerning the behaviours of individuals. This data was collected through a field sociological survey that interviewed 88 persons of the studied district (Beck, Colbeau-Justin, Cartier and Saikali, 2011). The questionnaire focused on several subjects including knowledge about earthquakes, risk perception, earthquake experience and associated behaviours, etc. The statistical analysis of the survey allowed us to define several categories of behaviours. Some extra behaviours, that were not observed in the survey but have been reported in similar cases, were also considered and implemented into the model (leader-follower behaviour, for example).

MODEL DESIGN

As a first step, the purpose of the simulation is to reproduce, based on the survey information, the mobility of people and their behaviours following an earthquake in Lebanon. In an earthquake, natural and artificial obstacles (e.g. escarpments and stairways), and obstacles induced by the earthquake itself, can injure people and constrain their movement towards safe zones. The simulator will provide information concerning the physical and human damage incurred (the number of damaged buildings, blocked streets, number of victims, number of people in dangerous zones, evacuation time, etc.).

Once we have reproduced mobility behaviours that align to the survey information, as a second step, we experiment with difference scenarios by changing the ways that people behave and analysing the results. We want to see what happens (damages, number of fatalities and exposure to danger) if people act in different ways. These results can serve to inform people about how to improve their behaviours, for example by forming the basis of an information campaign.

Entities and their attributes

The model consists of six entities: Human, Street, Building, Green space, Quake, and Obstacle (figure 2). Their attributes are divided into three types: position and visualization (location, colour, etc.), entity characteristics (e.g. Human entity has attributes about age, sex, profession, etc., and the Building entity has attributes about height, capacity, etc.). Finally there are special attributes concerning the entity's behaviours (e.g. target to reach, street_knowledge, etc. of the Human entity).

Human agent

Since the Human agent is the most important entity we focus on describing this agent in detail. Human agents represent people in the district. After an earthquake, agents can decide to move, on roads that are not blocked, or stay at their current location. They can become informal leaders who know the safe areas and can guide other agents to their targets. Agents can also become followers by searching for a leader in their perception zones and then following that agent. If a follower cannot find a leader, it wanders around.

Human agents can observe obstacles and other agents inside its perception radius (attribute radius_perception). After meeting a big obstacle, the Human agent can act in different ways. If it is a leader, it will choose another way to reach target. If it is a wanderer, it will choose another target to go to. Each Agent update their street knowledge, remembering which streets are blocked so that they may be avoided in future.

Obstacle agent

Obstacle agents can affect the behaviour of Human agents. In this work, we focus on the obstacles induced by an earthquake itself, particularly by the buildings damaged by the earthquake. Based on the damage level of buildings, we can group obstacles into three levels: big, medium and small. The big obstacles destroy all the streets and people in those zones, so no agent can pass these obstacles. In the zone of medium and small obstacles, people are exposed to the danger (attribute exposed); this means that the agent is in the danger zone and there is a high risk of being injured. Medium obstacles can constrain the movement of agents in their zone, making them move more slowly.

Other agents

The Quake agent has attributes for the intensity and time of the earthquake. Building agents have attributes for damage level of the building. Street agents are considered as a weighted graph where each street is an edge and each intersection is a vertex. Each unaffected street has a small weight (e.g. the length of street), but blocked streets have a much higher weight than normal streets (e.g. a billion). By using this weighted graph approach each agent has their own representation (beliefs) about the state of the streets.

The spatial scale in our model is in meters and each simulation cycle lasts one second. The simulation begins after all human agents have decided upon their target, and ends when all the leaders have reached their targets.

Agent behaviours

Human behaviours:

How people behave in reality is very complex. In our model, we just focus on the mobility of people and do not consider more complex behaviours. From the survey, we synthesize 6 behaviours:

1. Move without changing activity. This covers people that had intended to move and do not change their activity as a result of the earthquake (e.g. people who are on their way to work) 2. Stay in place without changing activity. This covers people who remain where they are and continue doing the same activity (e.g. those who are at home and do not have employment). 3. Change activity in order to go to a safe place. This covers people who intentionally move to a safe place just after an earthquake. 4. Change activity in order to go to an unsafe place. This covers peoples who intentionally move to an unsafe place (e.g. someone who goes to a relative's home even though it may be unsafe). 5. Change activity and stay in place. This covers people who had planned to go somewhere, but because of the earthquake they remain where they are. 6. Imitate others. This addresses leader-follower behaviours.

The simulation begins when all of the agents have determined where to go. Thus, from a mobility point of view, behaviours 1, 3 and 4 are equivalent (the agent moves), and behaviours 2 and 5 are also equivalent (the agent does not move).

Human behaviours are synthesized into three groups: movement, obstacle perception and imitation. While movement behaviour concerns the mobility of people, obstacle perception and imitation relate to determining and updating target locations and knowledge of the blocked streets. The streets are considered as a weighted graph on which the Human agent moves and determines the shortest path to reach its target. Moreover, the street knowledge of agent is considered as list of weights applied for the graph. Each agent has a different knowledge of the streets, so the weights applied for the graph are different according to each agent. Thus the shortest path to the target depends on two factors: the agent's target and the agent's street knowledge (figure 3).

• Movement behaviour

• Obstacle perception behaviour

Agents within the zone of medium and small obstacles are considered to be exposed to danger. In addition, medium obstacles constrain the movement of agents, making them move along their route more slowly.

Agents can observe obstacles in their zone of perception and when they encounter a big obstacle they cannot proceed. However, they remember blocked streets so that they may be avoided in the future when finding another path to the target. The difference between the behaviour of leaders and followers is that if they meet an obstacle, leaders will find other ways to reach their target, while followers will find another target to go to. In our model, we consider 5 blocked streets configurations (Fig. 4). In figure 4 the red circle is an obstacle and a blue square is the Human agent. In (1) the Agent is in the area with no exit. In (2) there is a long street section before the agent moves towards the obstacle. In (3) there is a short street section before the agent moves towards the obstacle. In (4) the Agent is in the street with two obstacles on both sides. In [START_REF] Demazeau | From interactions to collective behaviour in agent based-systems[END_REF] the Agent sees an obstacle on another street.

The simplest method to implement this behaviour is to update an agent's knowledge of blocked streets and target path when the agent encounters an obstacle. However, this leads to several problems. First, with configuration (2), the agent is still in the blocked street, it must go back to choose another street, but the street section to go back to is longer than the rest of the street, so the agent often chooses to enter into the obstacle. To overcome this problem we make all agents update their street knowledge and go back to the previous intersection after meeting an obstacle (Figure 5). After that, the weights of blocked streets are updated, leader agents build other paths to reach their targets, and wandering agents choose other targets and build paths to go. If an agent encounters an obstacle as in the configuration (5), figure 4, then it will normally update the weights.

Figure 5. Agent moves to the previous intersection to update its street knowledge

Secondly, with configurations ( 1) and ( 4), the agent has no choice but to pass over obstacles to get out. In order to fix this problem we add a parameter concerning the maximum number of times that the agent encounters an obstacle (2-3 times). After this, the agent stops and stays where it is. This feature relates to the psychological aspect of people in the evacuation. If a person tries to find a way to its target but repeatedly doesn't success it will eventually give up. Figure 6 shows the obstacle perception behaviour.

• Imitation behaviour

For the followers and wanderers, at the beginning, they search randomly around them for a leader. If they find a leader, then they become followers. They update their targets to the position of leader and update paths to follow it. Figure 7 shows the process of imitation. slowest follower (figure 8). The Obstacle agent has two behaviours: street damage and human damage. With street damage, each agent finds all the streets around it that it overlaps and updates the "blocked" attribute of that Street agent to "true". With human damage, each Obstacle agent finds all of the Human agents with which it overlaps and destroys them.

DEVELOPMENT

The AMEL model has been implemented on the GAMA platform (Grignard, Taillandier, Gaudou, Huynh and Drogoul, 2013); GAMA is an open-source generic agent-based modelling and simulation platform. It provides an intuitive modelling language with high-level primitives to define agents and their environment. GAMA includes a powerful Integrated Development Environment to help non-computer scientists to develop complex models with powerful features in terms of Geographical Information Systems (GIS) integration and high-level tools (e.g. decision-making or clustering algorithms). In addition, both the language and the software have been designed to allow the development of big models with a huge number of agents (with various architectures from reflex-based to belief-desire-intention architectures). GAMA also allows modellers to manage various kinds of complex environments such as square, hexagonal or irregular grids, networks, or continuous environments linked to GIS data.

GAMA enables a hierarchical and dynamic organization of agents. The platform is also easily extensible to add new features to models (e.g. to give the possibility of integrating equation-based (ODE) models into agents) or new agent architectures or features. GAMA has been successfully used to develop various large-scale applications that share the need for a tight integration and management of huge GIS data, and for strong interactions between a complex environment and the agents. It was used for instance in the MAELIA platform for simulating water management problems (Gaudou, Sibertin-Blanc, Thérond, Amblard, Arcangeli, Balestrat, Charron-Moirez, Gondet, Hong, Mayor, Panzoli, Sauvage, Sanchez-Perez, Taillandier, Nguyen, Vavasseur and Mazzega, 2013).

Figure 9 shows a screenshot of the simulation. The main display represents the district with the streets in black, green space in green, and the buildings in yellow. Human agents are represented by small circles of different colours: hot colours (red, pink) for agents who stay where they are; cold colours (blue, purple) for agents who move; and the colour cyan for imitating agents. We also created other displays: charts of the different behaviours; charts about the number of victims, people exposed, leaders and followers, etc; and a chart about the total exposition time. On the right, we can see the parameters of the simulator; these concern the location of the human agents and their behaviours, as well as some global simulator parameters. By varying the value of these parameters we can create different scenarios which may help to see the effect of the different information campaigns. For example, if we want to see the effect on the number of victims if an earthquake occurs during the night, we can set the value for the "people in buildings" parameter (probability of people who are in the buildings) equal to 100 (percent), indicating that at night everyone is inside a building.

VALIDATION

Once the model has been implemented, it must be validated. This is current work in progress. Every newly implemented "brick" of the model is tested individually and all the associated scenarios are evaluated in order to understand the impact of the new brick. As we developed the simulator incrementally it eventually becomes increasingly complex. The challenge is to understand the impact of any new changes on the global model. For example, once human exposure to hazards is implemented, we have to verify that the fatalities figures are realistic. For these kinds of validations we ask the experts who collected the data and those that have real life experience of the consequences of earthquakes. This process of validation involves all other steps of the methodology and is central to the development process. In practice we have found this to be a highly iterative step where we frequently request validation from experts when changing our model and after implementation. This validation constitutes a major step that should not be underestimated.

FIRST RESULTS

Although we are currently in the validation stage, we provide some first results so that the reader may see the type of results that may be achieved with the simulator. The results concern two scenarios: 1. The survey scenario, this is the real situation where most people stay where they are and do not go to the safe place; 2. An optimistic scenario where everyone goes to a safe place (e.g. green space, schools, outside of the city, or even just into the streets).

We can see in the graphs below that the number of exposed people in the optimistic scenario is less than that in the survey scenario. Also the number of 'others' (i.e. those who are not exposed and injured, in other words they are in the safe places) and those who go to outside of the city (considered as a safe place) in the optimistic scenario is more than in the survey scenario.

SUMMARY

Lebanon is an earthquake prone area, lying on the Mount Lebanon Thrust Fold and the Yamouneh faults, which have the capacity to generate a 7.5 magnitude earthquake. Although major earthquakes are rare, minor ones occur frequently and there is a deep concern that a large earthquake could destroy the city, as it did in 551 AD. Despite this threat the study conducted as part of this project showed that, possibly because people are constantly exposed to many small earthquakes, a significant number of them do not change their behaviours during an earthquake and continue to do their originally planned activities. Whilst in small earthquakes this is not problematic, the consequences of being unprepared in the case of a large one could be devastating.

We have developed a model and first prototype of a simulator that shows pedestrians' movements following different human behaviours after an earthquake. Using GAMA we were able to incorporate the actual map of the neigbourhood in Beirut, including current streets and actual buildings. This allowed us to develop a realistic simulator with agents moving from known buildings along current streets. One particular feature of how we have designed the model is that each agent has its own personal view of the state of the road infrastructure. This replicates the real situation where people are only aware of the routes open to them after an earthquake as they move around the neigbourhood towards their target. As agents explore their neigbourhood, their knowledge of the streets is incrementally updated.

We are currently in the process of validating the model and simulator. In particular we are addressing the problem of an agent's exposure to danger and how fatalities may be quantified given the different severities of building damage. Validation is a long process since we need to interact with local experts from different domains and the data required is not always available. Nevertheless, despite some difficulties with validation, the results of the simulator are encouraging.

Figure 1 .

 1 Figure 1. Methodology adopted for AMEL Following figure 1, the first step covers performing detailed field studies of the real situation in order to assess the human behaviours and their underlying motivations. The second step, Model design, involves developing the formal model, for example by using UML, of what has been obtained through the analysis of field studies data. Validation has been put at the heart of the methodology and relies heavily on the data obtained via field studies in the analysis phase. The final step, Simulation, covers experimentation. The bi-directional arrows ensure that iteration plays a major role. The above methodology was first described in 1999 (Dugdale, Pavard and Soubie, 1999; Dugdale, Pavard and Soubie 2000) and has been used for designing and developing agent-based simulators in several works over the years, for example in (Bellamine-Ben Saoud, Ben Mena, Dugdale, Pavard and Ahmed, 2006; Kashif, Binh, Dugdale and Ploix, 2001; Kasif, Ploix, Dugdale and Binh, 2013).

  Developing a model of evacuation after an earthquake in Lebanon Proceedings of the ISCRAM Vietnam 2013 Conference -Institut de la Francophonie pour l'informatique (IFI), Ha Noi, Vietnam, October 30 th -November 1 st 2013 Serge Stinckwich, Ho Tuong Vinh, Julie Dugdale and Nguyen Hong Phuong, eds.

Figure 2 .

 2 Figure 2. Class diagram

  Developing a model of evacuation after an earthquake in LebanonProceedings of the ISCRAM Vietnam 2013 Conference -Institut de la Francophonie pour l'informatique (IFI), Ha Noi, Vietnam, October 30 th -November 1 st 2013 Serge Stinckwich, Ho Tuong Vinh, Julie Dugdale and Nguyen Hong Phuong, eds.

Figure 3 .

 3 Figure 3. The movement process

Figure 4 .

 4 Figure 4. Blocked street configurations (with damages)

  Developing a model of evacuation after an earthquake in Lebanon Proceedings of the ISCRAM Vietnam 2013 Conference -Institut de la Francophonie pour l'informatique (IFI), Ha Noi, Vietnam, October 30 th -November 1 st 2013 Serge Stinckwich, Ho Tuong Vinh, Julie Dugdale and Nguyen Hong Phuong, eds.

Figure 6 .

 6 Figure 6. Obstacle perception process

Figure 7 .

 7 Figure 7. Imitation process for followers and wanderersFor the leaders, every 10 cycles of the simulation, they perceive their followers and adjust their speeds to the

Figure 8 .

 8 Figure 8. Speed adjustment process for leader

  Developing a model of evacuation after an earthquake in Lebanon Proceedings of the ISCRAM Vietnam 2013 Conference -Institut de la Francophonie pour l'informatique (IFI), Ha Noi, Vietnam, October 30 th -November 1 st 2013 Serge Stinckwich, Ho Tuong Vinh, Julie Dugdale and Nguyen Hong Phuong, eds.

Figure 9 .

 9 Figure 9. Screenshot of the simulation

Figure 10 .

 10 Figure 10. Fatality and exposure graphs. On the left, the survey scenario. On the right an optimistic scenarioConcerning the total exposure time, we can see that in the current situation (survey scenario) there are more people who are exposed for a long time, whereas in the optimistic scenario there are more people who are exposed for the short time.

  Developing a model of evacuation after an earthquake in LebanonProceedings of the ISCRAM Vietnam 2013 Conference -Institut de la Francophonie pour l'informatique (IFI), Ha Noi, Vietnam, October 30 th -November 1 st 2013 Serge Stinckwich, Ho Tuong Vinh, Julie Dugdale and Nguyen Hong Phuong, eds.

Figure 11 .

 11 Figure 11. Number of agents exposed. On top, the survey scenario. On the bottom the optimistic scenario. The horizontal axis shows time slot periods and the vertical axis shows the of exposed people in a corresponding time period.

ACKNOWLEDGMENTS

This study is funded by the French National Research Agency (ANR RISKNAT 2009 LIBRIS). The authors would like to thank the different persons involved in the collection of data: Jocelyne Adjizian-Gérard, Rita Zaarour, Nada Saliba, Pierre-Charles Gérard and the students involved in the survey for the building inventory (Geography Department of the Saint Joseph University, Beirut, Lebanon); Jacques Harb for the damages modelling (Notre Dame University, Beirut, Lebanon); Maud Saikali and Carine Azzam (Saint Joseph University, Beirut, Lebanon), Ludvina Colbeau-Justin and Stéphane Cartier (CNRS-PACTE Laboratory, Grenoble, France) for the risk perception survey. The authors also thank Patrick Taillandier (Rouen University, France), Nicolas Marilleau (Institute for Research and Development, France) and Benoit Gaudou (Toulouse University, France) for their precious help on GAMA. Julie Dugdale would like to acknowledge the support of Developing a model of evacuation after an earthquake in Lebanon Proceedings of the ISCRAM Vietnam 2013 Conference -Institut de la Francophonie pour l'informatique (IFI), Ha Noi, Vietnam, October 30 th -November 1 st 2013 Serge Stinckwich, Ho Tuong Vinh, Julie Dugdale and Nguyen Hong Phuong, eds.

the University of Agder, to which she is affiliated.


Z E A Fellah

S Berger

W Lauriks

C Depollier

C Ariste

J.-Y Chapelon

Measuring the porosity and the tortuosity of porous materials via reflected waves at oblique incidence

Keywords: numbers: 43, 20, Bi, 43, 20, Hq ͓ANN͔

An ultrasonic reflectivity method is proposed for measuring porosity and tortuosity of porous materials having a rigid frame. Porosity is the relative fraction by volume of the air contained within a material. Tortuosity is a geometrical parameter which intervenes in the description of the inertial effects between the fluid filled the porous material and its structure at high frequency range. It is generally easy to evaluate the tortuosity from transmitted waves, this is not the case for porosity because of its weak sensitivity in transmitted mode. The proposed method is based on measurement of reflected wave by the first interface of a slab of rigid porous material. This method is obtained from a temporal model of the direct and inverse scattering problems for the propagation of transient ultrasonic waves in a homogeneous isotropic slab of porous material having a rigid frame ͓Z. E.

I. INTRODUCTION

The ultrasonic characterization of porous materials saturated by airlike plastic foams, fibrous, or granular materials is of a great interest for a wide range of industrial applications. These materials are frequently used in the automotive and aeronautics industries or in the building trade. The determination of the properties of a medium from waves that have been reflected by or transmitted through the medium is a classical inverse scattering problem. Such problems are often approached by taking a physical model of the scattering process, generating a synthetic response for some assumed values of the parameters, adjusting these parameters until reasonable agreement is obtained between the synthetic response and the observed data. Two important parameters which appear in theories of sound propagation in porous materials [START_REF] Fellah | Direct and inverse scattering of transient acoustic waves by a slab of rigid porous material[END_REF][START_REF] Fellah | Transient acoustic wave propagation in rigid porous media: A time-domain approach[END_REF][START_REF] Fellah | An approach to direct and inverse time-domain scattering of acoustic waves from rigid porous materials by a fractional calculus based method[END_REF][START_REF] Fellah | Application of fractional calculus to the sound waves propagation in rigid porous materials: Validation via ultrasonic measurements[END_REF][START_REF] Johnson | Theory of dynamic permeability and tortuosity in fluid-saturated porous media[END_REF][START_REF] Allard | Propagation of Sound in Porous Media: Modeling Sound Absorbing Materials ͑Chapman and Hall[END_REF][START_REF] Biot | The theory of propagation of elastic waves in fluid-saturated porous solid. I. Low frequency range[END_REF][START_REF] Biot | The theory of propagation of elastic waves in fluid-saturated porous solid. II. Higher frequency range[END_REF][START_REF] Attenborough | Acoustical characteristics of porous materials[END_REF][START_REF] Beranek | Acoustic impedance of porous materials[END_REF][START_REF] Zwikker | Sound Absorbing Materials ͑Elsevier[END_REF][START_REF] Champoux | New empirical equations for sound propagation in rigid frame fibrous materials[END_REF][START_REF] Lafarge | Dynamic compressibility of air in porous structures at audible frequencies[END_REF] are the porosity and the tortuosity. Porosity is the relative fraction, by volume, of air contained in the connected pores in the material. Unlike other parameters included in the description of different various phenomena occurring in the acoustical propagation of porous material at high frequency range such as tortuosity, [START_REF] Zwikker | Sound Absorbing Materials ͑Elsevier[END_REF] viscous characteristic length, [START_REF] Johnson | Theory of dynamic permeability and tortuosity in fluid-saturated porous media[END_REF] and thermal characteristic length, [START_REF] Champoux | New empirical equations for sound propagation in rigid frame fibrous materials[END_REF] or in the low frequency range such as flow resistivity [START_REF] Zwikker | Sound Absorbing Materials ͑Elsevier[END_REF] and thermal permeability, [START_REF] Lafarge | Dynamic compressibility of air in porous structures at audible frequencies[END_REF] the porosity is a key parameter playing an important role in the propagation at all frequencies. As such, in studies of acoustical properties of porous materials, it is highly desirable to be able to measure this parameter.

Beranek [START_REF] Beranek | Acoustic impedance of porous materials[END_REF] described an apparatus ͑porosimeter͒ used to measure the porosity of porous materials. This device was based on the equation of state for ideal gases at constant temperature ͑i.e., Boyle's law͒. Porosity can be determined by measuring the change in air pressure occurring with a known change in volume of the chamber containing the sample. In the Beranek apparatus, both pressure change and volume change are monitored using a U-shaped fluid-filled manometer. An alternate technique for measuring porosity is a dynamic method proposed by Leonard. [START_REF] Leonard | Simplified porosity measurements[END_REF] Techniques that use water as the pore-filling fluid, rather than air, are common in geophysical studies. [START_REF] Guyon | Transport properties in sintered porous media composed of two particles sizes[END_REF][START_REF] Johnson | Tortuosity and acoustic slow waves[END_REF] Mercury has been used as the pore-filling fluid in other applications. [START_REF] Van Brakel | Mercury porosimetry: State of the art[END_REF] However, for many materials, the introduction of liquids into the material is not appropriate. Recently a similar device to that of Beranek, involving the use of an electronic pressure transducer, was introduced by Champoux et al. [START_REF] Champoux | Air-based system for the measurement of porosity[END_REF] This device can be used to measure very slight changes in pressure accurately, and the output can be recorded by a computer.

The tortuosity ␣ ϱ , namely the structure factor k s by Zwikker and Kosten, [START_REF] Zwikker | Sound Absorbing Materials ͑Elsevier[END_REF] or the parameter q 2 by Attenborough, [START_REF] Attenborough | Acoustical characteristics of porous materials[END_REF] is an important parameter which intervenes in the description of the inertial interaction between the fluid and the structure in the porous material at high frequency range. In the case of cylindrical pores making an angle with the direction of propagation, ␣ ϱ ϭ1/cos 2 . The tortuosity can be evaluated by electrical measurements, [START_REF] Johnson | Tortuosity and acoustic slow waves[END_REF] or by using a superfluid [START_REF] Fellah | Application of fractional calculus to the sound waves propagation in rigid porous materials: Validation via ultrasonic measurements[END_REF] He as the pore fluid. [START_REF] Johnson | Theory of dynamic permeability and tortuosity in fluid-saturated porous media[END_REF] It can also be evaluated by acoustical techniques as an ultrasonics measurement of transmitted waves. [START_REF] Fellah | Direct and inverse scattering of transient acoustic waves by a slab of rigid porous material[END_REF][START_REF] Fellah | An approach to direct and inverse time-domain scattering of acoustic waves from rigid porous materials by a fractional calculus based method[END_REF][START_REF] Fellah | Application of fractional calculus to the sound waves propagation in rigid porous materials: Validation via ultrasonic measurements[END_REF][START_REF] Leclaire | Determination of viscous and thermal characteristics lengths of plastic foams by ultrasonic measurements in helium and air[END_REF] In this work, we present a method of measuring the porosity and the tortuosity from the measurement of acoustic wave reflected by a slab of porous material at oblique incidence. This method is based on a temporal model of the direct and inverse scattering problems for the propagation of transient ultrasonic waves in a homogeneous isotropic slab of porous material having a rigid frame initially introduced by two of the authors ͑Z.E.A.F and C.D͒ in Refs. 1-4. The viscous and thermal losses of the medium are described by the Johnson et al. [START_REF] Johnson | Theory of dynamic permeability and tortuosity in fluid-saturated porous media[END_REF] and Allard 6 model modified by a fractional calculus based method in order to be used in time domain. Reflection and transmission scattering operators of a slab of porous material are derived for an oblique incidence and thus the responses of the medium to an incident acoustic pulse are obtained.

The outline of this paper is as follows. In Sec. II, a time domain model is given, and the basic equations of wave propagation in porous material are written for an oblique incidence. Section III is devoted to the direct problem and to the expressions of the reflection and transmission kernels in the time domain at oblique incidence. Finally in Sec. IV, an experimental validation with ultrasonic measurement is treated for two plastic foams having different values of flow resistivity. The porosity and tortuosity are measured via reflected waves for different incidence angles.

II. MODEL

In the acoustics of porous materials, a distinction can be made between two situations depending on whether the frame is moving or not. In the first case, the dynamics of the waves due to the coupling between the solid frame and the fluid are clearly described by the Biot theory. [START_REF] Biot | The theory of propagation of elastic waves in fluid-saturated porous solid. I. Low frequency range[END_REF][START_REF] Biot | The theory of propagation of elastic waves in fluid-saturated porous solid. II. Higher frequency range[END_REF] In airsaturated porous media the structure is generally motionless and the waves propagate only in the fluid. This case is described by the model of equivalent fluid which is a particular case in the Biot model, in which the interactions between the fluid and the structure are taken into account in two frequency response factors: the dynamic tortuosity of the medium ␣͑͒ given by Johnson et al. [START_REF] Johnson | Theory of dynamic permeability and tortuosity in fluid-saturated porous media[END_REF] and the dynamic compressibility of the air included in the porous material ␤͑͒ given by Allard. [START_REF] Allard | Propagation of Sound in Porous Media: Modeling Sound Absorbing Materials ͑Chapman and Hall[END_REF] In the frequency domain, these factors multiply the density of the fluid and its compressibility respectively and represent the deviation from the behavior of the fluid in free space as the frequency increases. In the time domain, they act as operators and in the asymptotic domain ͑high frequency approximation͒ their expressions are given [START_REF] Fellah | Direct and inverse scattering of transient acoustic waves by a slab of rigid porous material[END_REF][START_REF] Fellah | Transient acoustic wave propagation in rigid porous media: A time-domain approach[END_REF][START_REF] Fellah | An approach to direct and inverse time-domain scattering of acoustic waves from rigid porous materials by a fractional calculus based method[END_REF][START_REF] Fellah | Application of fractional calculus to the sound waves propagation in rigid porous materials: Validation via ultrasonic measurements[END_REF] by

␣ ˜͑t ͒ϭ␣ ϱ ͩ ␦͑t͒ϩ 2 ⌳ ͩ f ͪ 1/2 t Ϫ1/2 ͪ , ͑1͒ ␤ ˜͑t ͒ϭ ͩ ␦͑t͒ϩ 2͑␥Ϫ1 ͒ ⌳Ј ͩ P r f ͪ 1/2 t Ϫ1/2 ͪ . ͑2͒
In these equations, ␦(t) is the Dirac function, Pr is the Prandtl number, and f are, respectively, the fluid viscosity and the fluid density, and ␥ is the adiabatic constant. The relevant physical parameters of the model are the tortuosity of the medium ␣ ϱ initially introduced by Zwikker and Kosten, [START_REF] Zwikker | Sound Absorbing Materials ͑Elsevier[END_REF] the viscous and the thermal characteristic lengths ⌳ and ⌳Ј introduced by Johnson et al. [START_REF] Johnson | Theory of dynamic permeability and tortuosity in fluid-saturated porous media[END_REF] and Allard. [START_REF] Allard | Propagation of Sound in Porous Media: Modeling Sound Absorbing Materials ͑Chapman and Hall[END_REF] In this model the time convolution of t Ϫ1/2 with a function is interpreted as a semi derivative operator following the definition of the fractional derivative of order given in Samko and colleagues 20

D ͓x͑ t ͔͒ϭ 1 ⌫͑Ϫ ͒ ͵ 0 t ͑ tϪu ͒ ϪϪ1 x͑u ͒du, ͑3͒
where ⌫(x) is the Gamma function.

In this framework, the basic equations of our model can be written as

f ␣ ˜͑t ͒ * ץv i ץt ϭϪٌ i p and ␤ ˜͑t ͒ K a * ץp ץt ϭϪٌ"v, ͑4͒
where * denotes the time convolution operation, p is the acoustic pressure, v is the particle velocity, and K a is the bulk modulus of air. The first equation is the Euler equation, the second one is the constitutive equation. For a wave propagating at oblique incidence in the plane ͑xoz͒ and making an angle along the x-axis, these equations become

f ␣ ϱ ץv x ͑ x,z,t ͒ ץt ϩ 2 f ␣ ϱ ⌳ ͩ f ͪ 1/2 ͵ 0 t ץv͑x,z,tЈ͒/ץtЈ ͱtϪt Ј dtЈ ϭϪ ץp͑x,z,t ͒ ץx , f ␣ ϱ ץv z ͑ x,z,t ͒ ץt ϩ ץ f ␣ ϱ ⌳ ͩ f ͪ 1/2 ͵ 0 t ץ v ͑ x,z,tЈ͒/ץtЈ ͱtϪt Ј dtЈ ϭϪ ץp͑x,z,t ͒ ץz , ͑5͒ 1 K a ץp͑x,z,t ͒ ץt ϩ 2͑␥Ϫ1 ͒ K a ⌳Ј ͩ f P r ͪ 1/2 ͵ 0 t ץp͑x,z,tЈ͒/ץtЈ ͱtϪt Ј dtЈ ϭϪ ץv͑x,z,t ͒ ץx Ϫ ץv͑x,z,t ͒ ץz ,
where v x , v z are the components of the particle velocity along the axes xЈ and z. In these equations, the convolutions express the dispersive nature of the porous material. They take into account the memory effects due to the fact that the response of the medium to the wave excitation is not instantaneous but needs more time to take effect.

III. DIRECT PROBLEM

The direct scattering problem is that of determining the scattered field as well as the internal field, that arises when a known incident field impinges on the porous material with known physical properties. To compute the solution of the direct problem one needs to know the Green's function [START_REF] Fellah | An approach to direct and inverse time-domain scattering of acoustic waves from rigid porous materials by a fractional calculus based method[END_REF] of the modified wave equation in the porous medium. In that case, the internal field is given by the time convolution of the Green's function with the incident wave, and the reflected and transmitted fields are deduced from the internal field and the boundary conditions. In this section some notation is introduced. The geometry of the problem is shown in Fig. 1. A homogeneous porous material occupies the region 0рx рL. This medium is assumed to be isotropic and to have a rigid frame. A short sound pulse impinges at oblique incidence on the medium from the left, it gives rise to an acoustic pressure field p(x,z,t) and an acoustic velocity field v(x,z,t) within the material, which satisfying the system of Eqs. ͑5͒, which can be written as

a ץv x ͑ x,z,t ͒ ץt ϩb ץ 1/2 v x ͑ x,z,t ͒ ץt 1/2 ϭϪ ץp͑x,z,t ͒ ץx , ͑6͒ a ץv z ͑ x,z,t ͒ ץt ϩb ץ 1/2 v z ͑ x,z,t ͒ ץt 1/2 ϭϪ ץp͑x,z,t ͒ ץz , ͑7͒ d ץp͑x,z,t ͒ ץt ϩ f ץ 1/2 p͑x,z,t ͒ ץt 1/2 ϭϪ ץv x ͑ x,z,t ͒ ץx Ϫ ץv z ͑ x,z,t ͒ ץz . ͑8͒
with

aϭ f ␣ ϱ , bϭ 2 f ␣ ϱ ⌳ ͩ f ͪ 1/2 , dϭ 1 K a and f ϭ 2͑␥Ϫ1 ͒ K a ⌳Ј ͩ f P r ͪ 1/2 .
In the region xр0, the incident pressure wave is given by

p i ͑ x,z,t ͒ϭ p i ͩ tϪ x cos c 0 Ϫ z sin c 0 ͪ , ͑9͒
where c 0 is the velocity of the free fluid (xр0); c 0 ϭͱK a / f . In the region 0рxрL, the pressure wave is given by

p i ͑ x,z,t ͒ϭ p i ͩ tϪ x cos Ј cЈ Ϫ z sin Ј cЈ ͪ , ͑10͒
where cЈ is the velocity in the porous material (0рxрL), and the refraction angle Ј is given by the Descartes-Snell law:

6 sin c 0 ϭ sin Ј cЈ . ͑11͒
To simplify the system of equations

ץv z ץz ϭϪ sin 2 c 0 2 ͩ aϩb ץ Ϫ1/2 ץt Ϫ1/2ͪ Ϫ1 * ץp ץt . ͑14͒
By using Eqs. ͑6͒, ͑8͒, and ͑14͒, the equation system ͑6͒, ͑7͒ and ͑8͒ can thus be simplified to

a ץv x ͑ x,z,t ͒ ץt ϩb ץ 1/2 v x ͑ x,z,t ͒ ץt 1/2 ϭϪ ץp͑x,z,t ͒ ץx , ͑15͒ d ץp͑x,z,t ͒ ץt ϩ f ץ 1/2 p͑x,z,t ͒ ץt 1/2 ϭϪ ץv x ͑ x,z,t ͒ ץx ϩ sin 2 c 0 2 ͩ aϩb ץ Ϫ1/2 ץt Ϫ1/2ͪ Ϫ1 * ץp͑x,z,t ͒ ץt .

͑16͒

From Eqs. ͑15͒ and ͑16͒, we derive the generalized lossy wave equation in the time domain along the x-axis

ץ 2 p͑x,z,t ͒ ץx 2 ϪA ץ 2 p͑x,z,t ͒ ץt 2 ϪB ץ 3/2 p͑x,z,t ͒ ץt 3/2 ϪC ץp͑x,z,t ͒ ץt ϭ0, ͑17͒
where the coefficients A, B and C are constants, respectively, given by

Aϭ 1 c 0 2 ͑ ␣ ϱ Ϫsin 2 ͒, Bϭ 2␣ ϱ K a ͱ f ͩ 1 ⌳ ϩ ␥Ϫ1
ͱPr⌳ Ј ͪ , and ͑18͒ The first one is related to the velocity of the projected wave along the x-axis cϭc 0 /ͱ␣ ϱ Ϫsin 2 . The other coefficients are essentially dependent on the characteristic lengths ⌳ and ⌳Ј, and express the viscous and thermal interactions between the fluid and the structure. The constant B governs the spreading of the signal while C is responsible of the attenuation of the wave. Obviously, a knowledge of these three coefficients allows the determination of the parameters ␣ ϱ , ⌳, and ⌳Ј. One way to solve Eq. ͑17͒ with suitable initial and boundary conditions is by using the Laplace transform. The approach is quite simple although the inverse Laplace transform requires tedious calculi. [START_REF] Fellah | An approach to direct and inverse time-domain scattering of acoustic waves from rigid porous materials by a fractional calculus based method[END_REF] A suitable setting for the introduction of the time domain solution of the modified wave propagation equation ͑17͒ is provided by the following model.

Cϭ 4␣ ϱ ͑ ␥Ϫ1͒ K a ⌳⌳ЈͱPr .
If the incident sound wave is launched in the region x р0, then the expression of the pressure field in the region at the left of the material is the sum of the incident and reflected fields

p 1 ͑ x,t ͒ϭ p i ͩ tϪ x cos c 0 ͪ ϩp r ͩ tϩ x cos c 0 ͪ , xϽ0.

͑19͒

Here, p 1 (x,t) is the field in the region xϽ0, p i and p r denote the incident and the reflected fields, respectively. In addition, a transmitted field is produced in the region at the right of the material. This has the form

p 3 ͑ x,t ͒ϭ p t ͩ tϪ L c Ϫ ͑ xϪL ͒cos c 0 ͪ , xϾL ͑20͒
͓ p 3 (x,t) is the field in the region xϾL and p t is the transmitted field͔.

The incident and scattered fields are related by the scattering operators ͑i.e., the reflection and transmission opera-tors͒ for the material. Theses are integral operators represented by

p r ͑ x,t ͒ϭ ͵ 0 t R ˜͑ ͒p i ͩ tϪϩ x c 0 ͪ d ϭR ˜͑t ͒ * p i ͑ t ͒ * ␦ ͩ tϩ x cos c 0 ͪ , ͑21͒ p t ͑ x,t ͒ϭ ͵ 0 t T ˜͑ ͒p i ͩ tϪϪ L c Ϫ ͑ xϪL ͒ c 0 ͪ d ϭT ˜͑t ͒ * p i ͑ t ͒ * ␦ ͩ tϪ L c Ϫ ͑ xϪL ͒cos c 0 ͪ . ͑22͒
In Eqs. ͑21͒ and ͑22͒ the functions R ˜and T ˜are the reflection and transmission kernels, respectively, for incidence from the left. Note that the lower limit of integration in Eqs. ͑21͒ and ͑22͒ is chosen to be 0, which is equivalent to assume that the incident wave front first impinges on the material at tϭ0.

The scattering operators given in Eqs. ͑21͒ and ͑22͒ are independent of the incident field used in the scattering experiment and depend only on the properties of the materials.

To derive the reflection and transmission coefficients, the boundary conditions of the flow velocity at the interfaces xϭ0 and xϭL are needed. [START_REF] Fellah | Direct and inverse scattering of transient acoustic waves by a slab of rigid porous material[END_REF] Using the same calculus steps given in Ref. 1, we obtain the reflection and transmission scattering operators for a slab of porous material at oblique incidence

R ˜͑t ͒ϭ ͩ 1ϪE 1ϩE ͚ͪ nу0 ͩ 1ϪE 1ϩE ͪ 2n ϫ ͫ F ͩ t,2n L c ͪ ϪF ͩ t,͑2nϩ2 ͒ L c ͪͬ , ͑23͒ T ˜͑t ͒ϭ 4E ͑ 1ϩE ͒ 2 ͚ nу0 ͩ 1ϪE 1ϩE ͪ 2n F ͩ tϩ L c 0 ,͑2nϩ1 ͒ L c ͪ , ͑24͒ with Eϭ ͱ1Ϫ sin 2 ␣ ϱ ͱ␣ ϱ cos . ͑25͒
Theses expressions take into account the n-multiple reflection in the material.

In most cases of porous materials saturated by air, the multiple reflection effects can be negligible because of the high attenuation of sound waves in these media. This depends on the thickness L and the flow resistivity of the porous material.

In the case where the multiple reflections can be neglected, the kernel of transmission is given by

T ˜͑t ͒ϭ 4E ͑ 1ϩE ͒ 2 F ͩ tϩ L c , L c ͪ ,
and the kernel of reflection by R ˜͑t ͒ϭr͑ t ͒ϩR͑ t ͒, ͑26͒

with

r͑t ͒ϭ ͩ 1ϪE 1ϩE ͪ ␦͑t͒ and ͑27͒ R͑t ͒ϭϪ 4E͑1ϪE ͒ ͑ 1ϩE ͒ 3 F ͩ t, 2L c ͪ .

r(t) is the instantaneously response of the porous material corresponding to the reflection contribution at the first interface (xϭ0).

In the case of a semi-infinite medium when L→ϱ, it can be shown [START_REF] Fellah | Direct and inverse scattering of transient acoustic waves by a slab of rigid porous material[END_REF]

that F(t,2L/c)→0 and R ˜(t)→r(t).

This means that R(t) is equivalent to the reflection at the rear interface (xϭL), which is the bulk contribution to the reflection. The part of the wave corresponding to r(t) is not subjected to the dispersion but is simply multiplied by the factor (1ϪE)/(1 ϩE). This shows that although tortuosity is a bulk parameter, it may be evaluated from the wave reflected at the first interface when the porosity is known and vice versa. Although generally speaking it is easy to evaluate the tortuosity from transmitted waves [START_REF] Fellah | Direct and inverse scattering of transient acoustic waves by a slab of rigid porous material[END_REF][START_REF] Fellah | An approach to direct and inverse time-domain scattering of acoustic waves from rigid porous materials by a fractional calculus based method[END_REF][START_REF] Fellah | Application of fractional calculus to the sound waves propagation in rigid porous materials: Validation via ultrasonic measurements[END_REF][START_REF] Leclaire | Determination of viscous and thermal characteristics lengths of plastic foams by ultrasonic measurements in helium and air[END_REF][START_REF] Brown | Evaluation of viscous characteristic length of air-saturated porous materials from the ultrasonic dispersion curve[END_REF] this is not the case for porosity because of its weak sensitivity in transmitted mode. [START_REF] Fellah | Direct and inverse scattering of transient acoustic waves by a slab of rigid porous material[END_REF] Figure 2 shows two simulated transmitted signals for a plastic foam M1. The first one ͑solid line͒ corresponds to the value of porosity 1 ϭ0.98 and the second one ͑dashed line͒ corresponds to 2 ϭ0.49. The parameters of the plastic foam M1 ͑thickness 1.1 cm, ϭ0.92, ␣ ϱ ϭ1.25, ⌳ϭ50 m and ⌳Јϭ150 m, flow resistivity ϭ38 000 Nm Ϫ4 s) have been determined using conventional methods. [START_REF] Fellah | Direct and inverse scattering of transient acoustic waves by a slab of rigid porous material[END_REF][START_REF] Fellah | An approach to direct and inverse time-domain scattering of acoustic waves from rigid porous materials by a fractional calculus based method[END_REF][START_REF] Fellah | Application of fractional calculus to the sound waves propagation in rigid porous materials: Validation via ultrasonic measurements[END_REF][START_REF] Leclaire | Determination of viscous and thermal characteristics lengths of plastic foams by ultrasonic measurements in helium and air[END_REF][START_REF] Brown | Evaluation of viscous characteristic length of air-saturated porous materials from the ultrasonic dispersion curve[END_REF] Readers can see the slight difference between the two curves for a 50% difference in porosity values, which is due to the dispersion phenomenon that is governed by viscous, thermal and inertial effects contributed by ␣ ϱ , ⌳ and ⌳Ј and plays a more important role in the Green function F(t,k) than .

In Fig. 3 we show by numerical simulation, the difference between the reflected wave at the first interface and the total reflected wave. The difference between the two curves is marginal. This means that the wave reflected by the slab may be approximated by the reflected wave by the first interface r(t) with a good level of accuracy. It is easy to verify that this approximation is valid for any incidence angle . In the next section, we will show, by a simple measurements, how it is the possible to determine the porosity and the tortuosity by measuring reflected wave by a slab of porous material.

IV. ULTRASONIC MEASUREMENT

In this section, we measure the tortuosity and the porosity knowing the reflection coefficient at the first interface for different values of the incidence angle . The expression of the reflection coefficient at the first interface is given by r͑t, ͒ϭ ␣ ϱ cos Ϫͱ␣ ϱ Ϫsin 2 ␣ ϱ cos ϩͱ␣ ϱ Ϫsin 2 ␦͑t͒.

͑28͒

For two values of the incidence angle 1 and 2 , it is easy to calculate the expression of the tortuosity function of the reflection coefficients r 1 ϭr(t)͉ 1 and r 2 ϭr(t)͉ 2 corresponding to the angles 1 and 2 , respectively,

␣ ϱ ϭ ͩ ͑ 1Ϫr 2 ͒͑ 1ϩr 1 ͒cos 2 ͑ 1ϩr 2 ͒͑ 1Ϫr 1 ͒cos 1 ͪ 2 sin 2 1 Ϫsin 2 2 ͩ ͑ 1Ϫr 2 ͒͑ 1ϩr 1 ͒cos 2 ͑ 1ϩr 2 ͒͑ 1Ϫr 1 ͒cos 1 ͪ 2

Ϫ1

.

͑29͒

Knowing the value of the tortuosity, we deduce the expression of the porosity function of i and r i by the expression

ϭ ␣ ϱ ͑ 1Ϫr i ͒cos i ͑ 1ϩr i ͒ ͱ ␣ ϱ Ϫsin 2 i , iϭ1,2. ͑30͒
As an application of this model, some numerical simulations are compared to experimental results. Experiments are carried out in air with two broadband Ultran NCT202 transducers having a 190 kHz central frequency in air and a bandwidth at 6 dB extending from 150 kHz to 230 kHz. A goniometer used in optic has been employed for the positioning of the transducers. Pulses of 400 V are provided by a 5052PR Panametrics pulser/receiver. The received signals are amplified up to 90 dB and filtered above 1 MHz to avoid high frequency noise. Electronic perturbations are removed by 1000 acquisition averages. The experimental setup is shown in Fig. 4. The distance between the transducers and the samples is 20 cm. The duration of the signal plays an important role, its spectrum must verify the condition of high frequency approximation referred to in the previous section.

The parameters of the first investigated plastic foam M2 are: resistivity ϭ2500 Nm Ϫ4 s, 4.1 cm, ⌳ ϭ230 m, and ⌳Јϭ690 m. Figure 5 shows the incident signal generated by the transducer ͑solid line͒ and the reflected signal by the plastic foam M2 ͑dashed line͒ for normal incidence. Figure 6 shows their spectra. From the spectra of the two signals, the reader can see that they have practically the same bandwidth which means that there is no dispersion. The reflected signal from the foam M2 is very small compared with the incident signal because of the low value of its flow resistivity .

Figure 7 represents the incident signal generated by the transducer ͑solid line͒ and the reflected signal by the plastic foam M2 ͑dashed line͒ for an incidence angle ϭ53 deg. Figure 8 depicts their spectra, showing that there is also no dispersion in the diffusion process. This concurs with the theory below which predicts that the reflected wave from the first interface (xϭ0) is measured and simply attenuated by the factor (␣ ϱ cos Ϫͱ␣ ϱ Ϫsin 2 )/(␣ ϱ cos ϩͱ␣ ϱ Ϫsin 2 ). Table I gives the experimental data of the reflection coefficient for different values of the incidence angle and Table II give the values of the tortuosity and porosity calculated using of Eqs. ͑29͒ and ͑30͒ for each pair of incidence angles. The average value of the tortuosity obtained from these measurements is ␣ ϱ ϭ1.07Ϯ0.09 and the average value of the porosity is ϭ0.97Ϯ0.01.

Figure 9 shows the experimental data of the reflection coefficient for different values of the incident angle , and the simulation of the variation of the reflection coefficient with the incident angle for the values of the tortuosity ␣ ϱ ϭ1.07 and porosity ϭ0.97. The value of the porosity of the plastic foam M2 given by the porosimeter 18 is ϭ0.98 Ϯ0.02, and the value of the tortuosity given by classical method 1,3,4,19,21 is ␣ ϱ ϭ1.06Ϯ0.08. The reader can see the slight difference between the values of the porosity and tortuosity measured using this method and the other classical methods. [START_REF] Fellah | Direct and inverse scattering of transient acoustic waves by a slab of rigid porous material[END_REF][START_REF] Fellah | An approach to direct and inverse time-domain scattering of acoustic waves from rigid porous materials by a fractional calculus based method[END_REF][START_REF] Fellah | Application of fractional calculus to the sound waves propagation in rigid porous materials: Validation via ultrasonic measurements[END_REF][START_REF] Champoux | Air-based system for the measurement of porosity[END_REF][START_REF] Leclaire | Determination of viscous and thermal characteristics lengths of plastic foams by ultrasonic measurements in helium and air[END_REF][START_REF] Brown | Evaluation of viscous characteristic length of air-saturated porous materials from the ultrasonic dispersion curve[END_REF] Figure 10 shows the comparison between the simulated reflected signal at the first interface calculated for ␣ ϱ ϭ1.07 and ϭ0.97, and the experimental reflected signal for ϭ17 deg. Figure 11 shows the same comparison for ϭ53 deg. The difference between the simulated reflected signal and experimental reflected signal is relatively weak which leads to the conclusion that the values of the tortuosity and porosity obtained are satisfactory.

Another plastic foam M3 ͑thickness 2.1 cm, ⌳ ϭ23 m, and ⌳Јϭ69 m) with a higher value of the flow resistivity ϭ125 000 Nm Ϫ4 s has been investigated. Figure 12 shows the incident and reflected signals for ϭ13 deg and Fig. 13 shows their spectra. Here also, we can conclude from the two spectra that there is no dispersion. Figure 14 shows the incident and reflected signals for ϭ49 deg. Figure 15 shows their spectra. For this plastic foam M3, the reflected signal is more important than the one of the foam M2 because of its high value of the flow resistivity. The approximation that the reflected signal is essentially due to the reflection at the first interface is well verified in this case.

Table III gives the experimental data of the reflection coefficient for different values of the incidence angle and Table IV gives the values of the tortuosity and porosity calculated with Eqs. ͑29͒ and ͑30͒ for each pair of the incidence angles. The average value of the tortuosity obtained is ␣ ϱ ϭ1.72Ϯ0.02 and the average value of the porosity obtained 16 shows the simulated variation of the reflection coefficient at the first interface for the value of the tortuosity ␣ ϱ ϭ1.72 and porosity ϭ0.82 and experimental data of the reflection coefficient. The value of the porosity of the plastic foam M3 given by the porosimeter 18 is ϭ0.8Ϯ0.1, and the value of the tortuosity obtained by classical methods 1,3,4,19,21 is ␣ ϱ ϭ1.7Ϯ0.1. Figure 17 shows the comparison between the simulated reflected signal at the first interface calculated for ␣ ϱ ϭ1.72 and ϭ0.8, and the experimental reflected signal for ϭ13 deg. The difference between the two results is small which leads to the conclusion that the values of the tortuosity and porosity obtained are satisfactory.

V. CONCLUSION

A method for the measurement of the porosity and the tortuosity simultaneously from the measurements of reflected waves at two oblique angles of incidence has been proposed. Generally, it is easy to evaluate the tortuosity from transmitted waves, this is not the case for porosity because of its weak sensitivity in transmitted mode. This method is an alternative to the usual method involving the use of a porosim-eter introduced by Beranek and improved by Champoux et al. or the other ultrasonic methods developed by Leclaire et al.

The method is based upon the propagation equation in time domain in a slab of porous material having a rigid frame at asymptotic domain ͑high frequency range͒. A time domain model of wave propagation in such material is worked out from the concept of fractional calculus. The kernels of the reflection and transmission scattering operators are derived giving simple relations between these operators and the parameters of the medium. Because of the high attenuation in such media, the reflected wave by the porous material can be approximated by the reflected wave at the first interface. As a consequence, it gives a simple relation between the porosity, tortuosity, incident angle and reflected wave, proving in such a way that all physical parameters of acoustic material can be evaluated from acoustic measurements.

An experimental validation of the theoretical expression of the reflection scattering operators illustrates the high level of correspondence between numerical and experimental results and show that this time domain model is well suited for the characterization of porous media via acoustic wave propagation.

FIG. 1 .

 1 FIG. 1. Geometry of the problem.

FIG. 2 .

 2 FIG. 2. Transmitted simulated signals for ϭ0.98 ͑solid line͒ and for ϭ0.49 ͑dashed line͒.

FIG. 6 .

 6 FIG.6. Spectrum of experimental incident signal ͑solid line͒ and spectrum of experimental reflected signal ͑dashed line͒ for ϭ0 deg.

FIG. 13 .

 13 FIG.13. Spectrum of experimental incident signal ͑solid line͒ and spectrum of experimental reflected signal ͑dashed line͒ for ϭ13 deg.

FIG. 16 .

 16 FIG. 16. Simulation of the variation of the reflection coefficient ͑solid line͒ with the incident angle and experimental data of the reflection coefficient ͑star͒ for the material M3.

  ͑5͒, we can then use the following property:

	ץ ץz	ϭϪ	sin c 0	ץ ץt	,					͑12͒
	which implies
	ץv z ץz	ϭϪ	sin c 0	ץv z ץt	and	ץp ץz	ϭϪ	sin c 0	ץp ץt	.	͑13͒

From Eqs. ͑7͒ and ͑13͒, we obtain then the relation

TABLE II .

 II Values of the tortuosity and porosity calculated for each pair of incidence angles 1 and 2 ͑foam M2͒.

TABLE III .

 III Experimental data of the reflection coefficient for different values of the incidence angle ͑foam M3͒.

	Incidence angle ͑deg͒	0	13.5	23	32.5	39.7	49.8	54	60.3
	Reflection coefficient r	0.2274	0.221	0.2106	0.1905	0.1668	0.1163	0.0864	0.024

TABLE IV .

 IV Values of the tortuosity and porosity calculated for each pair of incidence angles 1 and 2 ͑foam M3͒.

		Incidence angle 1 / 2 ͑deg͒	Tortuosity	Porosity
	͑a͒
		0/13.5	1.9	0.86
		0/23	1.680	0.816
		0/32.5	1.700	0.820
		0/39.7	1.715	0.824
		0/49.8	1.712	0.823
		0/54	1.709	0.822
		0/60.6	1.729	0.827
		13.5/23	1.592	0.796
		13.5/32.5	1.665	0.814
		13.5/39.7	1.695	0.821
		13.5/49.8	1.701	0.822
		13.5/54	1.704	0.823
		13.5/60.3	1.722	0.827
		23/32.5	1.714	0.823
	͑b͒
		23/39.7	1.729	0.826
	is ϭ0.82Ϯ0.04. Figure	23/49.8 23/54	1.719 1.719	0.824 0.824
		23/60.3	1.735	0.828
		32.5/39.7	1.743	0.829
		32.5/49.8	1.72	0.824
		32.5/54	1.721	0.824
		32.5/60	1.739	0.828
		39.7/49.8	1.709	0.823
		39.7/54	1.713	0.824
		39.7/60.3	1.738	0.821
		49.8/54	1.695	0.829
		49.8/60.3	1.761	0.829
		54/60.3	1.801	0.832

J. Acoust. Soc. Am. 113 (5), May 2003 0001-4966/2003/113(5)/2424/10/$19.00 © 2003 Acoustical Society of America

J. Acoust. Soc. Am., Vol. 113, No. 5, May 2003 Fellah et al.: Pulse propagation in porous media

J. Acoust. Soc. Am., Vol. 113, No. 5, May 2003 Fellah et al.: Pulse propagation in porous media

J. Acoust. Soc. Am., Vol. 113, No. 5, May 2003 Fellah et al.: Pulse propagation in porous media

J. Acoust. Soc. Am., Vol. 113, No. 5, May 2003 Fellah et al.: Pulse propagation in porous media

J. Acoust. Soc. Am., Vol. 113, No. 5, May 2003 Fellah et al.: Pulse propagation in porous media

ACKNOWLEDGMENT

The authors would like to express all thanks to the reviewers of the paper for their critical remarks.


Z E A Fellah
email: zine.fellah@fys.kukeuven.ac.be

) S Berger

W Lauriks

C Depollier

M Fellah

Chapman F Allard

Measuring the porosity of porous materials having a rigid frame via reflected waves: A time domain analysis with fractional derivatives

An ultrasonic reflectivity method is proposed for measuring porosity of porous materials having a rigid frame. Porosity is the relative fraction by volume of the air contained within a material. It is important as one of the several parameters required by acoustical theory to characterize porous materials like plastic foams and fibrous or granular materials. The proposed method is based on a temporal model of the direct and inverse scattering problem for the propagation of transient ultrasonic waves in a homogeneous isotropic slab of porous material having a rigid frame. This time domain model of wave propagation was initially introduced by the authors ͓Z.

I. INTRODUCTION

The ultrasonic characterization of porous materials saturated by air such as plastic foams, fibrous or granular materials is of great interest for a wide range of industrial applications. These materials are frequently used in the automotive and aeronautics industries and in the building trade. The determination of the properties of a medium using waves that have been reflected by or transmitted through the medium is a classical inverse scattering problem. Such problems are often approached by taking a physical model of the scattering process, generating a synthetic response for a number of assumed values for the parameters, and adjusting these parameters until a reasonable level of correspondence is attained between the synthetic response and the data observed. One important parameter which appears in theories of sound propagation in porous materials [1][2][3][4][5][START_REF] Allard | Propagation of Sound in Porous Media ͑Chapman and Hall[END_REF][START_REF] Biot | [END_REF][8][9][10][START_REF] Zwikker | Sound Absorbing Materials ͑Elsevier[END_REF][START_REF] Champoux | [END_REF][13] is porosity. Porosity is the relative fraction, by volume, of the air contained within the material. Unlike other parameters included in the description of different various phenomena occurring in the acoustical propagation of porous media at high frequency range such as tortuosity, [START_REF] Zwikker | Sound Absorbing Materials ͑Elsevier[END_REF] viscous characteristic length, 5 and thermal characteristic length, [START_REF] Champoux | [END_REF] or at low frequency range such as flow resistivity [START_REF] Zwikker | Sound Absorbing Materials ͑Elsevier[END_REF] and thermal permeability, 13 porosity is a key parameter that plays an important role in propagation at all frequencies. As such, in studies of acoustical properties of porous materials, it is extremely useful to be able to measure this parameter.

Beranek 10 described an apparatus ͑porosimeter͒ used to measure the porosity of porous materials. This device was based on the equation of state for ideal gases at constant temperature ͑i.e., Boyle's law͒. Porosity can be determined by measuring the change in air pressure occurring with a known change in volume of the chamber containing the sample. In the Beranek apparatus, both pressure change and volume change are monitored using a U-shaped fluid-filled manometer. An alternate technique for measuring porosity is a dynamic method proposed by Leonard. 14 Techniques that use water as the pore-filling fluid, rather than air, are common in geophysical studies. 15,16 Mercury has been used as the pore-filling fluid in other applications. 17 However, for many materials, the introduction of liquids into the material is not appropriate. Recently a similar device to that of Beranek, involving the use of an electronic pressure transducer, was introduced by Champoux et al. 18 This device can be used to measure very slight changes in pressure accurately, and the output can be recorded by a computer.

In this work, we present a simple method of measuring porosity by measuring a wave reflected by a slab of porous material. This method is based on a temporal model of the direct and inverse scattering problem for the propagation of transient ultrasonic waves in a homogeneous isotropic slab of porous material having a rigid frame. This work was initially introduced by the authors of Refs. 1-4. The viscous and thermal losses of the medium are described by the Johnson et al. 5 and Allard 6 model modified by a fractional calculus based method in order to be used in the time domain. Reflection and transmission scattering operators of a slab of porous material are derived and thus the responses of the medium to an incident acoustic pulse are obtained.

The outline of this work is as follows. In Sec. II, a time domain model is given, the connection between the fractional derivatives and wave propagation in rigid porous media in high frequency range is established and the basic equations are written in the time domain. Section III is devoted to the direct problem and to the expressions of the reflection and transmission kernels in the time domain and numerical results are given. Section IV deals with the inverse problem and an expression of porosity based on the incident and reflected waves from a slab of porous materials is established. Finally in Sec. V, the experimental results are given to validate this method.

II. MODEL

In the acoustics of porous materials, a distinction can be made between two situations depending on whether the frame is moving or not. In the first case, the dynamics of the waves due to the coupling between the solid frame and the fluid are clearly described by the Biot theory. [START_REF] Biot | [END_REF]8 In airsaturated porous media the structure is generally motionless and the waves propagate only in the fluid. This case is described by the model of equivalent fluid which is a particular case in the Biot model, in which the interactions between the fluid and the structure are taken into account in two frequency response factors: the dynamic tortuosity of the medium ␣() given by Johnson et al. 5 and the dynamic compressibility of the air included in the porous material ␤() given by Allard. [START_REF] Allard | Propagation of Sound in Porous Media ͑Chapman and Hall[END_REF] In the frequency domain, these factors multiply the density of the fluid and its compressibility respectively and represent the deviation from the behavior of the fluid in free space as the frequency increases. In the time domain, they act as operators and in the high frequency approximation their expressions are given by Fellah and Depollier 1,2 and Fellah et al. 3,4 as

␣ ˜͑t ͒ϭ␣ ϱͫ ␦͑t͒ϩ 2 ⌳ ͩ f ͪ 1/2 t Ϫ1/2 ͬ , ͑1͒ ␤ ˜͑t ͒ϭ ͫ ␦͑t͒ϩ 2͑␥Ϫ1 ͒ ⌳Ј ͩ Pr f ͪ 1/2 t Ϫ1/2 ͬ .

͑2͒

In these equations, ␦(t) is the Dirac function, Pr is the Prandtl number, and f the fluid viscosity and the fluid density, respectively, and ␥ is the adiabatic constant. The relevant physical parameters of the model are the tortuosity of the medium ␣ ϱ initially introduced by Zwikker and Kosten [START_REF] Zwikker | Sound Absorbing Materials ͑Elsevier[END_REF] and the viscous and the thermal characteristic lengths ⌳ and ⌳Ј introduced by Johnson et al. 5 and Allard. [START_REF] Allard | Propagation of Sound in Porous Media ͑Chapman and Hall[END_REF] In this model, the time convolution of t Ϫ1/2 with a function is interpreted as a semi-derivative operator following the definition of the fractional derivative of order given in Samko et al. [START_REF] Samko | Fractional Integrals and Derivative: Theory and Applications ͑Gordon and Breach Science[END_REF] D ͓x͑ t ͔͒ϭ 1

⌫͑Ϫ ͒ ͵ 0 t ͑ tϪu ͒ ϪϪ1 x͑u ͒du, ͑3͒
where ⌫(x) is the gamma function.

In this framework, the basic equations of our model can be expressed as follows:

f ␣ ˜͑t ͒ * ץv i ץt ϭϪٌ i p and ␤ ˜͑t ͒ K a * ץp ץt ϭϪٌv, ͑4͒
where * denotes the time convolution operation, p is the acoustic pressure, v is the particle velocity, and K a is the bulk modulus of the air. The first equation is the Euler equation; the second is a constitutive equation obtained from the equation of mass conservation associated with the behavior ͑or adiabatic͒ equation. For a wave propagating along the x axis, these equations become

f ␣ ϱ ץv ץt ϩ2 f ␣ ϱ ⌳ ͩ f ͪ 1/2 ͵ 0 t ץv/ץtЈ ͱtϪt Ј dtЈϭϪ ץp ץx , ͑5͒ 1 K a ץp ץt ϩ2 ␥Ϫ1 K a ⌳Ј ͩ Pr f ͪ 1/2 ͵ 0 t ץp/ץtЈ ͱtϪt Ј dtЈϭϪ ץv ץx . ͑6͒
In these equations, the convolutions express the dispersive nature of the porous material. They take into account the memory effects due to the fact that the response of the medium to the wave excitation is not instantaneous but needs more time to take effect. The retarding force is no longer proportional to the time derivative of the acoustic velocity but is found to be proportional to the fractional derivative of order 1/2 of this quantity. This occurs because the volume of fluid involved in the motion is not the throughout the whole length of the signal as it is in the case of a fully developed steady flow. The phenomenon may be understood by considering such a volume of fluid in a pore to be in harmonic motion. At high frequencies, only a thin layer of fluid is excited and so the average shear stress is high. At a lower frequency, the same amplitude of fluid motion allows a thicker layer of fluid to participate in the motion and consequently the shear stress is lower. The penetration distance of the viscous forces and therefore the excitation of the fluid depends on frequency. In the time domain, such dependence is associated with a fractional derivative.

III. DIRECT PROBLEM

The direct scattering problem is that of determining the scattered field, as well as the internal field, that arises when a known incident field impinges on the porous material with known physical properties. To compute the solution of the direct problem, it is necessary to know the reflection and/or the transmission scattering operators which depends on the Green 3 function of the porous medium. In that case, the reflected field is given by the convolution of the reflection operator with the incident field and the transmitted field is given by the convolution of the transmission operator by the incident field.

Consider a homogeneous slab of porous material which occupies the region 0рxрL. The geometry of the problem is shown in Fig. 1. This medium is assumed to be isotropic and to have a rigid frame. A short sound pulse impinges normally on the medium from the left. It gives rise to an acoustic pressure field p(x,t) and an acoustic velocity field v(x,t) within the material, which satisfies the propagation equation 1

ץ 2 p ץx 2 ϪA ץ 2 p ץt 2 ϪB ͵ 0 t ץ 2 p/ץtЈ 2 ͱtϪt Ј dtЈϪC ץp ץt ϭ0, ͑7͒
where the coefficients A, B, and C are constants expressed, respectively, by

Aϭ 1 c 2 ϭ f ␣ ϱ K a , Bϭ 2␣ ϱ K a ͱ f ͩ 1 ⌳ ϩ ␥Ϫ1 ͱPr⌳ Ј ͪ , ͑8͒ Cϭ 4␣ ϱ ͑ ␥Ϫ1͒ K a ⌳⌳ЈͱPr .
The first one is related to the velocity cϭ1/ͱ f ␣ ϱ /K a of the wave in the air included in the porous material. ␣ ϱ is the refractive index of the medium which changes the wave velocity from c 0 ϭͱK a / f in free space to cϭc 0 /ͱ␣ ϱ in the porous medium. The other coefficients are essentially dependent on the characteristic lengths ⌳ and ⌳Ј and express the viscous and thermal interactions between the fluid and the structure. The constant B governs the spreading of the signal while C is responsible for the attenuation of the wave.

To derive the reflection and transmission scattering operators, it is assumed that the pressure field and the flow velocity are continuous at the boundary of the material

p͑0 ϩ ,t ͒ϭ p͑0 Ϫ ,t ͒, p͑L Ϫ ,t ͒ϭ p͑L ϩ ,t ͒, ͑9͒ v͑ 0 Ϫ ,t ͒ϭv͑ 0 ϩ ,t ͒, v͑ L ϩ ,t ͒ϭv͑ L Ϫ ,t ͒,
where is the porosity of the medium and Ϯ superscript denotes the limit from the left and the right, respectively. The initials conditions are given by p͑x,t ͉͒ tϭ0 ϭ0, ͑10͒ ץp ץt ͯ tϭ0 ϭ0, which means that the medium is idle for tϭ0.

If the incident sound wave is launched in the region x р0, then the general solution of Eq. ͑7͒ in the region to the left of the material is the sum of the incident and reflected field

p 1 ͑ x,t ͒ϭ p i ͩ tϪ x c 0 ͪ ϩ p r ͩ tϩ x c 0 ͪ , xϽ0. ͑11͒
Here, p 1 (x,t) is the field in the region xϽ0, p i and p r denote the incident and the reflected field, respectively. In addition, a transmitted field is produced in the region to the right of the material. This has the form

p 3 ͑ x,t ͒ϭ p t ͩ tϪ L c Ϫ ͑ xϪL ͒ c 0 ͪ , xϾL,

͑12͒

where p 3 (x,t) is the field in the region xϾL and p t is the transmitted field.

The incident and scattered fields are related by the scattering operators ͑i.e., reflection and transmission operators͒ for the material. These are integral operators represented by

p r ͑ x,t ͒ϭ ͵ 0 t R ˜͑ ͒p i ͩ tϪϩ x c 0 ͪ d ϭR ˜͑t ͒ * p i ͑ t ͒ * ␦ ͩ tϩ x c 0 ͪ , ͑13͒ p t ͑ x,t ͒ϭ ͵ 0 t T ˜͑ ͒p i ͩ tϪϪ L c Ϫ ͑ xϪL ͒ c 0 ͪ d ϭT ˜͑t ͒ * p i ͑ t ͒ * ␦ ͩ tϪ L c Ϫ ͑ xϪL ͒ c 0 ͪ .
In Eqs. ͑13͒ the functions R ˜and T ˜are the reflection and the transmission kernels, respectively, for incidence of waves from the left. Note that the lower limit of integration in Eq. ͑13͒ is set 0, which is equivalent to assuming that the incident wave front first impinges on the material at tϭ0. The scattering operators given in Eqs. ͑13͒ are independent of the incident field used in scattering experiment and depend only on the properties of the materials. Using Eqs. ͑4͒, ͑7͒, ͑9͒, ͑13͒, and Laplace transform calculus, we can derive the reflection and transmission scattering operators given by In most cases, for porous materials saturated by air, the multiple reflection effects are negligible because of the high attenuation of sound waves in this kind of media. So, by taking into account only the reflections at the interfaces x ϭ0 and xϭL, the kernel of transmission is given by

R ˜͑t ͒ϭ ͩ Ϫϩͱ␣ ϱ ϩͱ␣ ϱ ͚ͪ nу0 ͩ Ϫͱ␣ ϱ ϩͱ␣ ϱ ͪ 2n ϫ ͫ F ͩ t,2n L c ͪ ϪF ͩ t,͑2nϩ2 ͒ L c ͪͬ , ͑14͒ T ˜͑t ͒ϭ 4ͱ␣ ϱ ͑ ͱ␣ ϱ ϩ ͒ 2 ͚ nу0 ͩ Ϫͱ␣ ϱ ϩͱ␣ ϱ ͪ 2n ϫF ͩ tϩ L c 0 ,͑2nϩ1 ͒ L c ͪ .
T ˜͑t ͒ϭ 4ͱ␣ ϱ ͑ ϩͱ␣ ϱ ͒ 2 F ͩ tϩ L c , L c ͪ ,
4ͱ␣ ϱ ͑ ͱ␣ ϱ Ϫ ͒ ͑ ͱ␣ ϱ ϩ ͒ 3 F ͩ t, 2L c ͪ ,

͑18͒

in the case of a semi-infinite medium when L→ϱ; F(t,2L/c)→0 and R ˜(t)→r(t) ͑see Appendix͒. This means that r(t) is equivalent to the reflection at the interface xϭ0 and that R(t) is equivalent to reflection at the interface x ϭL, which is the bulk contribution to the reflection. The part of the wave corresponding to r(t) is not subjected to the dispersion but it is simply multiplied by the factor (ͱ␣ ϱ Ϫ)/(ͱ␣ ϱ ϩ). This shows that although tortuosity is a bulk parameter, it may be evaluated from the wave reflected at the first interface when the porosity is known and vice versa. Although generally speaking it is easy to evaluate the tortuosity from transmitted waves 3,4,[START_REF] Leclaire | [END_REF]21 this is not the case for porosity because of its weak sensitivity in transmitted mode. Figure 2 shows two simulated transmitted signals for a plastic foam M 1, the first one ͑solid line͒ corresponds to the value of porosity 1ϭ0.98 and the second one ͑dashed line͒ corresponds to 2ϭ0.49. The parameters used in the simulation namely thickness: 1.1 cm, ϭ0.92, ␣ ϱ ϭ1.25, ⌳ϭ50 m, and ⌳Јϭ150 m have been determined using conventional methods. 3,4,[START_REF] Samko | Fractional Integrals and Derivative: Theory and Applications ͑Gordon and Breach Science[END_REF][START_REF] Leclaire | [END_REF] Readers can see the slight difference between the two curves for a 50% difference in porosity values, which is due to the dispersion phenomenon that is governed by viscous, thermal, and inertial effects contributed by ␣ ϱ , ⌳, and ⌳Ј, and plays a more important role in the Green function F(t,k) than .

Figure 3 shows the two contributions to the reflection: r(t) and R(t) for a plastic foam M 1. The contribution of R(t) to the reflection is negligible when it is compared to the contribution of r(t). In Fig. 4 we show by numerical simulation, the difference between the reflected wave at the first interface and the total reflected wave. The difference between the two curves is slight. This means that the wave reflected by the slab may be approximated by the reflected wave by the first interface r(t) with a good level of accuracy.

IV. INVERSE PROBLEM

A porous material having a rigid frame is characterized by four parameters in the high frequency range, namely, porosity , tortuosity ␣ ϱ , the viscous characteristic length ⌳, and the thermal characteristic length ⌳Ј. These values are crucial for the behavior of the sound waves in such materials. So, it is therefore fairly important to work out experimental methods and efficient tools to assess them. Therefore, a basic inverse problem associated with the slab may be stated as follows: to find the values of the parameters of the medium from the measurements of the transmitted and/or reflected signals outside the slab. As shown in Sec. II, the solution to the direct problem is the system of two operators expressed as functions on , ␣ ϱ , ⌳, and ⌳Ј. The inversion algorithm for finding the values of the parameters of the slab is based on a fitting procedure. To find the values of the parameters , ␣ ϱ , ⌳, and ⌳Ј such that the transmitted and reflected signal describes the scattering problem in the best possible way ͑e.g., in the least-squares sense͒. Solving the inverse problem using experimental data from transmitted waves is treated 3,4,22 and an estimation of ␣ ϱ , ⌳, and ⌳Ј has been given with a good level of correspondence from ultrasonic measurements. Here we will try to solve the inverse problem using experimental reflected wave data, with the aim of returning to the porosity which, because of its weak sensitivity, cannot be characterized in transmitted mode.

The inverse problem is to find the porosity value which minimizes the function

Uϭ ͚ iϭ1 iϭN ͑ p r ͑ t i ͒Ϫ p exp r ͑ t i ͒͒ 2 , ͑19͒
where p r (t i ) iϭ1,2,...N represents the discrete set of values of the simulated reflected signal given by Eq. ͑13͒ and p exp r (t i ) iϭ1,2,...N is the discrete set of values of the experimental signal. In the approximation that the reflected wave is essentially due to the reflection at the first interface xϭ0 R ˜͑t ͒ϭr͑ t ͒ϭ aϪ1 aϩ1

␦͑t͒.

͑20͒

where aϭͱ␣ ϱ / and p r ͑ t i ͒ϭ aϪ1 aϩ1 p i ͑ t i ͒. ͑21͒

p i (t i ) is the discrete set of values of the experimental incident signal given out by the transducer. To minimize U(a) we must have

ץU ץa ϭ0,
and

ץ 2 U ץa 2 Ͼ0, ͑22͒ ץU ץa ϭ ץ ץa ͚ iϭ1 iϭN ͩ aϪ1 aϩ1 p i ͑ t i ͒Ϫ p exp r ͑ t i ͒ ͪ 2 ϭ0, ͑23͒ then 2 ͚ iϭ1 iϭN ͑ aϩ1 ͒Ϫ͑ aϪ1 ͒ ͑ aϩ1 ͒ 2 p i ͑ t i ͒ ͩ aϪ1 aϩ1 p i ͑ t i ͒Ϫ p exp r ͑ t i ͒ ͪ ϭ0, ͑24͒ and ͚ iϭ1 iϭN p i ͑ t i ͒ ͩ aϪ1 aϩ1 p i ͑ t i ͒Ϫ p exp r ͑ t i ͒ ͪ ϭ0,

͑25͒

which is equivalent to

aϭ ͚ iϭ1 iϭN p i ͑ t i ͒͑ p i ͑ t i ͒ϩ p exp r ͑ t i ͒͒ ͚ iϭ1 iϭN p i ͑ t i ͒͑ p i ͑ t i ͒Ϫ p exp r ͑ t i ͒͒ ϭ ͱ␣ ϱ , ͑26͒
the expression of the porosity is

ϭ ͚ iϭ1 iϭN p i ͑ t i ͒͑ p i ͑ t i ͒Ϫ p exp r ͑ t i ͒͒ ͚ iϭ1 iϭN p i ͑ t i ͒͑ p i ͑ t i ͒ϩ p exp r ͑ t i ͒͒ ͱa ϱ ,

͑27͒

the second derivative term

ץ 2 U ץa 2 ϭϪ 8 ͑ aϩ1 ͒ 3 ͚ iϭ1 iϭN ͩ aϪ1 aϩ1 ͑ p i ͑ t i ͒͒ 2 Ϫ p exp r ͑ t i ͒. p i ͑ t i ͒ ͪ ϩ 8 ͑ aϩ1 ͒ 4 ͚ iϭ1 iϭN ͑ p i ͑ t i ͒͒ 2 , ͑28͒
from Eq. ͑24͒ we have

ץ 2 U ץa 2 ϭϪ 2 ͑ aϩ1 ͒ ץU ץa ϩ 8 ͑ aϩ1 ͒ 4 ͚ iϭ1 iϭN ͑ p i ͑ t i ͒͒ 2 , ͑29͒
when ץU/ץaϭ0, the condition ץ 2 U/ץa 2 Ͼ0 is always verified because aϾ1, (␣ ϱ Ͼ1 and Ͻ1).

Relation ͑27͒ is the expression of the porosity as a function of the tortuosity, a discrete set of values of the incident and reflected signal, we will use this relation in Sec. V to estimate the value of the porosity from experimental incident and reflected data.

V. ULTRASONIC MEASUREMENTS

In application of this model, some numerical simulations are compared to experimental results. Experiments are carried out in air with a broadband Ultran NCT202 transducer having a 190 kHz central frequency in air and a bandwidth at 6 dB extending from 150 to 230 kHz. This transducer is used simultaneously as a transmitter and receiver in order to detect the reflected wave. Pulses of 400 V are provided by a 5052PR Panametrics pulser/receiver. The received signals are amplified up to 90 dB and filtered above 1 MHz to avoid high frequency noise ͑energy is totally filtered by the sample in this upper frequency domain in transmitted mode but not in reflected mode͒. Electronic perturbations are removed by 1000 acquisition averages. The experimental setup is showed in Fig. 5. The duration of the signal plays an important role in solving the inverse problem, its spectrum must verify the condition of high frequency approximation referred to in the previous section.

Figure 6 shows the incident signal given out by the transducer and the reflected signal by the plastic foam M 2 ͑thickness 4.1 cm, ␣ ϱ ϭ1.07Ϯ0.01, ⌳ϭ230 m, and ⌳Ј ϭ690 m). Figure 7 shows their spectrums. From the spectrums of the two signals, the reader can see that they have the same bandwidths which means that there is no dispersion. This concurs the theory below which predicts that the reflected wave from the first interface xϭ0 is measured and simply attenuated by the factor (ͱ␣ ϱ Ϫ)/(ͱ␣ ϱ ϩ).

Using the incident and reflected experimental data for the foam M 2 and relation ͑27͒, the value of the porosity optimized by solving the inverse problem is ϭ0.968 Ϯ0.005. This value of the porosity is obtained on the basis of an average of a set of values, and its accuracy is calculated knowing the accuracy of the tortuosity using the relation ␦/ϭ(1/2)(␦␣ ϱ /␣ ϱ ), where ␦ is the accuracy of the porosity and ␦␣ ϱ is the accuracy of the tortuosity.

Figure 8 shows the comparison between experimental reflected signals and simulated reflected signals for the optimized value of the porosity, the difference between the two curves is slight which leads us to conclude that the optimized value of the porosity is good. The value of the porosity of the foam M 2 given by the porosimeter 10,18 is ϭ0.97Ϯ0.02. The flow resistivity of the plastic foam M 2 is ϭ3000 N m Ϫ4 s. This value is low compared with other plastics foams, and is the reason why the value of the porosity is very close to 1.

Another plastic foam M 3 is studied ͑thickness 0.6 cm͒. This has a high level of flow resistivity, ϭ125 000 N m Ϫ4 s and the ultrasonic parameters are ␣ ϱ ϭ1.70Ϯ0.01, ⌳ϭ23 m, ⌳Јϭ69 m. Figure 9 shows the incident and reflected signal and Fig. 10 shows their spectra respectively. Here again, we can conclude from the two spectra that there is no dispersion. By solving the inverse problem from the incident and reflected data and relation ͑27͒ the optimized value of the porosity is ϭ0.821Ϯ0.005. Figure 11 shows the comparison between simulated reflected signal obtained by the optimized value of the porosity and the experimental signal. The readers can see the slight difference between the two curves that proves the optimized value of the porosity.

The main advantage of this method is its simplicity compared with the porosimeter, 18 which involves more complicated measurement steps and a more expensive experimental setup.

VI. CONCLUSION

A simple method for the measurement of the porosity has been proposed. This method is based on a concept of The method is based upon the propagation equation in the time domain in a slab of porous material having a rigid frame at high frequency range. A time domain model of wave propagation in such material is worked out from the concept of fractional calculus. The kernels of the reflection and transmission scattering operators are derived giving simple relations between these operators and the parameters of the medium. Because of the high attenuation in such media, the wave reflected by the porous material can be approximated by the wave reflected at the first interface. As a consequence, it gives a simple relation between the porosity, tortuosity and reflected wave, thereby proving that all physical parameters of acoustic material can be evaluated using acoustic measurements.

An experimental validation of the theoretical expression of the reflection scattering operators illustrates the high level of correspondence between numerical and experimental results and shows that this time domain model is well suited for the characterization of porous media via acoustic wave propagation. As part of our future plan, we hope to deal with the problem of oblique incidence that can give us more data for the reflected mode of the inverse problem.

͑15͒FIG. 1 .

 1 FIG. 1. Geometry of the problem.

͑16͒ and the kernel of reflection by R ˜͑t ͒ϭr͑ t ͒ϩR͑ t

FIG. 2 .

 2 FIG. 2. Transmitted simulated signals for ϭ0.98 ͑solid line͒ and for ϭ0.49 ͑dashed line͒.

FIG. 4 .

 4 FIG. 4. Reflected wave at the interface xϭ0 ͑solid line͒ and the total reflected wave ͑dashed line͒.

FIG. 6 .

 6 FIG. 6. Experimental incident signal ͑dashed line͒ and experimental reflected signal ͑solid line͒.

FIG. 8 .

 8 FIG. 8. Comparison between experimental reflected signal ͑solid line͒ and simulated reflected signal for ϭ0.968 ͑dashed line͒.

FIG. 10 .

 10 FIG. 10. Spectrum of experimental incident signal ͑dashed line͒ and spectrum of experimental reflected signal ͑solid line͒.

J. Appl. Phys., Vol. 93, No. 1, 1 January 2003 Fellah et al. Downloaded 18 Jul 2003 to 194.57.163.19. Redistribution subject to AIP license or copyright, see http://ojps.aip.org/japo/japcr.jsp

ACKNOWLEDGMENTS

The authors would like to thank Professor Nora Berrah ͑Physics Department, Western Michigan University͒ and Dr. Smaine Zeroug ͑Ph.D. Program Manager, RT Modeling and Inversion Applications Schlumberger-Doll Research, Ridge-field͒ for their helpful comments.

APPENDIX: GREEN FUNCTION OF THE MEDIUM

The Green function of the propagation equation ͑7͒ is given in Ref. 3 by

where h(,) has the following form:

and where ͑,, ͒ϭ͓⌬ͱ͑ Ϫ͒ 2 Ϫk 2 ϩbЈ͑Ϫ ͔͒ 2 /8, bЈϭBc 0 2 ͱ,c ЈϭC.c 0 2 and ⌬ϭbЈ 2 Ϫ4cЈ when k→ϱ, ⌶, and h(,) tends to zero, then the Green function F(t,k) tends also to zero.