Toward multi-granularity spatiotemporal simulation modeling of crowd movement for dynamic assessment of tourist carrying capacity

ABSTRACT Dynamic process simulation and prediction of crowd movement are effective approaches to understanding the complex human behavior system in GIScience. At present, obtaining full-sample individual trajectory data still faces challenges because of privacy and cost constraints, thereby resulting in difficulty solving geographic modeling problems that require full-sample individual data. In this paper, a general model for crowd movement simulation is proposed by taking the dynamic evaluation of tourist carrying capacity as an example. Such method is a multi-granularity coupling model, which considers behavioral process and spatiotemporal heterogeneity of tourists. First, a secrete event-based logic model of tourist behavior is proposed. Second, a social force-based inference method of tourist path is designed. Finally, the simulation and evaluation model of remaining spatial carrying capacity of tourists based on a behavioral dynamic system is achieved. In addition, the correctness and applicability of the model are demonstrated through a case study. The proposed model will positively affect time- and space-sharing analysis and assessment of crowd flow within a specific area of activity.


Introduction
The simulation of crowd flow is an important means to reproduce the real scenario of crowd movement. Compared with the crowd flow modeling based on real sampling data, simulation modeling can reflect the reality more comprehensively, and the data of whole sample can be obtained. This method can greatly reduce the cost of data acquisition and solve some GIS modeling problems that failed to be solved only by sample data, such as the evaluation of tourism spatiotemporal carrying capacity. Tourist carrying capacity of a tourist area refers to the maximum number of tourists that can be received in a particular region, and it is an important indicator in the assessment of tourism carrying capacity. This factor affects the development of tourist areas and experience. The tourist carrying capacity plays an increasingly important role in the development and management of tourist areas with the increase of population and the number of transnational tourists. Many scholars have used this indicator in the fields of tourism management and tourism geography (Canestrelli and Costa 1991;Cupul-Magaña and Rodríguez-Troncoso 2017;Naranjo-Arriola 2021).
Problems caused by the crowded flow of people in scenic spots, such as the decline of tourist experience and destruction of scenic resources and ecological environment, have become increasingly serious and exerted a strong negative influence on the reputation of scenic spots in recent years (Popp 2012;Szuster and Peng 2021). Although a large number of scholars have attempted to calculate the overall tourist carrying capacity of tourist areas from various angles (Saveriades 2000;González-Guerrero et al. 2016), the nonuniformity of people flow leads to important changes over time and important heterogeneity in space in real-time carrying capacity of such areas (Y. Chen, Chen, and Mu 2021). This dynamic change process is difficult to describe using the static overall statistics of tourist carrying capacity. Therefore, further investigation on the tourist carrying capacity considering spatiotemporal heterogeneity is necessary for the dynamic prediction and evaluation of tourist flow capacity.
The modeling and analysis of tourist carrying capacity are hindered by two main factors: (1) difficulty in obtaining tracking data of individual tourists of the whole sample and (2) typical exploration of tourist carrying capacity from a multidisciplinary perspective (Navarro Jurado et al. 2012). On the one hand, practical factors, such as privacy protection and data acquisition cost, remarkably reduce the possibility of obtaining full-sample individual behavior data, although data acquisition is technically and objectively feasible. Hence, challenges also exist in further prediction and analysis of dynamic tourist carrying capacity. A large number of studies based on smalland medium-size individual trajectory samples of tourists rather than full samples also present this difficulty (Zhang et al. 2017;Muler Gonzalez, Coromina, and Galí 2018;Huang et al. 2020). On the other hand, the modeling of tourist carrying capacity emphasizes the spatial and temporal heterogeneity of carrying capacity from the geo-spatiotemporal perspective and regards the carrying capacity as a dynamic spatiotemporal process. The modeling method based on traditional GIS analysis is limited in the evaluation of tourist carrying capacity, and dynamically modeling tens of thousands of individuals is challenging (Ren-jun 2005). However, the construction of effective prediction and evaluation methods can only be achieved with the support of individual behavior data of nearly full samples and multigranularity modeling methods. Therefore, the individual tourist behavior trajectory data generated by the simulation and agent-based modeling approach with multiple granularities can be used to achieve this goal.

Research perspectives of scenic carrying capacity
The earliest proposal of the tourism environment capacity related to the carrying capacity of tourism can be traced back to the 1960s (Coccossis and Mexa 2004). Although relevant concepts have been proposed, an in-depth investigation is still lacking. Notably, the number of cases and empirical studies on tourism carrying capacity has grown after a period of development (O'Reilly 1986;Martin and Uysal 1990). Quantitative research on tourism carrying capacity has attracted considerable attention since the beginning of the 21st century. Researchers have begun to build models from different perspectives to evaluate the tourism carrying capacity of various regions from the economic perspective while considering the impact of tourism development and management on economic development as well as calculate a suitable value for economic development (Sowman 1987). Studies on tourism carrying capacity from the social perspective typically consider the impact of social factors, such as the number of tourists and tourism policies; determine the tourism carrying capacity without affecting social development; and formulate a reasonable tourism area management plan (Saveriades 2000;C.-L. Chen and Teng 2016).
Research on tourism carrying capacity from the ecological perspective has attracted considerable attention, given that environmental and resource issues have intensified (Ríos-Jara et al. 2013;Lobo 2015). These studies generally consider the impact of tourism development on the ecological environment and evaluate the tourism carrying capacity of ecological sustainable development. Multi-angle research on tourism carrying capacity has obtained increasingly comprehensive results in recent years with the improvement of tourism carrying capacity theory and mathematical models (Guo and Chung 2019;Wang et al. 2020). The outbreak of COVID-19 in 2020 has exerted a serious negative effect on the tourism industry. Therefore, research on tourism carrying capacity considering the epidemic background has been performed from new perspectives (Bustos et al. 2021), such as emotions (Tokarchuk, Barr, and Cozzio 2022). Hence, the research angle of tourism carrying capacity has become increasingly comprehensive; mathematical models have improved, and studies in conjunction with other disciplines have begun to emerge (H. Zhang et al. 2022). However, these studies are commonly performed at the macrolevel, and only few studies are conducted at the microlevel.
In general, previous studies on tourist carrying capacity in scenic spots relied on tourist behavior sample data, which were obtained through real-time monitoring, mobile phone signaling, and social networks. However, Lin (2011) proposed the use of the Internet of Things technology to build an intelligent management system for tourist areas as well as monitor and evaluate the carrying capacity of tourists in real time. Researchers can analyze behavioral patterns and preferences of tourists because the use of major social media platforms by tourists and the development of technology can provide massive data that can indirectly reflect their tourism behavior (Shao, Chang, and Morrison 2017). Meanwhile, numerous call detail record and mobile phone data are also used to monitor the distribution and flow of tourists in scenic spots, and they can be used in the operation and management of scenic spots (Y. Luo and Zhang 2018;Qin et al. 2019). Some prediction models, such as random forest, have also been applied to predict movement patterns of pedestrians, but the high prediction requirements on sample data and large prediction granularity highly limit the generalizability of methods (Wu and Chen 2021). Methods can only obtain effective evaluation and prediction results in specific scenic spots or within a specified period of time because of the high cost of data acquisition, maintenance, and update as well as personal privacy issues. The requirements for the predictability of tourism management and formulation of specific measures in different periods and spaces have increased because of the rapid development of the tourism industry, increasing the scale of tourists in various places, and constant change of behavioral patterns of tourists. Therefore, a model for the simulation, evaluation, and prediction of tourist carrying capacity in scenic spots with strong generality, multigranularity, and low dependence on tourist sample data is urgently needed.

Analysis methods of scenic carrying capacity
Models of tourism carrying capacity have been widely used in the planning and management of tourist areas in recent years. Some scholars used the time-varying parameter model and empirical mode decomposition to examine natural scenic spots, and back-propagation artificial neural network was utilized to predict the daily tourism carrying capacity (Y. Wang, Xie, and Jie 2019). Other scholars also used theme parks as an example to present the conceptual framework of tourism carrying capacity; analyze the impact of different factors on the tourism carrying capacity using the path model, variance analysis, linear regression test, and other methods; and suggest management strategies (Y. Zhang et al. 2017). Some researchers used an archipelago in Indonesia as an example to calculate the optimal tourism carrying capacity by establishing a coupling model of social and ecological carrying capacity, and the results were applied to the sustainable development of island tourism management and the formulation of island tourism policies (Adrianto et al. 2021). A recent study analyzed the risk of over-tourism in Tunisia from two aspects of tourism intensity and density indexes using the tourism carrying capacity model by calculating and predicting the developmental stage of tourism in the country (Widz and Brzezińska-Wójcik 2020). Researchers developed a tourism carrying capacity index for a Greek island, which focuses on environmental and human factors, to assess island sustainability in the Mediterranean region (Leka et al. 2022).
At present, the majority of tourism recommendations ignore the individual needs of tourists. Thus, some scholars use a multi-agent system to help tourists meet different tourism needs and present personalized suggestions and tourism planning for various needs and preferences of tourists (Chiu and Leung 2005;Sebastia, Giret, and Garcia 2010). Routes can be provided to tourists using a multi-agent system in accordance with their preferences and locations of scenic spots (Lee, Chang, and Wang 2009;Stoyanova-Doycheva et al. 2020). Some scholars proposed a dynamic visit coordination model based on multiple agents for large-scale crowd activities and designed a tour coordination algorithm that comprehensively considers tourist preference, scenic spot congestion, and road impedance (PCI) to improve congestion and social services (Li 2014).
System dynamics methods have been applied to tourism behavior. Some scholars adopted the method of system dynamics to simulate and evaluate the impact of different policies on tourist behavior and provided suggestions for the recovery of tourism after COVID-19 in Maldives (Gu et al. 2022). A quantitative decision support system based on the system dynamics model was proposed to examine the tourism management of the Galapagos Islands and illustrate the impact of tourism and population growth (Pizzitutti et al. 2017). An effective model with system dynamics can be established to analyze the dynamics of tourism activities in a specific region, and the tourism demand of the region can be predicted by simulating different scenarios (Vetitnev, Kopyirin, and Kiseleva 2016;Žibert et al. 2020).
The discrete event model can be used to simulate the diffusion phenomenon of human flow. Some scholars examined the operation of air taxis through discrete event simulation to determine the number of taxis required to meet the urban air traffic demand of New York City (Rajendran and Shulman 2020). Other scholars utilized the discrete event model based on airport saturation caused by the rapid development of tourism and commerce in Mexico in recent years to simulate and evaluate the saturation degree of the country's airspace (Ann, E, and A 2021). The discrete event model can be applied to find ways to reduce the impact of port road traffic on local communities and the environment (Preston et al. 2020). Stewart et al. (2017) explored the impact of tourist location on the evacuation exit. The researchers designed the exit and evacuation model on a unity platform that can provide risk assessment and an emergency evacuation plan for future work by observing the typical time of three tourist groups passing through the tunnel.

General framework of the model
A tourist area contains the two basic elements of tourists and scenic spots. Individual tourists can travel to different scenic spots and wander between them. The tourist area in this scenario contains several scenic spots and individual tourists. These multigranularity geographic elements with nested relationships constitute a typical medium/microgeographic scenario. Therefore, the following problems exist during modeling of the tourist carrying capacity of tourist areas: (1) traditional methods in the macroscopic scale for evaluating the carrying capacity of tourist areas often ignore the heterogeneity of the spatial structure within the region, and (2) the simulation model that portrays the behavior of tourists in the microscopic scale, such as the multi-agent simulation model (Shirzadi Babakan and Alimohammadi 2016; Heppenstall et al. 2021), exponentially increases the computer arithmetic requirements with the increase of the regional scale, thereby leading to difficulties in simulating large tourist areas and high tourist traffic. Therefore, this section will introduce a model framework that considers macrosystem statistics and microbehaviorial analysis for simulating the evaluation of tourist carrying capacity in tourist areas ( Figure 1).
This study aims to simulate individual tourists within and between scenic spots through multiagents at the minimum granularity, treating individual tourists as independent agents to simulate their spatial behavior, and focus on their generated spatial behavioral trajectories. It also models tourist stock and spatial flows among scenic spots using spatial system dynamics at a high granularity and represents the behavioral logic of tourists using discrete events. First, the scenic space and its main tourism resource area distribution are identified, and the Figure 1. Framework of the spatiotemporal carrying capacity evaluation model for a scenic region. The geographical elements in the scenic region (a) are partitioned into sub-regions as shown in (d). Then, a multi-granularity simulation model of crowd movement is constructed by coupling the discrete event model, social force model, and system dynamics model. cognitive spatial scope, perceptual spatial scope of tourists in mesoscopic and microscopic views, and operational spatial scope in the microscopic view are sequentially set as basic construction parameters to describe the behavioral flow of tourists within the tourist area ( Figure 1a). Notably, spatial perception and sensing from the meso-and microperspectives mean the same, except in the setting of the specific range size based on the needs of the model. The former is primarily used for microscopic individual behavior simulation and generally expressed using specific values, whereas the latter is utilized for describing the behavioral flow of the tourist area as a whole and generally defined in accordance with the topological relationship of the behavioral dynamic system. Second, the geographic environment is partitioned and abstracted, and the scenic scene is abstracted into a network structure that considers spatial relationships ( Figure 1b). Finally, a system dynamics model can be established by combining the discrete event model of tourists constructed in accordance with the spatial cognition and perception of tourists ( Figure 1c) and the simulation of tourist flow characteristics based on the social force modeling of the micro-partitioning scenario ( Figure 1d) to describe the tourist flow inside the tourist area in the full spatiotemporal scenario. Thus, a complete evolution of tourist carrying capacity within the tourist area is established (Figure 1e).

Discrete event-based logical model of individual tourism behavior
Although behavioral patterns of individual tourist in tourist areas have been extensively investigated, research on tourist carrying capacity inside tourist areas based on micro-individuals is limited. Notably, obtaining individual track data close to the full sample is difficult, and the calculation of tourist carrying capacity of tourist areas focusing on spatiotemporal heterogeneity within tourist areas must rely on large or even close to the full sample of individual track data. Therefore, a logical model of individual tourist behavior based on discrete events is introduced in this section. The proposed model will lay the foundation for the subsequent simulation of the overall tourist flow and assessment of the tourist carrying capacity of the tourist area.
Tourists' behavior is driven by various motives, such as pursuit of novelty and socialization, which are primarily manifested in six categories: visiting, eating, entertainment, shopping, moving, and staying (Kock, Josiassen, and Assaf 2018). We can simply decompose and abstract the behavioral logic of tourists into multiple discrete tourism behaviors, given that the ability of individuals to simultaneously undertake more than one task is limited ( Figure 2a) (Hägerstrand 1975). Discrete events of tourists' behaviors determine the number and order of their visits to specific scenic spots as well as generate certain behavioral expectations, such as stay and departure times. Behavioral discrete events (behavioral plan) may change as the space and time in which an individual is located change. Therefore, behavioral discrete events of tourists should also present certain end conditions, such as reaching target areas and a specified time.
The end of a discrete event also indicates that a new discrete event is created or that the individual has ended all tourism behavior in the area. The Heppenstall graph provides an effective theoretical framework for expressing behavioral processes of individual tourists under the spatiotemporal background ( Figure 2b) (Hägerstrand 1989). Combining the spatiotemporal path and discrete-event conceptual models can express the typical tourist behavioral pattern (Crease and Reichenbacher 2013). The time period of the experiment is set from 9:00 a.m. to 1:00 p.m. under the assumption that scenic spot S will close at 1:00 p.m. Tourists A and B plan to visit scenic spot S, while tourists C and D plan to visit the lake together. Tourist A enters first through entrance B and leaves the tour area close to 1 p.m. after completing the plan to move to scenic spot A and join the tour. Tourists C and D enter separately shortly after the entrance of tourist A, and the two individuals show similar behavioral plans to move around the lake and join the tour because of similar motivations (their relationship may be similar to that of companions or relatives). The two tourists reach entrance B at 1 p.m., and their current plan ends while generating a new behavioral plan that will begin from entrance B. Tourist B enters the tour area at the latest since the time is close to 1 pm and the site is about to close. Tourist B's plan is forced to end early as the individual leaves the tourist area. Table 1 explains the six types of tourist behavior and the sub-behaviors they contain in detail to facilitate the following quantitative representation of various types of basic tourists' behavior and describe the behavior from a functionalist perspective using three indicators, namely, distance of active, speed of active, and range of active (Huang et al. 2016), where the set of differential equations of the social force model is primarily used to calculate the distance of active and the speed of active of the tourist.

Social force-based simulation algorithm for the behavior of tourist groups
Behavioral motivation of tourists inside a tourist area is influenced not only by individual attributes but also by the dynamic spatiotemporal environment. Therefore, a microscopic individual tourist behavior simulation framework based on the logical model of individual behavior is proposed in this section to achieve the behavior simulation and analysis of tourist groups within the region and provide the data basis for the tourist carrying capacity assessment model while enhancing its credibility and interpretability.
The microscopic behavior simulation architecture for tourists is primarily divided into three layers: One of the most common tourists' behavior in scenic spots. Affected by the psychology of exploration (Schaller et al. 2017), tourists generally visit various attractive things.

Stand-by view
Personal territory (Helbing 1991 No change in one's spatial position microscopic spatial cognition, microscopic spatial perception, and tourist behavior state layers ( Figure 3). This architecture is a microrealization of the discrete-event model constructed in Section 2.2, which is used to describe the specific behavior of tourists. We introduce concepts of spatial cognition and perception at the microscopic scale in this section. Spatial cognition refers to tourist's mental representation of the overall structure of the region, distance, orientation, and size of each specific tourism resource at a specific time and place. The individual's behavioral discrete event model constructed in the previous section is expressed in this layer as a spatial cognitive map in the microcosm, which is generated in the individual's psyche based on his spatial cognition and a comprehensive representation of the current tourist area. The spatial cognitive map includes the order of their expected behavioral goals and information on the general direction of travel, psychological distance to the goal, expected arrival time, and its relationship to drive the tourists to produce a specific tourism behavior. Figure 3a shows the spatial cognitive map constructed by the red tourist at time 1 and located at the beginning. Dark blue areas A, B, and C represent regions with a high degree of awareness for the tourist, which are generally referred to as famous scenic spots. Corresponding light blue areas D, E, and F represent regions with a low degree of awareness for the tourist, which generally rely on visual and other senses for perception. Light yellow areas G and H are impassable areas, whereas the remaining gray regions are road areas. The tourist will generate a behavior plan based on the psychological distance from each dark blue area, location relationship among areas, and time expectations. Dotted red arrows point to targets excluded from the plan. The tourist first calculates the degree of awareness with the tourist resource area A within his/her cognitive range as follows: where is the inherent degree of region awareness by tourist α; D αA represents the distance from tourist α to region A, which is normalized; f αA t ð Þ is a monotonic function from 0 to 1 over time, which indicates the degree of decay of tourist α's perception of region A over time, and h αA is equal to 1 when region A is a highly perceived region and 0 if otherwise.
The higher the cognitive degree, the easier it is to be selected as the target region for the next behavior. Based on the order from largest to smallest to determine, after confirming the first target area, the calculated location will be changed to the location of that area. Then, the abovementioned calculation process Figure 3. The process simulation of tourists' behavior. First, the spatial cognitive map of tourists is initialized in the spatial cognitive layer as shown in (a). Then, the cognitive map of tourists is modified in the spatial perception layer (b and c). Finally, the movement state of the tourists at the next moment is calculated in the behavioral state layer (d and e).
is repeated until the plan size reaches the preset behavior plan size. The unselected region will no longer participate in the comparison. Thus, Figure 3a produces a sequential order of 1, 2, and 3 three behavioral targets. Solid red arrows indicate the general direction of behavior, and the thickness of the arrow indicates the expected value of the direction, which can be calculated as follows: where Spd t αA is the degree of perception of tourist α to region A, and S αA is the position in the behavioral plan. Notably, a high expectation value corresponds to a low probability for change. Figure 3 shows that arrows point to the three planned target regions a, b, and c. The size of the target region is calculated using Equation (3). The default completion of this behavioral event occurs when the individual enters the region. The value is generally set to match the actual size of the area.
( 3) where A α is the actual area of region A. Here, we temporarily disregard the irregularity of its performance shape and regard the regular graph as its contour. A high degree of cognition corresponds to a small ratio of its area to the actual area of the region. The second spatial perception layer is used to correct the cognitive map generated by the spatial perception layer (Figure 3b). The tourist will perceive things within a predefined spatial perception range and then calculate the tourist interest value for all regions separately within the range using Equation (4). A large interest value denotes a high likelihood that the region will be added to the current behavioral plan. Figure 3c presents that a new spatial cognitive map is generated because the tourist adds region D to the current behavioral plan.
where I αA refers to the intrinsic interest value of tourist α to region A. The default is equal to 1 when region A is a high-awareness region. We equate the intrinsic interest value to the intrinsic awareness level, in which is a monotonic function that varies from 0 to 1 over time and D αA is the normalized distance of tourist α to region A. Finally, the current spatiotemporal position and behavioral state of the tourist directly changes because of the third behavioral state layer. This layer is based on the social force model shown in Figure 3d, which can be used to calculate the specific forces on the individual tourist at the current time by confirming obstacles, other individuals, and the target direction in the preset operation space, thereby obtaining the speed and position at the next moment  pedestrian movement toward the target and interaction terms between pedestrians and various entities that are quantified as different forces (Helbing 1991;Helbing and Molnár 1995;Helbing, Farkas, and Vicsek 2000). Equations of the social force model are presented as follows: where m α is the mass of tourist α; ω α is the velocity in the current environment; F α t ð Þ is the combined force on tourist α, and � is a random variable that represents the occasional uncertain behavior of the tourist. The first, second, third, and fourth terms in the expression of F α t ð Þ refer to the driving force, combined force between tourists and combined force of forces acting among other tourists, combined force of forces acting between tourists and various obstacles, and combined force of attraction, respectively.
If external interference is absent, then the tourist will move with velocity v 0 α in the desired direction ẽ α t ð Þ. However, a certain deviation exists in the actual velocity ṽ α t ð Þ from the desired velocity v 0 αẽ α t ð Þ À � because of the necessary avoidance behavior. This deviation can be corrected by the relaxation time τ α , which can be described in the following form of acceleration: The tourist wants to reach the destination r 0 α with a path that can usually be abstracted as the edge of a polygon. Assuming that r k α is the next destination to be reached, then the desired direction of the tourist's motion is expressed as follows: where r α t ð Þ is the actual position of the tourist at moment t.
A repulsive force is exerted on tourist β when tourist α is close to tourist β. This effect can be expressed as follows: where A α and B α are constants representing the intensity and scope of interaction between tourist α and other tourists, respectively; r αβ is the sum of radius of the two tourists interacting; d αβ is the distance between the two tourists, and ñ αβ is the unit vector of tourist β pointing to tourist α. Physical forces of mutual contact are generated when bodies of tourists come into contact with one another, and the sum of the radii of interacting tourists is greater than the distance between them. The physical force is composed of two parts: body force and sliding friction. The body force is a compressive force that resists other tourists and avoids the destruction of the body, as shown in the first term on the right side of Equation 10. Sliding friction force comes from a tendency, that is, a tangential force, to avoid tourists passing others at a large speed and a small distance, as shown in the second term on the right side of Equation 10.
where k and K are constants that denote human elasticity and sliding friction coefficients, respectively, denotes the tangential vector, and is the tangential relative velocity. θ x ð Þ is a segmented function that can be expressed as follows: When x > 0, the distance between interacting tourists is less than their sum of radius, and θ(x) is taken as x, in other cases it is taken as 0, which indicates the absence of contact between tourists. Tourists will maintain a certain distance from the boundary of buildings, walls, and streets that can be expressed in the form of force. The repulsive effect caused by the boundary B can be used in a similar manner to the repulsive force among tourists based on the following equation: Tourists are sometimes attracted to other entities, and the attraction induced at location r i can be calculated using the following monotonically increasing potential function W αirαi ; t ð Þ, which is similar to the repulsive force: where r αi ¼r α Àr i is the position vector difference between the tourist and the attraction point. Attraction decreases with time and interest until it becomes equal to 0. The tourist-object attraction is generally expressed using Equation 9, whereas the tourist-tourist attraction is expressed as the product of the constant C αβ and the vector ñ αβ using Equation 14.

Behavioral dynamics system-based analysis of spatiotemporal tourist carrying capacity
We propose the concept of remaining tourist carrying capacity (RTCC), which is the ratio of the residual receivable tourists to the maximum number of tourists in a specific spatial range, to quantify the tourist carrying capacity of the tourist area from the spatial and temporal perspectives. We construct a dynamic system considering the complexity of the spatial structure and spatiotemporal correlation, which is also known as the behavioral dynamical system, to describe the tourist flow inside the tourist area and calculate the RTCC inside the tourist area dynamically in this section ( Figure 5). The proposed system can enhance the time-space continuity of the calculation results and their interpretability sufficiently. As shown in Figure 5a, the cognitive area range of the unit scenic spot under the mesoperspective set is used as an example. First, the tourist area is spatially partitioned into tourism resource areas a1, a2, a3, a4, and a5 and road network areas s1, s2, s3, s4, and s5 based on the type and boundary of the area and other attributes. Tourism resource areas a1 and a2 are merged because of their adjacent locations and similar characteristics (Figure 5b). Then, discrete geographic regions are abstracted into regional nodes to form a network structure based on geographic location relationships (Figure 5c).
Parameters are set in accordance with the behavioral decision information of the tourist population and the discrete event model of individual tourists constructed in Section 2.2 to determine the flow direction of tourists between regional nodes and their basic flow intensity (Figure 5d). The arrow direction represents the overall behavioral direction of the population between two regions, and the thickness represents the basic flow between the intensity of regions. For example, tourism resource region a4 is regarded as a typical shopping service resource, and road network regions s3, s2, and s5 show two-way circulation and low basic intensity. Thus, the structure of the directed graph based on tourists' behavioral decisions is constituted (Figure 5e) and then used to describe the basic flow trend of tourists within the tourist area.

Figure 5.
Integrating multiple methods to simulate crowd flow. From spatial relationship modeling based on networks (b, c) to spatial interaction modeling of networks with spatial discrete events (d, e) and then to complex spatial process modeling of crowd movement integrating multi-agents, spatial discrete events, and system dynamics (f, g).
Finally, further microscopic behavior simulation experiments using the social force model will be conducted for representative areas to consider the changes in tourist flows at the microscopic scale with changes in the spatiotemporal structure of the area and the total number of people (Figure 5f). The experiments will fit the empirical curves used to describe the variation of pedestrian flow intensity under different scenarios and obtain a table of system properties for each area, which will dynamically affect the variation of the crowd flow rate among areas. Thus, the behavioral dynamics system that considers the spatial heterogeneity and microscopic manifestation of tourists' group behavior is constructed (Figure 5g). The dynamics system can be used to calculate the tourist flow condition in the entire time and space, and its time slice is the expression of the RTCC of each region. The behavioral dynamics system shown in Figure 6 consists of four main components: stock element representing the total number of instantaneous tourists within the region, flow element representing the flow of tourists, attribute table of system variables (attribute links) affecting the dynamics of the flow, and system boundary. The magnitude of each stock depends on the stream to which it is linked, and the stock is calculated as follows: where s denotes the stock value, t denotes the time, f in denotes the flow to the stock, and f out denotes the flow out of the stock.
The flow element describes the behavioral flow of tourists within the tourist area, and the size of the flow rate f is related to the geographical semantics it represents, instantaneous flow rate, and size of the connected stock. F can be expressed as follows: where s in denotes the stock size pointed by nonflow; d denotes the geographical distance of the element connected by the flow; v s in ð Þ denotes the average rate of the crowd obtained from the crowd flow function of s in ; t wait denotes the average dwell time of the crowd, and p denotes the probability of selecting this flow element.
The maximum carrying capacity b imax within each tourism resource region is calculated on the basis of the flow of tourists in each region using the following behavioral dynamics system: where A i denotes the actual area of the region, and b is a constant that denotes the maximum tourist carrying capacity per unit in the region. The saturation degree of each region is then calculated to express RTCC B i of the region. A saturation degree less than 1 indicates that the region still presents a certain RTCC, and that more than 1 indicates that the number of tourists carried in the region exceeds the preset maximum tourist carrying capacity. RTCC B i of the region is calculated as follows: Figure 6. Logical model of the crowd flow based on the spatial dynamics system. A basic behavioral dynamic system consists of regional stocks (s1, a1, and a3), system variable tables, flows (F1, F2, F3, F4, and F5), attribute links, and system boundary.
where s i is the current stock of tourists in the region, and b imax is the maximum tourist carrying capacity of the current region.

Study area and parameter description
Xuanwu Lake Park is a typical comprehensive tourist area with various kinds of tourism resources. In general, in fine-grained attraction modeling, scholars tend to model indoor attractions and outdoor attractions as separate geographical units (Khatibi et al. 2020;Y. Zhang et al. 2017). Moreover, indoor attractions can be modeled independently from outdoor attractions because the attractions in Xuanwu Lake Park have relatively clear geographical boundaries and a certain distance between them. However, in some scenic areas with a mixture of indoor and outdoor attractions and dense scenic areas (for example, natural scenic areas), the interaction among attractions will be difficult to generalize using a systematic network, in which case combining attractions and treating them as a special kind of hybrid attractions for modeling and simulation may achieve better results. The experiment abstracts and classifies tourist areas based on tourists' behavior semantics and characteristics in the tourist area, given that the article primarily models the behavior of tourists (Figure 7). The main scenic spots are divided into indoor and outdoor scenic spots (Figure 7a and Figure 7b, respectively). Indoor scenic spots include museums and memorials. Homogeneous tourist behavioral logic with a high degree of similarity in the relative visiting sequence of various visiting areas is usually generated, given that the interior of the building is a typical regular structure with clear visiting areas and consistent visiting sequence. Meanwhile, outdoor scenic spots are primarily composed of some natural or man-made landscapes. The spatial structure of outdoor spaces is often complex, and viewing areas commonly present a fuzzy boundary because such regions are primarily generated around the landscape. Therefore, tourists' visiting behavior in outdoor spaces usually changes around the spatial location of the main landscape and presents stronger behavioral heterogeneity compared with that in indoor spaces. In addition, restaurants, shops, toilets, and other areas with specific service functions are considered as types of service facilities. These areas are simplified into regions where tourists will stay for a specific time in the experiment given their minimal impact on the overall tourist flow (Figure 7c). The areas connecting scenic spots and service facilities are collectively referred to as types of road network. Tourists will follow the transportation area restricted by the road network and go to their current target location (Figure 7d). The remaining areas where tourists cannot pass are the types of obstacle area. In this study, Ring and Sakura Islands in Nanjing Xuanwu Lake Park are used as sample areas for simulation and analysis. Ring Island covers an area of 127,700 m 2 and contains numerous historic sites, scenic spots, and service facilities. Sakura Island covers an area of 65,900 m 2 with a variety of cherry trees, which is connected to Ring Island by the White Bridge. The experimental area comprises 16 service areas, 10 outdoor scenic spots, and seven indoor scenic spots (Figure 8). Notably, the simulations presented below are under general conditions, and they do not consider specific scenarios (for example, heavy rainfall and fire).
The maximum spatial carrying capacity of main scenic spots of Nanjing Xuanwu Lake is calculated on the basis of the carrying capacity index of the relevant tourist area in accordance with the guidelines for measurement of the carrying capacity of tourist areas issued by the China National Tourism Administration in 2014 (LB/T 034-2014). Among them, the maximum spatial carrying capacity index of the types of road network is 1.5 m 2 /person. The maximum space carrying capacity index of the types of service facilities and indoor scenic spots is 1.9 m 2 /person and that of the types of outdoor scenic spots is 2 m 2 /person. The specific result data of scenic spots are listed in Table 2.

Simulation of the individual behavior process of tourists in scenic spots from the microperspective
The flow of tourists inside each indoor and outdoor scenic spot plays a decisive role in the overall flow of tourists in the tourist area, given that the major tourism resource area is located inside the tourist area. Therefore, individual behavioral processes of tourists in major indoor and outdoor scenic spots are simulated at the microscopic scale to dynamically describe the flow distribution of tourist areas based on the simulated logical model constructed in Section 2.3 and obtain the crowd flow characteristics in specific areas, that is, the crowd flow curve.
Lunar Garden is used as an example of an outdoor scenic spot (Figure 9b). First, the discrete event model of tourist behavior in the scenic spot is built. The process begins when the tourist arrives at the scenic spot. Second, tourists will randomly seek and stay in viewing spots located in transportation areas throughout the scenic spot. Finally, tourists will head for the exit, leave the scenic spot after the departure time, and continue with the remainder of their behavioral discrete event (Figure 9a). Figure 9c presents that tourists enter the region one after another in 60s when most of the region is still unmanned during simulation. Tourists gradually fill the region at 200 s. Small-scale congestion occurs in the area at 600 s. Large-scale congestion exists in the scenic spot, particularly at some intersections, at 1800 s because of the rapid increase of the population.
Simulation processes of indoor scenic spots are similar to those of outdoor scenic spots. Guo Pu Memorial Hall is used as an example in Figure 10b. Tourists will begin by entering the building and then visit viewing areas A, B, and C in a relative order. Finally, tourists will head for the exit, leave the building after visiting, and continue with the remainder of their behavioral discrete event ( Figure 10a). As shown in Figure 10c, tourists enter the region one after another at 60s during simulation. Additional tourists have gathered in viewing area A, whereas other viewing areas temporarily idle. In addition to viewing area A, visiting areas B and C also begin to fill with tourists at 200 s. The intersection of the entrance and viewing area A presents a certain degree of congestion at 600 s, and tourists also begin to gather in other viewing areas. Large-scale congestion is observed in the entrance area and the intersection between viewing areas A and B at 1800 s.  Finally, the average velocity of the tourist flow within the region in unit time is calculated during the simulation (Figure 11). The simulation data clearly showed that the average crowd flow velocity within the area initially decreases and then flattens with the increase of the number of people inside, although indoor and outdoor scenic spots are significantly different with regard to spatial structure and individual behavioral pattern. The ExpDec1 model is used to fit the two groups of experimental data; the fitting model is a one-phase exponential decay function with time-constant parameter, and its full equation is as follows: The results are listed in Table 3. The crowd flow curve of the indoor scenic spot Guo Pu Memorial Hall is estimated using the function y = 0.559*exp(−x/208.323) + 0.106, with a reduced Chi-SQR of only 9.17266E-4 and a maximum R-Square and Adjusted R-Square of 0.92955 and 0.92951, respectively. The fitting curve of Guo Pu Memorial Hall is presented in Figure 11a. The crowd flow curve of the outdoor scenic spot Lunar Garden is estimated using the function y = 0.646*exp(−x/614.979) + 0.098, with a Reduced Chi-SQR of only 3.46903E-4 and a maximum R-square and Adjusted R-Square of 0.97632 and 0.97631, respectively. The fitting curve of Lunar Garden is illustrated in Figure 11b. The results demonstrated that the satisfactory fitting effect of the two crowd flow curves can reflect the characteristics of crowd flow inside scenic spots to a certain extent.

Simulation of the flow process of tourist groups in scenic spots from the mesoscopic perspective
Geographical elements are further abstracted in this study after simulating individual behavioral process within the tourism resource region from the microscopic perspective to obtain the crowd flow in major scenic spots. Moreover, the overall behavioral dynamic system of the Xuanwu Lake Park tourist area is constructed from the medium perspective ( Figure 12). The core part of the Xuanwu Lake tourist area is selected for modeling and then abstracted into a network structure with spatial relations based on the regional classification results in this experiment. Types of scenic spots and road network areas are considered high cognitive areas, and the spatial cognition range of tourists is set as two regional stocks from the medium scale in this experiment.
In addition, types of service areas are considered low cognitive areas, and the spatial perception range of tourists is set as the single regional stock from the  medium scale in this experiment. The discrete model built for tourists' individual behavior is shown in Figure 12a. The expected flow probability among regions is set on the basis of the estimation of the relative attractiveness of each scenic spot and the area of the tourist area. The structure of the network is adjusted on the basis of the discrete event model to form the directed graph structure. Finally, the average crowd flow velocity of each regional stock is provided by the simulated crowd flow curve of each region. The behavioral dynamic system of the Xuanwu Lake tourist area is illustrated in Figure 12b. The hourly average inflow velocity of tourists at the entrance on 3 May 2022, can be obtained using a step interpolation of the hourly crowd count monitoring data provided on the official website of the Xuanwu Lake tourist area (Table 4, http://www.xuanwuhu.net/ ). Crowd flow processes of the major part of Xuanwu Lake tourist area are simulated on the basis of these data. The simulation starts from 7:00, and the number of tourists inside the tourist area is assumed to be equal to 0 at 7:00. Figure 13 shows the comparison of the monitoring data and the simulation results after sample interpolation. The range of crowd flow is primarily within the simulation area of this experiment from the beginning of the experiment until 240 min. Hence, the simulation results are very close to the real monitoring data. However, monitoring data are gradually larger than the simulation results after 240 min, and the difference in values between the two gradually stabilize at about two times because the statistical  scope of the tourist area is larger than the experimental area, including the peripheral walking path around the lake and other regions. Trends of the two curves are nearly the same, and the inflection point of the simulation result curve is slightly earlier than that of the monitoring data curve. This finding is consistent with the influence of monitoring data, including some regions excluded from this experimental simulation. The experimental simulation results are basically consistent with the overall monitoring crowd flow changes. Hence, the proposed model in this article is verified to a certain extent. The simulated spatiotemporal slice of tourist flow processes within the tourist area is shown in Figure 14. Figure 14a and Figure 14b present that the filling speed of tourists in scenic spots near the entrance, exit, and Sakura Island is significantly faster than that in other areas. Two peak periods of crowd flowing into the park are observed at 240 min, and the number of tourists in Ring Island begins to increase. The number of tourists in some popular scenic spots with large areas also generally increases to a high level (Figures 14 c-). The number of tourists in each scenic spot gradually decreases with the decrease of crowd flowing into the park at 650 min (Figure 14f). The road network from the entrance through Sakura Island to the exit uses a relatively simple road structure and densely connects scenic spots. Therefore, the number of tourists carried by the road network area is generally larger than that by the road network area in the left through Ring Island. Notably, congestion is evident in the vicinity of some large scenic spots.

Validation of simulation results
The credibility of the simulation results can be verified by analyzing the correlation between the population monitoring data in the scenic spots and the simulation data. Real-time tourist data cannot be obtained within each of the attractions in the tourism area; thus, Baidu's population heat map data was used as a proxy. Baidu heat map data is a big data product from Baidu Inc., which calculates and generates a population aggregation density map based on the location information of mobile phone users when they visit Baidu products (for example, Baidu search and Baidu map). This paper also analyzes the correlation between the tourist aggregation density extracted from the Baidu heat map and the tourist density generated by the simulation using the geographically weighted regression (GWR) method. Two moments of data were selected for the experiment, 7:00 and 18:00 on 21 August, when the scale of tourists was relatively large and the fitting results were more stable (Zhu et al. 2018). Table 5 shows the results of GWR analysis of the simulation results for these two moments with the population heat data at the corresponding moments, where Sigma is the square root of the regularized residual sum of squares, and R-Square(R2) is the same as Adjusted R-Square (R2-Adjusted), which is a measure of fitness. The results shown in Table 5 indicate that the simulated results correlate well with the real visitor data. Figure 15 (a.1) and Figure 15 (b.) show the simulated tourist densities for the two moments, and Figure 15 (a.2) and Figure 15 (b.) show the extracted tourist aggregation densities for the two moments. Most of the simulated results of these two peak population moments are consistent with the Baidu statistics. However, some areas show that the simulated results differ significantly from the real data, such as the entrance and some of the indoor attractions on the west side of the ring continent, where the density is significantly higher than the simulated results, probably because of the general premise of the simulation, which does not consider the fact that the entrance/exit area generates tourist aggregation (Y. Zhang et al. 2017) and that the hot weather (Brandenburg, Matzarakis, and Arnberger 2007) has a significant impact on tourist behavioral preferences.

Analysis result
The spatiotemporal slice of changes in the RTCC of regional tourists in each region of the tourist area is shown in Figure 16. The size of the RTCC is represented by regional saturation, that is, the number of tourists in the region is divided by the maximum spatial carrying capacity of the region. High saturation corresponds to a small RTCC. The RTCC in major scenic spots (Sakura Island and near the entrance of Ring Island, including Lunar Garden, Wind Lotus Garden, and Lotus Square) gradually reduces to 50% of the maximum spatial carrying capacity and that of the remaining scenic spots idles at the stage where people gradually fill the scenic spots or at the beginning of the process simulation (Figures 16a and 16b, respectively).
The spatial carrying capacity of attractive scenic spots, including Lotus Square and Lawn, nearly reaches the maximum level for a long time, and their RTCC continues to decrease after the peak period of entry. Meanwhile, the remaining part of scenic spots, such as Rockery Waterfall, is in the idle period because its distance from the entrance and exit is large and the scenic spot is small (Figures 16 c-). The level of RTCC in the majority of scenic spots gradually returns to the initial stage after the decrease of crowd flowing into the park (Figure 16b and Figure 16f). The RTCC of processes in the majority of road network areas is at a high level, although some paths connecting scenic spots are saturated, and only a small part is lower than 50%. Notably, roads in the road network structure of the western part of Ring Island typically idle because of their complexity. This finding indicated that the construction of road network near some popular scenic spots should be strengthened, such as increasing parallel roads and widening roads.
The variation of the RTCC in each scenic spot is important in evaluating the tourist carrying capacity of the whole tourist area. The changing curves of outdoor scenic spots are shown in Figure 17. The majority of scenic spots can maintain their RTCC at more than half of the maximum spatial carrying capacity. However, the number of tourists carried by the Mi Fu Stone, Nona Tower, and Lama Temple gradually exceeds their maximum spatial carrying capacity at 100 min, and Nona Tower and Lama Temple even reach 3.5 times at one time.
On the one hand, Mi Fu Stone is located at the exit, and Nona Tower and Lama Temple are located at the connecting part of Ring and Sakura Islands. In addition, the cultural attributes of both tourist spots attract a large number of tourists. On the other hand, their internal microscopic space is relatively narrow. Their few channels, single structure, and concentrated viewing points result in crowd congestion, which is detrimental to crowd flow, particularly in Lama Temple. By contrast, Lunar Garden and Tower of Friendship are located at the entrance with a large flow of people. However, its looped microscopic spatial structure and scattered viewing points allow for satisfactory internal mobility and maintain its RTCC at a high level. Some areas, such as Huanliu Pavilion, are unattractive, and they exhibit low instantaneous flow of people. The high maximum carrying capacity of other areas, such as Sakura King, allows the region to maintain large RTCC at about 70% of the maximum carrying capacity.
The changing curves of indoor scenic spots are shown in Figure 18. Compared with those of outdoor scenic spots, changing curves of indoor scenic spots are more similar because of their matching microscopic internal structure and internal behavioral logic. The RTCC of each region is closely related to its geographical location and attraction. Attractions of Lianhu Academy, Guo Pu Memorial Hall, and Gunley Art Appreciation Center are similar. However, the RTCC of Lianhu Academy is relatively lower than the other two scenic spots while maintaining consistency in the changing trend because of its shorter distance from the entrance and larger number of tourists. The three other scenic spots are located in the middle part of the whole region, and their change is flatter than the first three regions because of their higher maximum carrying capacity. The RTCC is high, and it can be maintained at approximately 80% of the maximum carrying capacity.
The RTCC of the tourist area can be maintained at more than 50% of the maximum carrying capacity, and roads typically idle after two evident peak periods of entry. Hence, the whole tour process of tourists is comfortable. Indoor scenic spots commonly demonstrate better recovery capability, and they can rapidly recover the RTCC, although a more significant peak exists during the peak period of crowd flow. By contrast, outdoor scenic spots demonstrate relatively flat changing curves. In addition, scenic spots exhibit intensive differences. Some regions, such as Lama Temple, show serious congestion, whereas other regions, such as Lianhu Academy, indicate a congestion period during the peak period. Therefore, certain guidance or diversion measures remain necessary for specific periods in some regions.

Potential applications of spatiotemporal simulation modeling of crowd movement
The proposed analytical framework has two main practical implications for tourist area management. First, the dynamic simulation and quantitative measurement of the tourist carrying capacity can intuitively allow managers to understand the changing trend of such capacity within a tourism area and its corresponding spatial variability and relevance. This understanding can provide experimental support for the initiative of the tourism area to control the flow of tourists and prevent the deterioration of tourist experiences or the waste of several tourism resources because of the crowding of a certain part of the area during a specific time. Second, tourism area managers can formulate relevant scheduling and control rules in advance and use the proposed methodological framework to carry out multivariable, multi-scenario simulations to fully verify their effectiveness. This simulation-based approach to scenic area management can be particularly useful in scenarios such as emergency tourist evacuation and dynamic tourism resource deployment, which are difficult to reproduce realistically because of cost and other constraints.

Shortcoming and future directions
In-depth knowledge mining and interpretation are necessary for the experimental results, given that this study focuses on the introduction of the proposed computational model. Moreover, the effective conversion of the simulation and analysis results into an actual plan for scenic spot management still needs further exploration from the practical application perspective. Furthermore, simulations for different special scenarios must be considered by introducing other scenario-specific tourist behaviors, such as visitor panic behavior in the event of a disaster (Helbing, Farkas, and Vicsek 2000) and following behavior during orientation (Morone and Samanidou 2008). For some sudden and dangerous events, an advanced spatiotemporal simulation of multiple scenarios can be an important and effective method of early warning and response.

Conclusion
This study proposes a multi-granularity simulation and evaluation model for the spatiotemporal change process of tourist carrying capacity in tourist areas. It also suggests the conception of RTCC to indicate the number of additional tourists that can be accommodated within a certain space in a scenic spot for a period. After the dynamic simulation and evaluation on 3 May 2022, the case study reveals that the RTCC within the Nanjing Xuanwu Lake tourist area can be maintained at a comfortable level for a long time from the global perspective. However, the western area of Ring Island can significantly waste tourism resources when management guidance measures are not implemented, as this part is clearly more vacant than the eastern area and Sakura Island from the local perspective. Moreover, some scenic spots still exhibit insufficient RTCC during certain periods, such as the peak period of a park. Therefore, scenic spots with insufficient RTCC and their adjacent areas should implement advanced targeted diversion and flow restriction measures during a specific period to avoid large-scale cross-regional congestion.
The experimental results indicate that the heterogeneity of the internal spatial structure in a scenic spot, continuity of time sequence, and difference in tourist behavior exerts a serious influence on the calculation and evaluation of RTCC. The proposed method not only provides a multi-granularity dynamic behavior simulation approach for the quantitative analysis of tourism geography on crowd movement but also plays an important role in resource allocation, emergency warning, and dynamic control in the management of crowded activity areas (e.g. scenic areas). This study can promote quantitative approaches in the field of tourism geography from analytical modeling to simulation modeling. The simulated data can also solve problems about quality and acquisition costs in conventional monitoring data.

Disclosure statement
No potential conflict of interest was reported by the author(s).