A refined simulation method of building earthquake evacuation processes considering multi-exits and time-variant velocities

This study proposes a refined method for effectively simulating the evacuation processes of multi-story buildings during earthquakes. The method is developed by combining the exit-option algorithm and multi-velocity field model with the previously developed probabilistic network-based floor-field cellular automaton (PN-FFCA) model. Besides, the formula is also proposed for calculating the time-variant moving velocities of evacuees in the buildings subjected to the earthquake loads. Consequently, the refined method can fully consider the effect of multi-exits and time-variant velocities on the random evacuation processes of multi-story buildings during earthquakes. The recorded evacuation process of a school classroom during the 2022 Luding earthquake in China is reproduced using the refined method, which verifies its validity. In addition, the evacuation processes of an actual three-story office building subjected to random earthquake loads are simulated using the refined method. The simulated results are compared with those obtained from the PN-FFCA to demonstrate the advantages and utility of the refined method.


Introduction
Accurate simulation of building earthquake evacuation processes is significant for modern assessments of structural seismic performance.The simulated results can provide information on the total earthquake evacuation times, pedestrian casualty rates, pedestrian evacuation paths, and so on, of buildings during earthquakes.This information is helpful for urban disaster prevention and performance-based seismic design of buildings.
A building earthquake evacuation simulation method mainly consists of the stochastic models of emergency evacuation simulation, interaction analysis of evacuees' behavior, and structural damage or collapse, as well as structural seismic damage or collapse process analysis.In the last few decades, numerous macroscopic and microscopic emergency evacuation simulation models have been developed, including the fluid dynamics model (Franjo et al., 2022;Hoogendoorn and Bovy, 2000;Wouter et al., 2021), social force model (Helbing et al., 2000;Johansson et al., 2015;Liu, 2018), cellular automaton (CA) model (Hartmann and Hasel, 2014;Ji et al., 2018;Yang et al., 2005;Yue et al., 2007), lattice-gas model (Huo et al., 2014;Khain et al., 2014;Li et al., 2008), game theory model (Lo et al., 2006), agent-based model (Bao and Huo, 2021;Cui et al., 2021;Delcea and Cotfas, 2019), and so on.The macroscopic model, that is, the fluid dynamics model, treats the crowd as a whole and has high computational efficiency.However, it cannot reflect the interaction and behavioral differences between individuals.In contrast, the microscopic models, including the social force, CA, lattice-gas, game theory, and agent-based models, evaluate each pedestrian individually and concentrate on the interaction between pedestrians.Therefore, they can give relatively accurate simulation results.Among them, due to the ability to analyze the complex behavior of pedestrians and high computational efficiency, the CA model has been widely used in practical engineering (Chen et al., 2021;Hulse et al., 2020;Zheng et al., 2019) and commercial software, such as STEPS, FDS + Evac, and Building EXODUS.To effectively analyze the interaction between pedestrians during the evacuation, a two-dimensional floor-field cellular automaton (FFCA) model was proposed by introducing the dynamical field computing local interaction into the conventional CA model (Burstedde et al., 2001;Kirchner and Schadschneider, 2002).Based on the FFCA model, a probabilistic network-based floor-field cellular automaton (PN-FFCA) model was recently developed for solving the three-dimensional evacuation simulation problem, especially for simulating the emergency evacuation processes of multi-story buildings (He, 2020).Although the PN-FFCA model exhibits relatively high flexibility and effectiveness in simulating the building earthquake evacuation processes, it does not deeply consider the problems of dynamic selection of the target exit and real-time change of evacuees' moving speed during the evacuation processes.
The interaction between evacuees' behavior and structural seismic damage or collapse substantially impacts the earthquake evacuation processes.Li et al. (2015Li et al. ( , 2018) ) proposed a numerical simulation method using an improved CA model to evaluate the casualties under earthquake.This method uses relative displacement to define casualty occurrence criteria.However, the simulation does not consider the impact of the non-structural and structural damage process on the evacuation process.Liu et al. (2016) and Cimellaro et al. (2017) used agent-based models to simulate pedestrian evacuation from damaged structures after earthquakes.They considered the interaction between the agent and the external environment damaged by the earthquake, especially the structural and non-structural damaged components in the building.However, they focused on the simulation of building evacuation after the earthquake.Poulos et al. (2017) introduced a conceptual framework and numerical algorithm to assess the seismic risk of building occupants during earthquake events.The framework combines probabilistic seismic hazard analysis, inelastic structural response analysis, damage assessment, and the response of evacuating agents.Lu et al. (2019Lu et al. ( , 2020) ) proposed a pedestrian evacuation simulation framework considering falling debris caused by earthquakes, focusing on the impact of the collapse of the innerfilled wall on the earthquake evacuation process.He (2020) proposed an entirely random evacuation model for simulating and analyzing the earthquake evacuation process of multi-story buildings during the earthquake.The model coupled the emergency evacuation and structural damage processes but did not consider the pedestrian exit selection and velocity variance.
The structural seismic damage or collapse process analysis is a key issue of the earthquake evacuation simulation and analysis since it provides the environmental information about the evacuation area that significantly influences pedestrians' behavior.Various effective analysis methods of structural seismic damage and collapse have been developed in the past decades.Han and Chopra (2006) and Moon et al. (2012) introduced the modal pushover analysis method into the incremental dynamic analysis (IDA) method for effectively analyzing structural seismic damage and collapse.Zarfam and Mofid (2011) assumed the material constitutive relationship as a three-line model and used the IDA method based on modal pushover analysis to conduct seismic damage analysis of reinforced concrete structures.The fine analysis of structural seismic damage and collapse is usually required by the finite element method.Kaewkulchai and Williamson (2004) and Lynn and Isobe (2007) established the finite element method for continuous collapse analysis of structures.Munjiza et al. (2004) proposed a hybrid discrete element method for structural failure and collapse analysis.Bao and Kunnath (2010) proposed a simplified finite element analysis method for the continuous collapse analysis of reinforced concrete frame shear wall structures.Li et al. (2015Li et al. ( , 2018) ) used the finite element software Abaqus to simulate building collapse.In their study, it is assumed that the infilled wall will collapse out of plane when the peak acceleration at the center of the wall reaches 1.0 g.In their earthquake evacuation simulation, Liu et al. (2016) used Multiphysics simulation software LS-DYNA to analyze the damage processes of non-structural components.Cimellaro et al. (2017) and Poulos et al. (2017) used OpenSEES software to obtain the structural seismic damage state during pedestrians' evacuation.In their earthquake evacuation study, Lu et al. (2019Lu et al. ( , 2020) ) used the nonlinear multi-degree-of-freedom models of the multi-story frame and masonry structures, and tall buildings to analyze the collapse processes of masonry-infilled walls.He (2020) simplified the multi-story reinforced concrete frame structure into a multi-degreeof-freedom shear model, from which the structural seismic damage and collapse state during the earthquake evacuation can be analyzed efficiently.It should be noted that effective analysis of structural seismic damage and collapse processes, and distribution of the debris of structural and non-structural components is extremely important for improving the efficiency of earthquake evacuation simulation.
As mentioned above, the building earthquake evacuation simulation is a relatively complex scientific and engineering problem that involves random earthquakes, and human and structural behavior.Although several methods, such as the PN-FFCA method, have been developed in the past decade, some critical issues still need to be solved to improve their effectiveness and flexibility.In this study, a refined simulation method is developed by combining the exit-option algorithm and multi-velocity field model with the PN-FFCA method, and by proposing the formula for calculating the time-variant moving velocities of evacuees in the buildings during earthquakes.The effectiveness, advantage, and use of the proposed method are demonstrated by two examples, that is, the scenario recurrence analysis of earthquake evacuation of a school classroom during the 2022 Luding earthquake in China and the comparative analysis of earthquake evacuation of an actual threestory office building subjected to random earthquake loads.

Review of PN-FFCA method
An essential technique in the PN-FFCA method is to model the whole network for the considered multi-story building according to the discrete grid of evacuation area and the type of neighborhood used by the CA model.The particles representing pedestrians can randomly move on the nodes of the constructed whole network such that the earthquake evacuation process in the considered building can be effectively simulated.Figure 1 shows the connection of the whole network between two adjacent floors when the type of Moore neighborhood (Burstedde et al., 2001) is adopted.In the PN-FFCA method, the movement of a particle is determined by a vector of transition probabilities P i = ½p i, 0 , p i, 1 , :::, p i, ki , where p i, 0 represents the probability that the particle will remain in the original node and p i, j , j = 1, :::, k i represents the probability that the particle will move into a neighbor node, and k i represents the number of neighbor nodes of the considered node i.The calculation formula of the transition probability p i, j can be derived from the FFCA model: where in which t s and t d indicate the static floor-field values and dynamic floor-field values, respectively.Table 1 summarizes the typical parameters used in Equation 1 as follows.The static floor-field value is inversely proportional to the distance between node i and the nodes denoting the exits.Therefore, the static floor field t s (k) of each node k can be calculated from: where L represents the total number of exit nodes, K represents the total number of nodes, D i!l is the shortest distance from a node i = 1, :::, K to an exit node l, and D k!l is the shortest distance from node k to an exit node l, where i runs over all nodes and l runs over all exit nodes.
The dynamic floor field is a virtual trace left by the pedestrians and has its own dynamics, that is, diffusion and decay (with decay constant a 2 ½0, 1).It is used to model an attractive interaction between the particles.The calculation for the dynamic floor field t d is as follows.At the beginning of the simulation, the dynamic floor-field value t d for each node is set to a given value.As the pedestrian moves from a node to its neighboring node, the dynamic floor-field value t d of the considered node is increased by one.After all particles' movements in a time step have been completed, the dynamic field value t d of each node is reduced by one with probability a if it was created at the previous update step or earlier, where a 2 ½0, 1.The parameter a is set to be a = 0:3 in this study (Kirchner and Schadschneider, 2002).
The parameter n j in Equation 1 reflects the number of obstacles on the target node j, with n j = 1 if the corresponding node j is occupied and n j = 0 otherwise.It should be noted that the parameter n j should be modified by the structural damage processes.The parameter N in Equation 1 serves as a normalizing factor to make sure P ki j = 0 p i, j = 1.The variables J s and J d determine the coupling strength between the pedestrian and the static and dynamic floor field.The value of J s can be viewed as a measure of the knowledge of the pedestrians about the location of the exit, and the parameter J d for the coupling to the dynamic field controls the tendency to follow the lead of others.The above two parameters can be determined through real-life video verification or earthquake evacuation tests.
The correction parameter d i, j in Equation 1 is an inertia factor to avoid particles being confused by their generated paths.If the present movement direction is the same as the previous one, this movement is strengthened; otherwise, it is lessened.The parameter d i, j The coupling strength between the pedestrian and static floor field The coupling strength between the pedestrian and dynamic floor field n j Debris identification parameter d i, j Inertia factor can be calculated as follows: If direction vector[(position(t + 1) 2 position(t)] = direction vector [(position(t) 2 position(t 2 1)], setting: If direction vector [(position(t + 1) 2 position(t)] = direction vector [(position(t 2 1) 2 position(t)], setting: In other cases, Based on the vector of transition probabilities P i = ½p i, 0 , p i, 1 , :::, p i, ki , the updating and collision rules of classical CA models can be adopted to simulate the particles' movement on the whole network.For more details about the PN-FFCA method, the readers are referred to Burstedde et al. (2001) and He (2020).

Effect of multi-exits and time-variant velocities
As mentioned above, the PN-FFCA method does not fully consider the effects of multiexits, time-variant velocities, and psychological behavior of evacuees on the earthquake evacuation processes.To overcome the limitations, the PN-FFCA method is refined in this study by considering the above effects, which are described as follows.

Effect of multi-exits
For the buildings with multi-exits, at each time step in the PN-FFCA method, the particles need to select the preferred exit first and then the static floor fields in Equation 1 are calculated.Without loss of generality, consider the situation shown in Figure 2. It is assumed that there are two exits on the kth floor of a multi-story building.To select the preferred exit at a time step, the pedestrians must analyze the desirability of all exits.To do so, the concepts of position attraction and density attraction can be introduced.The position attraction is determined by the distance from the pedestrian to the exit, that is, the shorter the distance, the greater the value of that exit.For typical Exit i (i = 1, Á Á Á N), the preferred probability defining position attraction can be given by (Yuan and Hai, 2007) where l i is the distance from the pedestrian to Exit i and L = The parameter a is a scalar constant adjusting the sensitivity of l i .It is obvious that From Equation 8, it can be seen that the greater the distance between the particle and exit i, the lower the probability of choosing exit i.
The density attraction is determined by the number of pedestrians in the exit area (EA), depicted in green in Figure 2. The lower the density around the exit, the more significant the amount of the density attraction.Similar to Equation 8, the preferred probability of density attraction can be given by (Yuan and Hai, 2007): where d i is the occupied density of Exit i and The parameter b is a scalar constant adjusting the sensitivity of d i .
Note that even though position and density attraction are computationally independent, occupants are evaluated simultaneously when selecting an exit.Consequently, based on the above equations, the preferred probability P i for Exit i is: and where c and d are two constant scalars.After the particle completes the exit selection according to Equation 10, the refined PN-FFCA method employs the static floor field of the selected exit to calculate the transition probability P i, j , thus completing the particle movement.
Adding exit selection to the transition probability improves the pedestrian's intelligence in evacuation.The selection of exits balances the number of pedestrians moving through each exit and prevents excessive pressure on a particular exit, enhancing evacuation efficiency.In addition, it can guide the reasonableness of the exit settings.

Effect of time-variant velocities
To consider the effect of time-variant velocities on the earthquake evacuation processes, the multi-velocity field model (Yuan and Hai, 2007) is introduced into the PN-FFCA method.The multi-velocity field model requires information about the desired velocities of the particles at each time step, which can be obtained by the following formula (Wirz et al., 2013): where v f (i) is the initial velocity of particle i, and v d (i) is the desired velocity of the particle i. r(i) is the density of the crowd around the particle, and r max = 6:25m À2 is the maximum crowd density that was reached in this method.It is worth noting that Equation 13 can be further modified by considering the effects related to structural damage, such as panic from the collapse of structural components or the passage of pedestrians through structural debris, as described in section ''Calculation of time-variant velocities.'' Assuming that v max is the maximum value that v d (i) can achieve, the velocity ratio of pedestrian i is defined as'' A time step can be calculated by: where Dl represents the length of a cell of the discrete grid of the evacuation area, for example, Dl = 0:4m.
The updating rule of the PN-FFCA method should be corrected as follows: first, generate a pseudo-random series u = (u 1 , u 2 :::u i :::u n ) at each time step, where u i 2 ½0, 1 is the pseudo-random number corresponding to particle i; second, compare the magnitude of h i and u i .If h i ø u i , the particle can move on the whole network according to the transition probabilities defined by Equation 1. Otherwise, the particle will stay in the original node of the network.
By introducing the multi-velocity field model, the time-variant velocity-based judgment is used to simulate the pedestrians' movement, and the pedestrians with extremely low velocities may stay in their original positions for multiple time steps.
Another effect of the time-variant velocity is the change of the dynamic floor fields with the time-variant velocities, that is, the pedestrians usually tend to move to an area with a higher movement velocity so as to get to the selected exit as soon as possible (You et al., 2016).Therefore, fast pedestrians should be more attractive.For the original PN-FFCA method, the dynamic floor-field values generated by different pedestrians are the same, which makes it difficult to simulate the willingness of pedestrians to move to areas with higher moving velocity.To solve this problem, a dynamic floor-field increment rule is proposed.In this rule, when a pedestrian moves from node i to its neighboring node j, the dynamic floor-field increment Dt d of node i is calculated by: while the decay value of the dynamic floor field À Dt d caused by the considered particle in the subsequent time steps is calculated by: where v d (i) represents the velocity of the particle at node i, and a represents the velocity sensitivity parameter.Using this rule, high-velocity pedestrians generate larger numerical increments of the dynamic field, and the dynamic field dissipates slowly, making it more attractive for other pedestrians to follow.
To solve the conflict resolution among particles when numerous individuals compete for the same node on the whole network, it is assumed that only one node is permitted to enter the conflicted node, while all other nodes must remain in their original nodes.When numerous particles compete for the same node, it is understood that the speedier individual has an advantage.Consequently, the probability that the particle on node i will prevail in conflict can be estimated by: where v i is the velocity of the particle on node i in the competitive time step and a is a velocity sensitive parameter that indicates the competition's average capacity.To simplify, the value of a can be set to 1 (Luo et al., 2018).

Calculation of time-variant velocities
Three main factors in real-time affecting pedestrians' velocity in the earthquake evacuation processes are considered in this study, including panic caused by collapsing components, pedestrian health conditions, and the debris inside the buildings.Therefore, the new desired velocity of the ith pedestrian (particle) in cell j can be modified using the following equation: judgment is made to update the velocity modification coefficients.If the non-structural component collapses in cell j, the pedestrian state will change from normal to panic and the f i, j r will vary from 1 to 1.2.If cell j is covered with debris, the pedestrian's velocity will be reduced as follows: f i, j d = 0 for the collapsed structural components, f i, j d = 0.5 for the partition walls, and f i, j d = 0.7 for the suspended ceilings.If the pedestrian is hit by the collapsed component in cell j, the velocity decrease is as follows: f i, j h = 0.75 for the suspended ceilings, f i, j h = 0.5 for the partition walls, f i, j h = 0 for the structural components, and if none of the above three cases, f i, j h = 1.The values' choice of the coefficients f i, j s f i, j r f i, j d and f i, j h is referred to Poulos et al. (2017).Of course, more reasonable values of the coefficients should be identified by the corresponding experiment results.

Coupling between pedestrians' evacuation and structural damage or collapse
There are two critical issues in the coupling of pedestrian evacuation and structural damage, one is the probabilistic assessment of damage to structures and the other is the occurrence times of damage to structures.
According to the performance-based seismic design principles of structures, the probability of non-structural damage, structural damage, and structural collapse of a multistory building is equivalent to the probability of the corresponding inter-story drift ratio response exceeding a specified threshold within the reference time.The thresholds of the reinforced concrete frame structures can be determined by some performance-based seismic design codes; specifically, the threshold levels for a three-story reinforced concrete frame structure can be determined by the standard CECS 160-2004 (2004), as shown in Table 2 in Appendix.The percentage of structural and non-structural component damage should be determined by the seismic fragility analysis.In this study, the percentages of damage to structural and non-structural components of a three-story reinforced concrete structure are determined according to the study by Cimellaro et al. (2017) and He, (2020) and listed in Tables 3 and 4 in Appendix.Note that it is assumed that the percentage of damaged partition walls is one-third of that of the suspended ceiling in Table 3.
The occurrence times of seismic damage to structures are uncertain due to the randomness of earthquakes and structural characteristics.This study assumes that the structural damage is brittle, and the occurrence times are defined as the time of the first passage of the inter-story drift ratios exceeding the threshold values corresponding to non-structural damage, structural damage, and structural collapse.Some non-structural components, such as suspended ceilings, may be acceleration-sensitive.However, considering the current research status of the seismic damage of non-structural components, only the inter-story drift ratio is considered in the non-structural component failure analysis in this study.The inter-story drift ratio samples can be obtained by the Monte Carlo simulation method (He and Gong, 2016).The procedure of the Monte Carlo method is as follows: first, a seismic excitation sample is generated based on the specified random earthquake load model; second, the deterministically dynamic analysis of the considered structure (or its simplified analysis model) subjected to the seismic excitation sample is conducted to obtain an interstory drift ratio response sample; third, the above two steps are repeated until the total number of N inter-story drift ratio samples have been obtained.
Due to the randomness of the pedestrians' evacuation processes and structural seismic damage and collapse processes, their coupling should be performed based on the sampling.Specifically, the coupling procedure is as follows: (1) establish a stochastic model for seismic loading and generate a seismic loading sample; (2) calculate an inter-story drift ratio response sample using a simplified structural analysis model; (3) generate a structural damage or collapse process sample (including the occurrence time and distribution) according to the above calculation method; (4) consider using the structural damage or collapse sample to simulate an evacuation process for pedestrians; and (5) repeat Steps (1)-( 4) N times and conduct statistical analysis on the N simulation results.

Illustrative examples
Two examples are presented to demonstrate the accuracy and effectiveness of the developed method.To perform the analysis, the mathematical software MATHEMATICA 11 is used for programming.The codes can be downloaded from the website (https://pan.baidu.com/s/1vd513oJaHue8fsdgEB9-iA,Extract code: 0bvf).The analyses are conducted on a computer with an Intel I7 CPU and 16GB RAM.

Classroom earthquake evacuation simulation
This example is used to demonstrate the validity of the proposed method by comparing the simulation results with the real-life video recording of the emergency evacuation process of a middle school classroom hit by an earthquake.At 12:52 on 5 September 2022, a magnitude 6.8 earthquake occurred in Luding County, Ganzi Prefecture, Sichuan Province, China.The monitor clearly recorded the pedestrians' evacuation process of a classroom during the earthquake.Both the personnel evacuation trajectory and evacuation time were obtained from the record.
The classroom layout and student evacuation information are shown in Figure 3.As illustrated in Figure 3, the dimension of the classroom is 6:4m38:4m and the width of both exits is 1.2 m.There were 38 students and 1 teacher in the classroom when the earthquake occurred.The distance between two adjacent tables is 0.8 m.In Figure 3a, each particle (student) is labeled and the number on the label represents the evacuation time itself.There are four evacuation routes in total, which are indicated by four different colored arrows in Figure 3b.The model parameters reflect the behavioral characteristics of personnel in the evacuation process, and whether the model parameters are set reasonably directly affects the evacuation simulation results.In this example, the parameter calibration is performed as follows.First, parametric modeling is carried out to explore the influence of the parameter change on the total evacuation time.Then, based on the real-life video recording of classroom evacuation, analyze the video frame by frame to obtain the corresponding evacuation curve.Finally, identify the parameter values so that the simulated evacuation process best fits that from the evacuation video.After calibration, the key parameters of the refined PN-FFCA method are identified to be Js = 2:0, Jd = 1:0, b = 10, a = 1:0, and b = 3:0.It is worth noting that the students are very familiar with the evacuation scenario and in good health and trained to evacuate; therefore, the parameters calibrated for this scenario may be not suitable for evacuation behavior in other building types with greater diverse occupants and a lack of training.
The discrete grid cell is set to 0:4m 3 0:4m.Based on the grid cell size, the whole network for the evacuation simulation is established, as shown in Figure 4.
In this example, the maximum pedestrian density r max = 6:25m À2 , which is obtained from theoretical calculations based on the CA model.The initial velocity v f (i) is set to be 3 m/s based on the running velocity of the students in the classroom.Figure 5 shows the variation of the desired velocities v d (i) of the students obtained from the original and refined PN-FFCA method.In Figure 5, each line represents the velocity variation of a student during the evacuation.From Figure 5a, it can be observed that the desired velocities obtained from the original PN-FFCA method remain constant throughout the whole evacuation process, which is inconsistent with the actual evacuation.While for the refined PN-FFCA method, the desired velocities change at all times during the evacuation process, as shown in Figure 5b.Moreover, during the evacuation process, the minimum value of desired velocities among the students decreases first and then increases.The main reason is that students are uniformly distributed in the classroom before the evacuation begins, and the density around them is low, resulting in fast initial velocities.After the evacuation begins, students gather toward the exit, causing a sudden increase in density and a sharp decrease in the minimum desired velocity.Subsequently, as the number of students in the classroom decreases, the density around them begins to fall, and the overall velocities begin to increase.From the above analysis, it is clear that the refined PN-FFCA method has advantages compared to the original PN-FFCA method for treating timevarying velocities.
Figure 6 shows the distribution of students' positions at t = 2s, t = 4s, and t = 6s obtained, respectively, from the real-life record, original PN-FFCA method, and refined PN-FFCA method.It can be found in Figure 6 that the simulation results obtained by the refined PN-FFCA method are better than the original PN-FFCA method.Figure 6a shows that in the early stage of evacuation (t = 2s), students obeyed the teacher's command and began to evacuate with their hands over their heads.Students near the desks quickly rushed into the aisle between desks.It can be found that the aisle near the exit attracted more students to choose due to its location advantage, which is reflected in both the real video and simulated evacuation process.Figure 6b shows that in the middle stage of evacuation (t = 4s), students spontaneously divided into two groups and evacuated toward the exit close to them.At the same time, it was clear that many students were still occupying the aisle next to the exit.Finally, as seen in Figure 6c, the number of students had significantly decreased, and the evacuation process was complete.From Figure 6a and b, it is evident that the evacuation process obtained by the original PN-FFCA method is significantly faster than those obtained from the recorded and refined PN-FFCA method.However, in the second half of the evacuation, see Figure 6c, the evacuation process obtained by the original PN-FFCA method significantly slowed down.The main reason for this phenomenon is that the algorithms of velocity in the original and refined PN- FFCA methods are different.In fact, in the first half of the real evacuation (t\ 4 s), due to students' gathering toward the exit, the density around the students continued to increase and the velocity decreased, resulting in a slow evacuation process.This can be well simulated in the refined PN-FFCA method because the time-variant velocities take into account the impact of personnel density.In the original PN-FFCA method, the velocity is always constant, so the evacuation process is significantly faster than those obtained from the recorded and refined PN-FFCA method.In the second half of the evacuation (t ø 4 s), as the number of students in the classroom decreases, the students' velocity starts to increase.However, the velocity remains constant in the original PN-FFCA method and, consequently, the original PN-FFCA method gives slower evacuation progress.
Figure 7 shows the comparison between the recorded evacuation curve and the simulating evacuation curves, which further confirms the differences above mentioned.It is also observed from Figure 7 that, for the original PN-FFCA method, the last several escaped students spent more time moving through the area close to the exit.The reason causing the phenomenon is that these students move around themselves in the last several time steps until their dynamic potential completely decays and then they reach the exit.The refined PN-FFCA method solves this problem due to the multi-exit selection algorithm and improved dynamic field rules.
The above comparative analysis demonstrates the validity of the refined PN-FFCA method for simulating the earthquake evacuation processes of buildings with multiple exits and time-variant velocities.It is also worth noting that this example does not consider the effect of structural damage or collapse on the evacuation process because there was no significant damage or collapse in the classroom during the earthquake.Regarding the calibration of model parameters, the parameter values should be different for different evacuation scenarios.The parameters calibrated in this example are only applicable to the evacuees who are physically fit, extremely familiar with the evacuation building layout, and have received relevant evacuation training.For the determination of parameters in new scenarios, it needs to determine the behavioral characteristics of the personnel through earthquake evacuation tests and questionnaires.

Official building earthquake evacuation simulation
This example is to demonstrate the effectiveness and use of the refined PN-FFCA method by simulating and comparatively analyzing the earthquake evacuation process of an actual office building subjected to a random earthquake load.
Case scenario.The considered building is a three-story office building in use, of which the first and second stories are occupied by personnel.On weekdays, there is 60 staff in the building, including 33 staff on the first story and 27 on the second.The type of construction is a reinforced concrete frame with the concrete strength grade of C30.The height of the building is 4.2, 3.6, and 3.6 m from the first to the third story, respectively.The general framework of the building is shown in Figure 16a in Appendix.Figure 8 shows the floor plans of the first and second floors, where the rectangular boxes with thin lines represent desks, bookcases, and so on.
The whole network and staff distribution of the first and second stories are shown in Figure 9, where the red circles represent office staff, the black area represents walls, tables, and other obstacles, and the yellow areas represent stairs.Two evacuation exits are considered, one on the X axis of the first floor and the other on the Y axis of the first floor.The width of both exits is 1.6 m.In the discrete grid modeling of the evacuation area, the size of each cell is (0:4 3 0:4) m 2 .The modeling process of the whole network is as follows: first, project the upper surfaces of the stairs and platform between adjacent floors vertically onto the plane of the upper floor slab; second, discrete the evacuation areas and the vertical projection of the stairs and platform according to the cell size; third, choose Moore neighborhood to determine the neighboring cells that allow the particles in the central cell to move into and modify the Moore neighborhood of the stairs to establish the connection between the first and second floors to realize the movement of pedestrians between different stories.Finally, set the centers of cells as vertices and the lines between the centers of the neighboring cells that allow each other to move as edges to construct the undirected network.
Simulation and comparative analysis.Because of the randomness of earthquake and pedestrian evacuation as well as the coupling of pedestrian evacuation and structural damage samples, it usually needs a large number of runs of structural seismic response and damage simulation in the earthquake evacuation analysis.Hence, the simplified structural models should be adopted to simulate the seismic response and damage samples to save computational costs.For the reinforced concrete frame structure considered in this example, a simplified structural model, that is, the shear model, has been established (He, 2020), as shown in Figure 16 in Appendix.For the shear model, Rayleigh damping with a 5% damping ratio in Modes 1 and 2 is assumed.The hysteretic restoring force of each story of the structure is defined as follows: where k 0 is stiffness; a represents the nonlinear parameter value of 0.1.DX and Z represent the elastic and hysteretic components of inter-story displacement, respectively.Z follows the well-known Bouc-Wen hysteresis law: in which parameter values n = 2, A = 1, and g = b = A=(2x n y ), where x y = 0:01m is the yield displacement.To consider the randomness of earthquakes, the time-modulated Gaussian process is used to model the earthquake load.Specifically, the earthquake excitation is written as , where m(t) = 0:093t 3 exp(20:5t) is the modulating function and € U (t) is a stationary Gaussian process with zero mean and the well-known Kanai-Tajimi power spectrum density: in which S 0 = 0.03 m 2 /s 3 , v f = 15.7 rad/s, and z f = 0.6.
A total of 100 seismic response samples of the structure are calculated, from which the first occurrence times of structural and non-structural damage can be determined according to the method described in section ''Coupling between pedestrians' evacuation and structural damage or collapse.''Subsequently, the collapse distribution of structural components and non-structural components is generated by random sampling.For multi-story reinforced concrete frame structures, the coupling between structural damage and non- structural component damage is not as obvious as in tall structures (International Code Council Evaluation Service, 2010;Lu et al., 2017).Therefore, to simplify the structural damage analysis, the collapse distributions of different components are generated independently in this study.Figure 10 shows the collapse distribution of structural and nonstructural components for one sample, where the orange squares represent the fall location of damaged suspended ceiling components, the blue squares represent the fall location of damaged partition wall components, and the red squares represent the fall locations of structural components.To generate the distribution of damaged non-structural components, first randomly select them based on the percentage of damaged non-structural components and then vertically project them on the floor.To generate the distribution of damaged structural components, it is assumed that the damaged floor slabs, beams, and columns of the building uniformly scatter on the floor just below them.Therefore, the distribution of the damaged structural components is generated by randomly selecting the cells on the corresponding floors according to the area percentage of the damaged structural components.Note that, considering the small proportion of collapse of structural components, the overall failure of structural components is not considered in this study.In the earthquake evacuation simulation, once a cell is covered by debris of damaged structural and non-structural components, then the cell in the discrete grid and the corresponding node in the whole network will be regarded as permanently occupied in the subsequent evacuation simulation.On the constructed whole network of the office building, the building earthquake evacuation processes to the sample earthquake excitation mentioned above can be simulated using the original PN-FFCA and refined PN-FFCA methods.In this study, it is assumed that the evacuation starts after the pedestrians have received the earthquake warning.In addition, the effects of ground motion on pedestrian movement actions are not considered in this study.Considering the earthquake warning system, it is assumed that pedestrians will begin leaving 10 s before the event.In this evacuation simulation, the main model parameters of the refined PN-FFCA method are set as J S = 2:0, J d = 1:0, b = 10, a = 1:0, and b = 3:0, which have been verified in Example 1. Figure 11 shows the evacuation states at some specific time steps provided by the original PN-FFCA method, while Figure 12 shows those provided by the refined PN-FFCA method.It can be seen that the primary trend of pedestrian evacuation is the same for both models.When evacuation begins, staff on the first and second floors move toward the exit and stairwell, respectively, due to the attraction of the static floor.Subsequently, staff on the second floor enter the first floor through the stairwell and move toward the exit.When all staff reach the exits, the evacuation is completed.It is worth noting that herding behavior among staff is observed throughout the evacuation process due to the dynamic floor field, especially on the second floor.
To show the difference between the original and refined PN-FFCA methods, Figure 13 shows the pedestrian trajectories on the first floor to the considered sample earthquake excitation that are provided, respectively, by the original PN-FFCA and refined PN- FFCA methods.In this figure, different color lines are used to distinguish the evacuation trajectories of the staff coming from the different floors.The red line represents the evacuation trajectory of staff coming from the second floor, and the blue and black lines represent the evacuation trajectory of staff on the first floor.It can be seen in Figure 13 that there is a significant difference between the simulated results using the original and refined PN-FFCA methods.For the original PN-FFCA method, because staff cannot choose exits independently based on the evacuation situation, all staff from the second floor choose Exit 1 to evacuate based on the static floor field, which increases the traffic pressure at Exit 1 and is not conducive to the earthquake evacuation.For the refined PN-FFCA method, when staff reach the first floor through the stairs, they initially select Exit 1 according to the concept of position attraction.Subsequently, with the increasing number of staff around Exit 1, some chose Exit 2 to evacuate according to the concept of density attraction, which can be effectively reflected by the bifurcation of the red trajectory line in Figure 13b.Moreover, staff backtracking is observed in the refined PN-FFCA method, denoted by green lines.These staff initially chose Exit 1 based on the concept of position attraction, but the degree of congestion near Exit 1 increased during the evacuation, so they turned to Exit 2. Based on the above phenomenon, it can be found that the refined PN-FFCA method alleviates the traffic pressure at Exit 1, and staff can adjust their evacuation strategies based on the evacuation situation during the entire evacuation process, making it more reasonable and intelligent.
To further compare the difference in evacuation efficiency between the two methods, the evacuation curves are calculated based on the simulating results from the original and refined PN-FFCA methods and plotted in Figure 14.It can be seen in Figure 14a that, for the original PN-FFCA method, the calculated number of staff who escaped from Exit 2 is only about 10, though more staff can escape from the exit after the time of 10 s.In fact, the 10 staff are those whose initial locations are near Exit 2. In the simulation given by the original PN-FFCA method, when the staff on the second floor reach the first floor, they will select Exit 1 because the stairwell (i.e. the yellow square in Figure 13) is close to Exit 1 and, consequently, Exit 2 is not fully used.Since considering the choice of exits, the number of staff escaped from Exit 2 simulated by the refined PN-FFCA method increases to 20, as shown in Figure 14b, which is about twice that obtained from the original PN-FFCA method.In addition, the slopes of the evacuation curves shown in Figure 14a and b are different.The slope of the curve related to Exit 1 shown in Figure 14a remains nearly constant after about 8 s, while the slope of the curve related to Exit 1 shown in Figure 14b gradually decreases after 10 s.The reason is that the refined PN-FFCA method can let some staff escape Exit 2 after that moment.Finally, Figure 14 shows that the average total evacuation time calculated from the refined PN-FFCA method is 25 s, less than the 36 s calculated from the original PN-FFCA method.
Time-variant velocities of pedestrians.Different from the original PN-FFCA method, the motion velocities of staff are not constant in the refined PN-FFCA method.Two velocity variation curves for the seven staff are obtained by the refined PN-FFCA method and plotted in Figure 15. Figure 15a shows the curves that only consider the effect of density, that is, vd(i), see Equation 13, and Figure 15b shows the curves that consider both density and structural damage, that is, vd(i)Ã, see Equation 19.For clarity, a representative seismic response sample is selected for the analysis.In this sample, the time steps of the occurrence of non-structural component damage on the first floor, non-structural component damage on the second floor, and structural component damage on the second floor are 105, 97, and 106, respectively, which are marked with black, red, and green dashed lines in Figure 15.Note that no structural component damage occurred on the first floor.From Figure 15a, it can be seen that the fluctuation of vd(i) values firstly increases and then decreases with the increase in time steps.When the time step is 45, the value of vd(i) for staff 58 reaches the minimum value of 0.97 m/s.The abovementioned variation is similar to that observed in Example 1 classroom earthquake evacuation simulation.The reason is that although the building types in the two examples are different, the initial distribution of pedestrians is relatively dispersed, and the whole evacuation process goes through scattering, gathering, and then scattering.From Figure 15b, it can be seen that the fluctuation of vd(i)Ã is significant all the time.When structural damage occurs (i.e.time steps = 97), the vd(i)Ã reaches the minimum value of 0.75 m/s.By comparing Figure 15a and b, it can be found that before structural damage occurs, the value of vd(i) is the same as the value of vd(i)Ã, while after structural damage occurs, the value of vd(i)Ã is significantly smaller than the value of vd(i).This indicates that the collapsed structural and non-structural components are the primary hindrance to pedestrian evacuation following structural damage.In addition, when structural damage occurs, the minimum value of vd(i) is approximately twice the minimum value of vd(i)Ã.Therefore, it can be inferred that if the maximum pedestrian density occurs before the occurrence of structural damage, it is extremely beneficial for pedestrian evacuation.It should be pointed out that the refined PN-FFCA method still has some limitations.First, the model parameters identified and used in this study are only suitable to the building earthquake evacuation type where the evacuees have been trained and are all able-bodied.More validation examples are needed to identify the model parameters for different evacuation types when more recorded data are available.Second, only the inter-story drift ratio is considered in this study for determining seismic damage to partition walls and suspended ceilings.In fact, the non-structural damage to the building is partially controlled by the inter-story drift ratio, and more importantly, it is controlled by floor acceleration.Therefore, floor acceleration also needs to be considered in the structural seismic damage analysis.Finally, the method does not consider the couplings between structural and nonstructural damage.In addition, simplified discrete grid sampling is used to determine the damaged non-structural and structural element locations.This may lead to errors to some degree in calculating the time-varying motion velocities of evacuees.Overall, more systematic and in-depth theoretical and experimental research is required to further improve the effectiveness of the refined PN-FFCA method.

Declaration of conflicting interests
The author(s) declared no potential conflicts of interest with respect to the research, authorship, and/ or publication of this article.

Figure 1 .
Figure 1.Connection of the whole network between two adjacent floors.

Figure 2 .
Figure 2. Selection of exit for a specific particle.

Figure 3 .
Figure 3. Video information summary.(a) The classroom layout and evacuation time of each student.(b) Students' evacuation trajectories.

Figure 4 .Figure 5 .
Figure 4.The whole network of the evacuation area of the classroom.

Figure 8 .
Figure 8. Plan layout of the office building.(a) The first floor.(b) The second floor.

Figure 9 .
Figure 9. Floor networks and staff distribution.(a) The first floor.(b) The second floor.

Figure 10 .
Figure 10.A realization of the random distributions of the damaged elements on the floors.(a) The first floor.(b) The second floor.(For interpretation of the references to color in this figure legend, the reader is referred to the web version of this article.)

Figure 11 .
Figure 11.A snapshot of evacuation in the original PN-FFCA method.(a) The first floor.(b) The second floor.

Figure 12 .
Figure 12.A snapshot of evacuation in the refined PN-FFCA method.(a) The first floor.(b) The second floor.

Figure 13 .
Figure 13.Comparison of pedestrian trajectories.(a) The original PN-FFCA method.(b) The refined PN-FFCA method.(For interpretation of the references to color in this figure, the reader is referred to the web version of this article.)

Table 1 .
The typical parameters in Equation1