The Problem of Effective Evacuation of the Population from Floodplains under Threat of Flooding: Algorithmic and Software Support with Shortage of Resources

Abstract

Extreme flooding of the floodplains of large lowland rivers poses a danger to the population due to the vastness of the flooded areas. This requires the organization of safe evacuation in conditions of a shortage of temporary and transport resources due to significant differences in the moments of flooding of different spatial parts. We consider the case of a shortage of evacuation vehicles, in which the safe evacuation of the entire population to permanent evacuation points is impossible. Therefore, the evacuation is divided into two stages with the organization of temporary evacuation points on evacuation routes. Our goal is to develop a method for analyzing the minimum resource requirement for the safe evacuation of the population of floodplain territories based on a mathematical model of flood dynamics and minimizing the number of vehicles on a set of safe evacuation schedules. The core of the approach is a numerical hydrodynamic model in shallow water approximation. Modeling the hydrological regime of a real water body requires a multi-layer geoinformation model of the territory with layers of relief, channel structure, and social infrastructure. High-performance computing is performed on GPUs using CUDA. The optimization problem is a variant of the resource investment problem of scheduling theory with deadlines for completing work and is solved on the basis of a heuristic algorithm. We use the results of numerical simulation of floods for the Northern part of the Volga-Akhtuba floodplain to plot the dependence of the minimum number of vehicles that ensure the safe evacuation of the population. The minimum transport resources depend on the water discharge in the Volga river, the start of the evacuation, and the localization of temporary evacuation points. The developed algorithm constructs a set of safe evacuation schedules for the minimum allowable number of vehicles in various flood scenarios. The population evacuation schedules constructed for the Volga-Akhtuba floodplain can be used in practice for various vast river valleys.

1. Introduction

Statistics show an increase in the number of floods and annual losses from them, models forecast an increase in these trends [1,2,3], and 35–45 percent of the total number of natural disasters are due to hydrological events. The Center for Research on the Epidemiology of Disasters (CRED) 2022 reports highlight floods as the greatest hazard for many regions [4]. The death toll for the year was 5193 in Pakistan, India, Nigeria, South Africa, and Brazil alone. Moreover, the number of people affected by the floods in Pakistan and Bangladesh alone exceeded 40 million people. The economic losses are $35 billion in total in Pakistan, Australia, China, Nigeria, and India. For example, global climate and hydrological forecasts indicate an increase in the area that will be affected by floods in the future [5]. Thus, urbanized regions, transport networks, and basic infrastructure are under increasing threat of destruction due to floods.

The situation with the frequency of floods and their consequences in Russia and the CIS as a whole is in line with global trends [6,7,8,9]. Severe floods can be distinguished in the basin of the Amur River [10], the drainage basin of Lake Baikal [11], the regions of the Sochi Black Sea Coast and Krasnodar Territory [12] and other regions. Catastrophic events with a flood depth of more than 1.5 m are among the most severe types of emergency situations. Cases of flooding over an area of 100–1000 km${}^2$ belong to the category of large floods, and the excess of 1000 km${}^2$ has catastrophic consequences. Catastrophic floods due to heavy rains occur regularly in many regions of the country, the consequences of which are intensified due to the destruction of dams or the overflow of reservoirs [13,14]. Heavy rains in China triggered floods across the country in July 2020 that affected several million people and killed more than 140. Bangladesh is an example of a huge densely populated mega-region with regular flooding, since most of the country is floodplain [15].

Long-term efforts to create the hydraulic structures complex made it possible to significantly reduce the risk of high water in the Neva River and flooding of St. Petersburg [16]. A separate problem is surge floods, which are typical, for example, for the Kaliningrad region, the Lower Don [17], the Kuban river delta [7], the Black Sea bays [18], etc. The causes of severe floods are due to tsunami [19] and ice jams [20]. The most severe consequences are, as a rule, floods in mountainous areas. Note the catastrophic events in Krymsk (2012, Caucasus Mountains) [21], the accident at the Vajont dam (Vajont, 1963, Italy) [22,23], etc.

The limited time and transport resources are the main problem in organizing the effective evacuation of the population during floods. However, this problem can be largely compensated for by the availability of a relatively accurate forecast of the development of the emergency for some regions. Such water objects are the floodplains of large lowland rivers regulated by hydroelectric dams. Therefore, floodplains pose a serious danger to the life of the population during spring floods or heavy rainfall [15,24,25,26,27]. The effects of flash floods in urban areas are more severe. The influence of the anthropogenic factor can increase the number of flood victims in urban areas, especially in the road network [28,29]. Therefore, evacuation plans must be well-prepared and include the road network and traffic conditions.

Dam hydrograph control systems stimulate intensive floodplain urbanization and reduce the threat of large floods, but do not eliminate them. A slight excess of the safe level of flooding of these territories occurs quite often [1,30].

The low population density in vast floodplain areas and the flat relief make it very difficult to organize effective evacuation by land transport due to flooding of evacuation routes. On the other hand, the flooding of territories occurs relatively slowly with large differences in the flooding moments of different parts. This makes it possible to organize safe land evacuation with a relatively small number of vehicles in the presence of a reliable forecast of the spatial dynamics of flood waters and timely notification of the population. Therefore, the algorithms for organizing the safe evacuation of the population of floodplain territories under the conditions of a fairly accurate forecast of the development of floods are of great practical importance [31,32]. The improvement of high-performance computing technologies provides such forecasting and the development of appropriate algorithms has become possible on the basis of interdisciplinary simulations and solving various high-dimensional optimization problems [17,33,34,35,36,37,38,39,40].

Perhaps the main difficulty in modeling and optimizing evacuation processes is the high uncertainty of the components of the evacuation situation, in particular, the number of the evacuated population in each of the human settlements (HS), the moments of flooding of the HS, sections of the road network, the start time of evacuation, the number of vehicles for evacuation, localization of evacuation points (EP). The criterion for the effectiveness of the organization of evacuation is also ambiguous.

The problems of organizing the evacuation of the population during floods and modeling the corresponding evacuation situations are actively discussed in modern scientific literature. The most popular areas of research are listed below.

(A)

The choice of evacuation routes and locations of evacuation points are traditional problems [41,42,43,44,45,46,47,48,49,50,51,52]. The evacuation time is often chosen as an efficiency criterion [42,44,46,47,48,49], and the evacuation cost is an additional criterion [45,46,49].

(B)

Self-evacuation transport models [31,32,53,54] belong to the class of traffic microsimulations, describing the spread of traffic jams in the transport network. Stochastic programming methods allow minimizing the risk and time spent on the road to develop evacuation recommendations [54].

(C)

Uncertainty in the behavior of self-evacuating people when choosing evacuation routes (ER) and EP is studied on the basis of game-theoretic and agent-based modeling [47,50,53,55,56,57]. Simulating individual behavior over large areas using complex behavior models is computationally intensive.

These problems are relevant when urban areas are flooded, when overcoming the transport collapse with unorganized spontaneous self-evacuation of the population is essential. A high rate of flooding or the lack of a good forecast, the remoteness of evacuation points, the late start of evacuation require the involvement of the maximum possible number of vehicles to minimize evacuation time, which can exacerbate the problem of transport collapse. Therefore, the optimal choice of evacuation routes is an effective means of reducing the risk to life. Modeling the uncertainty in people’s behavior when choosing self-evacuation routes makes it possible to most accurately predict the development of an evacuation situation and make the choice of optimal organized evacuation routes. Maximizing the number of evacuation vehicles is not relevant in the case of low population density over a large area, low flooding rate, and a sufficiently accurate flood forecast. Thus, the traffic flows of spontaneous and organized evacuation do not create life risks. On the other hand, the lack of evacuation transports can be a major problem in such a situation.

Probabilistic modeling of the flooding occurrence and assessment of its risk for different zones of the territory are carried out in [58,59,60,61]. Flood risk assessment in the middle Chao Phraya River Basin (CPRB) was simultaneously analyzed and mapped as the product of flood hazard, and social vulnerability maps generated by fuzzy Analytic Hierarchy Process (AHP) and fuzzy logic [62]. This was calculated based on a physical status of evacuees, safe evacuation condition, shortest evacuation path, flood shelter, and road capacity. The aim of [63] is to identify urban areas in the Hail region that would be affected by simulated worst-case flash floods. The data sources for the soil type, infiltration, and initial moisture were utilized to create the coverage and index maps. To generate virtual floods, the authors ran the GSSHA model within the Watershed Modeling System (WMS) program to create the hazard map for flash flooding.

The authors of [64] performed the flood hazard assessment for a simulated 500-year flood event in the downstream floodplain of Onondaga Creek within Syracuse City, New York State, USA. Communities near Onondaga Creek were assigned designated evacuation shelters based on accessibility and distance to the shelters. The shortest possible evacuation routes were calculated. The impact of the lack of transport and time resources during evacuation on the mortality rate during floods is studied in [65]. The use of temporary engineering structures in areas of possible catastrophic flooding can significantly increase the efficiency of evacuation. The possibility of reducing flood damage through the construction of anti-flood dams is considered in [26,66]. An important direction is to identify the negative impact of floods on the health of the evacuated population [67,68]. It should be emphasized that such studies are carried out for specific areas, since each real water body has characteristic features [26,31,41,42,43,48,50,53,55,56,57,58,65].

General issues in [46,69,70] are the closest to our study. The above article is devoted to the problem of finding optimal evacuation routes using temporary shelters to minimize the total evacuation time and the risk of a possible secondary evacuation from several temporary shelters to permanent ones with restrictions on the total evacuation time, the number of vehicles and the capacity of the shelters. The authors proposed an approximation algorithm of polynomial complexity for solving the formulated two-criteria problem of integer programming. In addition, three model relationships between speed and traffic density are investigated. However, a specific territory is not considered in [46], and an example of model calculations is presented only for three NPs and three EPs.

Finally, since the most important component of evacuation modeling is a hydrodynamic model for forecasting the hydrological regime of a flooded area, we point out the current results of such computer modeling for various water systems [12,16,17,21,23,33,38,71,72]. Modeling of floods near rivers has shown its effectiveness for constructing cadastral maps of water supply [24,36,73]. Numerical models of water dynamics make it possible to better understand the mechanisms of occurrence of catastrophic consequences of floods, for example, for the tragedy in Krymsk (Krasnodar region, Russia) in 2012 [21]. Tsunamis or sea storm surges in estuaries, especially in the presence of wide deltas, are a source of flooding. Hydrodynamic models of such events describe well the dynamics of water level rise [8,17,37,74]. Separately, we single out a series of works devoted to the study of the hydrological regime of the Volga-Akhtuba Floodplain (VAF) during the spring flood, since the VAF is the testing ground for our study [24,25,75,76,77,78].

The minimization of risk to life is the general task of effectively organizing the evacuation of the population during floods in the literature cited above. Special cases are the minimization of the evacuation time of the population, and the maximization of the number of evacuated residents under given conditions (evacuation time and number of vehicles). The result is an optimal traffic schedule. Solving this problem is associated with several difficulties. First, there are no efficient solution algorithms due to the large dimension. Second, the constructed schedules are very sensitive to small changes in the number of vehicles. This instability hinders the practical implementation of optimal schedules under conditions of real transport uncertainty. Third, such optimal schedules may not be safe because they do not cover the entire population.

Our approach is based on solving the adjoint problem of finding the minimum number of vehicles that ensure the safe evacuation of the entire population. The proposed algorithm constructs a sequence of safe schedules corresponding to a sequence of decreasing number of vehicles. It is important that each of the schedules differs slightly from the two adjacent ones. Thus, the set of safe schedules is stable to changes in the number of vehicles and therefore can be used in real floods under conditions of transport uncertainty. The proposed approach is most effective in conditions of moderate floods in river floodplains with low population density and relative proximity to evacuation points. These conditions provide a sufficient margin of time for the repeated movement of each evacuation transport both along the same and along several different evacuation routes at a low average traffic flow density of the self-evacuating population.

The aim of our work is to develop a method for analyzing the minimum resource support for the safe evacuation of the population of floodplain areas based on a mathematical model for the development of floods and solving the problem of minimizing the number of vehicles on a set of safe evacuation schedules. The work [79] is devoted to finding the dependence of the vehicles minimum number on the parameters of the flood and the start of the evacuation, ensuring the safe evacuation of the entire population to permanent evacuation points. This work is aimed at solving the problem of safe evacuation in the event of a shortage of vehicles. In this case, the safety of evacuation is ensured by its division into two stages with the deployment of temporary evacuation points along evacuation routes.

2. Formulation of the Problem

Our study is aimed at calculating the functional dependence of the acceptable risk area boundary on the parameters of the evacuation situation: the hydrograph, the number of evacuation vehicles, the moment the evacuation begins, etc. We consider the evacuation of the entire population with its stay on ERs and in temporary EPs for an acceptable time as an acceptable risk. Here we do not study the problems of unorganized evacuation and the presence of the population in flooded HSs. Thus, our goal is to study the boundary of the domain of solution existence to the conditions system

$$R_h(E,Sh)=0,R_e(E,Sh)\le R^{\left(lim\right)},$$

(1)

where $R_h$ is the function of hydrological risk (damage to health and the threat of life loss due to being in a flooded area), $R_e={\displaystyle \underset{i,j}{max}}\left({\vartheta }_{ij}\right)$ ($i=1,\dots ,N_j$, $j=1,\dots ,n$) is the risk function due to a long stay on the ER and/or in the temporary EP under adverse external conditions (damage to health), $N_j$ is the number of inhabitants of the settlement j, ${\vartheta }_{ij}$ is the duration of stay of the i-th inhabitant of the j-th HS on the ER in the temporary EP, n is the number of human settlements, $R^{\left(lim\right)}$ is the maximum admissible value of the function $R_e$. The evacuation situation tuple E in (1) is defined as follows:

$$E\left(t\right)=\langle K_{fl}(t,Q),\overrightarrow{N},G,{\overrightarrow{t}}^{\left(0\right)},M\rangle ,$$

where t is the time after the start of the flood, $K_{fl}(t,Q)$ is the digital map of the flooded area, $Q=Q\left(t\right)$ is the river discharge (hydrograph), $\overrightarrow{N}=\left\{N_1\dots N_n\right\}$ is the vector of the total number of human settlements, G is the road network graph, ${\overrightarrow{t}}^{\left(0\right)}=\left\{t_1^{\left(0\right)}\dots t_n^{\left(0\right)}\right\}$ is the evacuation start moment vector, M is the number of evacuation vehicles.

The construction of the parameters ${\overrightarrow{t}}^{\left(0\right)}$ and M is the result of the preliminary stage of organizing the evacuation. These values form the time and transport resources of the evacuation schedule $Sh=\langle \tilde{G},\left|\right|m_{i,j}\left|\right|\rangle$, where $\tilde{G}$ is the evacuation route graph, $\parallel m_{i,j}\parallel$ is the evacuation transport schedule matrix ($t=t^{\left(0\right)},\cdots T$; $j=1,\cdots n$) along directed graph edges evacuation routes $\tilde{G}$. The elements of the matrix determine the number of vehicles on the j-th ER at the time t, where T is the end of the evacuation, not exceeding the flooding time of HS and ER according to the maps $K_{fl}(t,Q)$, $t^{\left(0\right)}=min(t_1^{\left(0\right)},\cdots ,t_n^{\left(0\right)})$.

The solution of the optimization problem

$$M(K_{fl}(t,Q),\overrightarrow{N},\tilde{G},{\overrightarrow{t}}^{\left(0\right)},S{h}^{sec},R^{\left(lim\right)})⟶\underset{S{h}^{sec}}{min}$$

(2)

belongs to the boundary of the existence domain of the system solution (1), where $S{h}^{sec}$ is the admissible (safe) schedule that ensures the fulfillment of constraints (1). This problem is a variant of the resource investment problem or a mirror copy of the classic RCPSP (Resource-constrained project scheduling with time windows) scheduling problem with deadlines [80,81,82,83,84].

We propose the algorithm for solving the problem (2) and the results of its implementation in simulation modeling of flood development in the Northern part of the Volga-Akhtuba floodplain (further abbreviations NPVAF and VAF are used).

3. Materials and Methods

3.1. Hydrodynamic Model of Flooding

This section is devoted to describing a model for forecasting the hydrological situation in an area that is under the threat of severe flooding. The model is based on shallow water equations [34,35,36,75] and has already been used to solve a wide range of problems related to modeling the hydrological regime as the Northern part of the Volga-Akhtuba Floodplain (NPVAF) and other river systems. This numerical model became the basis for the creation of flood cadastral maps [24], decision-making systems for the socio-economic development of floodplain territories [25], and the implementation expertise of various engineering projects [76].

The dynamics of surface waters on a non-uniform terrain $b(x,y)$ are based on the system of Saint-Venant equations in the following form [33,37,85,86]:

$$\frac{\mathsf{\partial }H}{\mathsf{\partial }t}+\frac{\mathsf{\partial }\left(H{u}_x\right)}{\mathsf{\partial }x}+\frac{\mathsf{\partial }\left(H{u}_y\right)}{\mathsf{\partial }y}=q(x,y,t),$$

(3)

$$\frac{\mathsf{\partial }\left(H{u}_x\right)}{\mathsf{\partial }t}+u_x\frac{\mathsf{\partial }\left(H{u}_x\right)}{\mathsf{\partial }x}+u_y\frac{\mathsf{\partial }\left(H{u}_x\right)}{\mathsf{\partial }y}=-gH\frac{\mathsf{\partial }}{\mathsf{\partial }x}(b+H)+H{f}_x^{\left(sum\right)}+H{f}_x^{\left(q\right)},$$

(4)

$$\frac{\mathsf{\partial }\left(H{u}_y\right)}{\mathsf{\partial }t}+u_x\frac{\mathsf{\partial }\left(H{u}_y\right)}{\mathsf{\partial }x}+u_y\frac{\mathsf{\partial }\left(H{u}_y\right)}{\mathsf{\partial }y}=-gH\frac{\mathsf{\partial }}{\mathsf{\partial }y}(b+H)+H{f}_y^{\left(sum\right)}+H{f}_y^{\left(q\right)},$$

(5)

where $H(x,y,t)$ is the surface water depth, $\overrightarrow{\mathbf{u}}(x,y,t)=\left\{u_x,u_y\right\}$ is the velocity of the vertically averaged flow determined by the x and y components ($u_x,u_y$), q is the functions sum of sources ($q>0$) and sinks ($q<0$) of water, g is the acceleration of gravity, $(f_x^{\left(q\right)},f_y^{\left(q\right)})$ is the specific force associated with the given source function q (because changes in mass due to q cause corresponding changes in momentum). The specific force ${\overrightarrow{\mathbf{f}}}^{\left(sum\right)}=\left(f_x^{\left(sum\right)},f_y^{\left(sum\right)}\right)$ depends on various factors [33,37,72], including the bottom friction ${\overrightarrow{\mathbf{f}}}^{\left(b\right)}$, the Coriolis forces, the wind action with speed of $\overrightarrow{\mathbf{W}}$, the internal friction due to turbulent viscosity, etc. (Figure 1). Equation (3) is the mass conservation law for an incompressible liquid layer with the source q. Equations (4) and (5) determine the law of change of two momentum components under the action of forces on the right side of these equations.

The quality of a hydrological forecast essentially depends on the force of bottom friction and the properties of the underlying surface

$${\overrightarrow{\mathbf{f}}}^{\left(b\right)}=-\frac{g{n}_M^2}{H^{4/3}}\left|\overrightarrow{\mathbf{u}}\right|\overrightarrow{\mathbf{u}},$$

(6)

where $n_M(x,y;H)$ is the roughness coefficient (or Manning coefficient), which determines the intensity of hydrological resistance to flow [87]. The value of $n_M(x,y)$ can strongly depend on coordinates ranging from 0.02 to 0.2, depending on the features of the underlying earth’s surface. The unit of $n_M$ is s/m${}^{1/3}$ and is traditionally not specified.

The source function is determined by two factors. First, the work of hydraulic structures, precipitation, river, and surface runoff are responsible for the ingress of water into the computational domain. Second, infiltration into the soil and evaporation removes water from the study area (see Figure 1) [85]. River discharge is the volume of water passing through the cross-section per unit time ($Q\left(t\right)$). The dependence $Q\left(t\right)$ is also called the river hydrograph, which is the most important characteristic of the hydrological regime of a water body and floodplain areas (Figure 2). The period of seasonal flooding in spring or due to prolonged rains is characterized by an increase in Q by several times compared to the low water level. Such a hydrograph may have a complex form (see Figure 2), but we use a rectangular hydrograph with a maximum value of $Q^{(max)}$ and a low water value of $Q^{(min)}$ (see the inset in Figure 2) when solving the problem of organizing safe evacuation in Section 3.5.

Numerical integration of Equations (3)–(5) is based on the Combined Smoothed Particle Hydrodynamics—Total Variation Diminishing, CSPH-TVD), which uses both the smoothed particle algorithm (SPH) and the TVD mesh algorithm at various computational stages [35,72].

The method provides conservatism and good balance. It is possible to carry out end-to-end stable calculations of the moving boundary between the dry bottom and the liquid on the topography, which contains arbitrary inhomogeneities, including waterfalls [87]. These properties of the numerical algorithm make it possible to simulate the dynamics of flooding with non-monotonous behavior of $Q\left(t\right)$, when water enters a certain area and then leaves it. The TVD part of the algorithm requires the use of a grid that covers the entire computational domain. We restrict ourselves to a homogeneous grid in the Cartesian coordinate system with cells $(x_i;y_j)$:

$$x_i=x_0+i\mathsf{\Delta }x,y_j=y_0+j\mathsf{\Delta }y,i=0,1,\dots ,N_x,j=0,1,\dots ,N_y,$$

where $\mathsf{\Delta }x$ and $\mathsf{\Delta }y$ are the cell sizes in the corresponding direction. The product $N_x·N_y$ gives the total number of cells in the numerical model.

3.2. Digital Hydrological Landscape Model as the Basis for Flood Modeling

We use the tested model of surface water dynamics for the Northern part of the Volga-Akhtuba Floodplain (NPVAF, Figure 3). The inset in the figure shows the NPVAF between the Volga and Akhtuba rivers below the dam of the Volga HPP (VHPP). The typical distance between the Volga river and the Akhtuba River in the VAF is 25–30 km. The VAF region is covered by a grid with the cell size of $\mathsf{\Delta }x=\mathsf{\Delta }y=15$ m as the base model.

The choice of a specific area for water regime modeling requires setting the Digital Hydrological Landscape Model (DHLM [77]), which includes a set of spatial matrixes on the $(x_i;y_j)$ grid. When studying the dynamics of surface waters, one can limit oneself to the matrixes $b_{ij}$, $n_{Mij}$. More complex models with groundwater and sediment dynamics require additional matrixes characterizing the properties of soil and water-resistant layer [35,75]. The digital hydrological landscape model defines the hydrological structure of the territory, consisting of surface flow, groundwater, and atmospheric water, and establishes the relationships between natural structures and hydrological processes [88]. The DHLM contains a set of digital geoinformation models required for hydrodynamic simulations in a given area.

The DHLM is based on a digital elevation model (DEM) and data on the hydrological network, which is a vector map of all watercourses. Small-scale inhomogeneities of the underlying surface provide hydrological resistance to the flow, which requires specifying in (6) the spatial distribution of the roughness coefficient $n_M$ on the grid $(x_i;y_j)$ (Figure 4). Different colors show different levels of $n_M(x,y)$), which are determined, among other things, by the nature of vegetation, soil properties, building density, etc. For the study of groundwater hydrology, the DHLM contains data on the depth of groundwater water, soil characteristics, infiltration coefficient, and porosity coefficient. In addition, water dynamics may depend on meteorological characteristics, such as the amount of precipitation, snowmelt patterns, and the rate of moisture evaporation.

The hydrological regime of NPVAF is determined by the discharge of water flow through the dam $Q\left(t\right)$ (m${}^3$/s). The territory of the entire VAF is characterized by a complex branched system of natural canals, the so-called eriks, which ensure the water flow from the Akhtuba and the Volga to the flat part of the interfluve when approximately Q > 17,000 m${}^3$/s [77]. Moreover, approximately 70–80% of water enters the floodplain through the Akhtuba River and 20–30% directly from the Volga through small canals on the left side. The contribution from the Akhtuba is the main one and the corresponding flooded area dominates compared to the direct contribution from the Volga (see Figure 3).

The value of water flow from the Volga river to Akhtuba plays a key role in the hydrological regime of the entire floodplain. The water discharge in the Akhtuba River $Q_A\left(t\right)$ directly depends on the discharge in the Volga river $Q_V\left(t\right)$. Figure 5 shows such a hydrological relationship based on the results of numerical simulation for 2022, which has a characteristic hysteresis shape. The ratio $Q_A/Q_V$ is in the range of 0.02–0.08, and only two percent of the Volga’s water enters the Akhtuba during the low season. The spring flooding stage is characterized by an increase in this proportion to 8 percent at its peak.

If the discharge level of the Volga river $Q^{(max)}$ exceeds approximately 30,000 m${}^3$/s, then an emergency occurs in the NPVAF. The catastrophic situation begins approximately with $Q^{(max)}$ > 35,000 m${}^3$/s. The level above 45 thousand m${}^3$/s in the Volga channel has not been observed in recent decades. Therefore, we believe that higher water discharge can only be caused by dam failure (Figure 6).

Figure 7 shows the result of the flood simulation at $Q^{(max)}$ = 45,000 m${}^3/$s at time $t=15$ h for two models $n_M$. Red color highlights areas that are flooded in the model $n_M\ne$ const, but remain without water when $n_M=$ const. Yellow color defines areas under water at $n_M=$ const, but not flooded at non-uniform $n_M$. Blue color indicates areas flooded in both cases. The upper left inset shows the result of flooding at time $t=16$ h: orange highlights the flooded areas in the model $n_M\ne$ const, blue is the result for $n_M=$ const. Some HSs and ERs (for example, No. 47) are flooded when $n_M=$ const. Similarly, the lower inset depicts calculations at the time $t=19$ h, when ERs No. 8, 41, 51 are flooded at $n_M=$ const; however, settlements and above escape routes are not flooded in the model $n_M\ne$ const. The upper right inset ($t=25$ h) demonstrates the opposite situation, when the calculation with $n_M\ne$ const gives the flooding of both the settlement and ER No. 72. Thus, the heterogeneity of the roughness coefficient seems to be an important factor in modeling the hydrological regime for solving the problem of organizing evacuation.

Thus, the emergency situation development model is a multi-layer geoinformation model of the floodplain, including layers of altitudinal topography, channel structure, road network, distribution of settlements with data on the number of inhabitants, as well as numerous layers of catastrophic floods maps, which are built based on the results of hydrodynamic simulations. Each such layer characterizes flooding at certain points in time (for example, hourly data). The topological structure of the road network is represented by the graph G.

3.3. Algorithm for Calculating the Flooding Moments of Settlements and Evacuation Routes

Real hydrographs of river floods $Q\left(t\right)$ have the form of a single-peak function or a sequence of peaks. The forecast $Q\left(t\right)$ becomes known shortly before the start of the flood. Therefore, direct hydrodynamic modeling for flood mapping can be difficult. In this regard, we calculate an approximate vector of flooding moments for settlements and evacuation routes based on <<base>> vectors of flooding moments, using the results of numerical simulations with constant hydrographs $Q_j(Q_j<Q_{j+1},j=1,...,n)$. The values of $Q_j$ and n are selected so that the approximation error does not exceed a given value (for example, 5 percent).

The predictive emergency hydrograph is stepped, forming a sequence of constant hydrographs $Q\left(t\right)=\left\{Q_i^c\right\}$$(i=1,\cdots ,m)$ of the same duration $\mathsf{\Delta }t$. At the first stage of the algorithm, we form a set of constant hydrographs ${\tilde{Q}}_i$$(i=1,\cdots ,m)$, given on the intervals $[t_0,t_o+i\mathsf{\Delta }t]$, whose water volume is equal to the volumes of the parts of the original hydrograph $Q\left(t\right)$ on the same intervals. The inequality $Q_j\le {\tilde{Q}}_i\le Q_{j+1}$ determines the search for the pair $Q_j,Q_{j+1}$ in the array $Q_j$$(j=1,\cdots ,n)$ for each $\widetilde{Q_i}$$(i=1,\cdots ,m)$. If the <<base>> vector for $Q_{j+1}$ contains numbers belonging to the time interval $[t_0+(i-1)\mathsf{\Delta }t,t_0+i\mathsf{\Delta }t]$, then they become elements of the moment vector of flooding settlements and road sections for $Q\left(t\right)$. These data allow us to calculate the vector of the maximum allowable duration of the safe evacuation of residents of each flooded settlement (

Table 1. The times of the flooding beginning of human settlements and roads during the first day of a large flood for the constant Manning coefficient ($n_M$ = 0.045) and the spatially distributed one ($n_M\ne$ const); $t_s$ is the flooding time of settlement; $t_r$ is the start time of flooding of the evacuation route from the corresponding settlement; L is the length of the evacuation route.

№	${\mathit{Q}}^{\left(\mathbf{max}\right)}$ = 35,000 m${}^3/$s				${\mathit{Q}}^{\left(\mathbf{max}\right)}$ = 60,000 m${}^3/$s				L, km
	${\mathit{t}}_{\mathit{s}}$		${\mathit{t}}_{\mathit{r}}$		${\mathit{t}}_{\mathit{s}}$		${\mathit{t}}_{\mathit{r}}$
	${\mathit{n}}_{\mathit{M}}$= 0.045	${\mathit{n}}_{\mathit{M}}\ne$  const	${\mathit{n}}_{\mathit{M}}$= 0.045	${\mathit{n}}_{\mathit{M}}\ne$  const	${\mathit{n}}_{\mathit{M}}$= 0.045	${\mathit{n}}_{\mathit{M}}\ne$  const	${\mathit{n}}_{\mathit{M}}$= 0.045	${\mathit{n}}_{\mathit{M}}\ne$  const
1	25	25	25	25	25	25	25	25	16.25
2	25	25	25	25	16	14	25	25	6.85
3	25	25	25	25	14	13	25	25	5.25
4	17	17	25	25	7	6	16	17	2.45
5	25	25	15	16	15	15	8	9	21.35
6	25	25	25	25	15	15	14	15	34.15
7	25	25	25	25	17	18	14	15	29.9
8	25	25	25	25	9	9	12	12	21.55
9	25	25	25	25	17	18	14	15	18.1
10	25	25	25	25	16	18	16	19	11.6
11	25	25	25	25	10	11	10	11	20
12	25	25	14	15	8	8	7	6	21.8
13	24	25	12	13	7	7	5	5	25.65
14	21	22	22	24	7	6	8	9	18.7
15	11	11	25	25	5	5	14	15	10.2
16	15	15	5	6	6	6	3	3	13.85
17	25	25	3	3	6	6	2	2	15.5
18	18	19	15	16	8	8	6	6	22.5
19	25	25	16	17	9	9	5	5	12.2
20	24	25	25	25	9	9	11	12	13.65
21	25	25	25	25	10	10	11	12	7.6
22	25	25	25	25	11	12	11	12	7.6
23	25	25	14	15	11	12	7	6	21.8
24	25	25	25	25	11	10	13	13	6.6
25	25	25	25	25	15	16	16	17	3.6
26	25	25	25	25	11	12	13	13	8.15
27	25	25	25	25	13	13	13	14	7.9
28	25	25	25	25	13	14	13	14	10.4
29	25	25	25	25	25	25	21	21	15.6
30	25	25	25	25	25	25	21	21	21.5
31	25	25	25	25	19	20	16	17	13.7
32	25	25	25	25	16	17	13	13	11.3
33	25	25	25	25	15	15	16	19	16.6
34	25	25	25	25	25	25	16	19	14.7
35	25	25	25	25	21	22	16	19	14.7
36	25	25	15	16	25	25	8	9	27.25
37	25	25	25	25	21	22	14	15	35.45
38	25	25	25	25	21	22	14	15	32.9
39	25	25	23	24	23	23	14	15	40.6
40	25	25	23	24	25	25	14	15	47.95
41	25	25	25	25	18	20	16	19	10.2
42	25	25	25	25	19	19	14	15	14.5
43	25	25	25	25	12	12	12	11	4.95
44	25	25	25	25	25	25	16	19	10.2
45	25	25	25	25	25	25	16	19	9.45
46	22	23	15	16	7	7	6	6	24.35
47	25	25	25	25	18	20	14	15	10.2
48	21	22	24	25	7	7	7	7	7.2
49	25	25	25	25	20	19	25	25	6.85
50	25	25	25	25	24	25	25	25	6.85
51	25	25	25	25	23	25	21	21	15.6
52	25	25	25	25	25	25	21	21	16.5
53	25	25	25	25	25	25	16	19	14
54	25	25	25	25	12	10	15	14	10.4
55	25	25	25	25	11	12	12	13	19.35
56	25	25	25	25	12	13	11	12	19.05

). This table contains the flooding moments of settlements and evacuation routes of the northern part of the floodplain. Distances to the nearest flooding points are also indicated. These data are computed from direct numerical hydrodynamic modeling for two maximum hydrograph values and two hydrological resistance models with constant and non-uniform Manning coefficients. These results of numerical simulations made it possible to formulate the problem of finding the minimum resource for the safe evacuation of the floodplain population.

3.4. Creation of Safe Schedule

An arbitrary evacuation schedule using the i-th STER is described by the evacuation matrix $\parallel m_{itj}\parallel$, ($i=1,\dots ,K_j$; $t=t^{\left(0\right)},\dots ,T$; $T=\underset{j}{max}{\tau }_{ij}$, $j=1,\dots ,n$), where $m_{itj}$ is the number of transport means that simultaneously evacuate the population from the j-th NP along its part of the i-th STER at the time step t, $t^{\left(0\right)}$ is the moment of evacuation start, ${\tau }_{ij}$ is the flood times of ERs, $K_j$ is the number of temporary EPs on the jth ER. An admissible (safe) schedule $S{h}^{sec}({\overrightarrow{t}}^{\left(0\right)},Q,R^{\left(lim\right)})$ is such that the elements of the evacuation matrix ensure the evacuation of the entire population of the flooded area before the onset of a dangerous situation. The safe schedule of the 1st stage satisfies the following condition

$$\begin{array}{cc}\hfill & {\tilde{M}}_{ij}^{\left(1\right)}=\sum _{p_{ij}=1}^{P_{ij}^{\left(1\right)}}{\tilde{m}}_{i{p}_{j}j},m_{itj}=m_{itj}^{\left(0\right)},m_{itj}^{\left(p_{ij}\right)}=m_{itj}^{(p_{ij}-1)}-{\tilde{m}}_{i{p}_{ij}j},\hfill \\ & {\tilde{m}}_{i{p}_{ij}j}=min\left(m_{i{p}_{ij}j}^{(p_{ij}-1)},\dots ,m_{i{p}_{ij}+{\theta }_{ij}j}^{(p_{ij}-1)}\right),\left[\frac{N_j}{a}\right]\le {\tilde{M}}_{ij}^{\left(1\right)},\hfill \\ & p_{ij}=1,\dots ,P_{ij}^{\left(1\right)},P_{ij}^{\left(1\right)}\le \widehat{{\tau }_j}-t_j^{\left(0\right)}-{\theta }_{ij}+1,\hfill \\ & i=1,\dots ,K,t=p_{ij},\dots ,p_{ij}+{\theta }_{ij}-1,j=1,\dots ,n,\hfill \end{array}$$

(7)

where a is the capacity of one vehicle, ${\theta }_{ij}\left(Q\right)=\theta (L_{ij}\left(Q\right)$ is the duration of one trip of an evacuation vehicle from the jth HS in its part ith STER with length $L_{ij}$ (forth-and-back). It is important that the values of $P_{ij}^{\left(1\right)}$ are defined ambiguously in (7). The minimum value (${\widehat{P}}_{ij}^{\left(1\right)}$) satisfies the inequalities

$$\left[\frac{N_j}{a}\right]\le {\tilde{M}}_{ij}^{\left(1\right)}\left({\widehat{P}}_{ij}^{\left(1\right)}\right)<1+\left[\frac{N_j}{a}\right],$$

which determines the actual duration of the 1st stage of evacuation from j-th HS by i-th STER.

The safe evacuation schedule of the second stage satisfies the condition

$$\begin{array}{cc}\hfill & {\tilde{M}}_{ij}^{\left(2\right)}=\sum _{p_{ij}=1}^{P_{ij}^{\left(2\right)}}{\tilde{m}}_{i{p}_{ij}j},m_{tj}=m_{itj}^{\left(0\right)},m_{itj}^{\left(p_{ij}\right)}=m_{itj}^{(p_{ij}-1)}-{\tilde{m}}_{i{p}_{ij}j},\hfill \\ \hfill & {\tilde{m}}_{i{p}_{ij}j}=min\left(m_{i{p}_{ij}j}^{(p_{ij}-1)},\dots ,m_{i{p}_{ij}+{\theta }_{ij}j}^{(p_{ij}-1)}\right),\left[\frac{N_j}{a}\right]\le {\tilde{M}}_{ij}^{\left(2\right)},\hfill \\ & p_{ij}=1,\dots ,P_{ij}^{\left(2\right)},P_{ij}^{\left(2\right)}\le R^{\left(lim\right)}-{\tilde{\theta }}_{ij}+1,\hfill \\ & i=1,\dots ,K;K=\prod _{j=1}^n{K}_j,\hfill \\ & t=p_{ij},\dots ,p_{ij}+{\tilde{\theta }}_{ij}-1;j=1,\dots ,n,\hfill \end{array}$$

(8)

where ${\tilde{\theta }}_{ij}=\theta \left({\tilde{L}}_{ij}\right)$ is the time of the evacuation of the population of the j-th NP from the temporary EP to the permanent EP along the second part of the i-th STER of length ${\tilde{L}}_{ij}$. The values of $\left(P_{ij}^{\left(2\right)}\right)$ in (8) are also not uniquely defined. The smallest of them $\left({\widehat{P}}_{ij}^{\left(2\right)}\right)$ satisfies the inequalities ${\displaystyle \left[\frac{N_j}{a}\right]\le {\tilde{M}}_{ij}^{\left(2\right)}\left({\widehat{P}}_{ij}^{\left(2\right)}\right)<1+\left[\frac{N_j}{a}\right]}$, which determines the actual duration of the 2nd stage of evacuation from the j-th HS.

We assume that the beginning of the 2nd stage of evacuation at each ER begins immediately after the end of the 1st stage. The set of beginning moments of the 2nd stage is described by the vector ${\overrightarrow{t}}^{\left(0\right)}$. A formal statement of the problem of minimizing the number of vehicles M on the sets of safe schedules (7), (8) and STER with parameters Q, ${\overrightarrow{t}}^{,\left(0\right)}$, $R^{\left(lim\right)}$ has the form:

$$\begin{array}{cc}\hfill & M_i^{min}({\overrightarrow{t}}^{\left(0\right)},Q,R^{\left(lim\right)})=\underset{t}{max}\left(M_{it}\right)⟶\underset{m_{itj}}{min},\hfill \\ \hfill & M_{it}=\sum _{j=1}^n{m}_{itj},\hfill \\ \hfill & t=t^{\left(0\right)},\dots ,T_i;T_i=\underset{j}{max}(T_{ij}^{\left(1\right)}+{\widehat{P}}_{ij}^{\left(2\right)});i=1,\dots ,K,\hfill \\ \hfill & M^{min}({\overrightarrow{t}}^{\left(0\right)},Q,R^{\left(lim\right)})=\underset{i}{min}{M}_i^{min}({\overrightarrow{t}}^{\left(0\right)},Q,R^{\left(lim\right)}).\hfill \end{array}$$

(9)

The moment of completion of the 1st stage of evacuation from the j-th HS along the i-th STER is equal to

$$T_{ij}^{\left(1\right)}=t_j^{\left(0\right)}+{\widehat{P}}_{ij}^{\left(1\right)}+{\theta }_{ij}-1.$$

(10)

The maximum duration of the 2nd stage of evacuation is equal to $R^{\left(lim\right)}$ for the sequential implementation of two stages of evacuation from the j-th HS in the case of $T_{ij}^{\left(1\right)}\le R^{\left(0\right)}$. The sequential implementation of the stages does not satisfy the system inequality (1) when $T_{ij}^{\left(1\right)}>R^{\left(lim\right)}$. Therefore, it is necessary either to implement both stages of evacuation in parallel, or to select the vector ${\overrightarrow{t}}^{\left(0\right)}$ so that the inequality $T_{ij}^{\left(1\right)}\le R^{\left(lim\right)}$ is satisfied automatically. This condition limits the number of admissible schedules and is necessary for the correct operation of the algorithm described below for solving the problem (7)–(9). We then consider the sequential implementation of two evacuation stages, with the second evacuation stage from each temporary EP starting immediately after the completion of the first stage. In this case, the matrixes (7) and (8) can be combined into one and taking $p_{ij}=T_{ij}^{\left(1\right)}+1,\dots ,T_{ij}^{\left(1\right)}+{\widehat{P}}_{ij}^{\left(2\right)}$ in (8).

The transport-pedestrian implementation of the 1st stage in the problem (7)–(9) requires the elimination of HSs located at a distance of less than 5 km from temporary EPs. The first stage is carried out on foot for the above HSs.

The result of solving problems (9) and (10) is the minimum allowable number of vehicles and the schedule of their movement, which provides a two-stage safe evacuation of the entire population. The safety of the 1st stage means the evacuation of the entire population of each human settlement to a temporary evacuation point before the flooding of the settlement and the evacuation route begins. The safety of the 2nd stage means the evacuation of the population from each temporary evacuation point to a permanent evacuation point in a time not exceeding the allowable value. The solution procedure is based on constructing a sequence of safe schedules and calculating the required number of transports. Each schedule in this sequence uses the entire allowed evacuation time, but the number of vehicles is reduced in each successive schedule until the limit is reached. Thus, the entire found sequence of schedules can be used for practical purposes.

The heuristic algorithm for solving the problem (7)–(9) is a multi-step process of redistributing vehicles between ERs throughout the evacuation. The initial admissible evacuation schedule $\parallel m_{itj}^{\left(0\right)}\parallel$ is the schedule in which each vehicle moves along a single ER. This schedule implements the 1st stage of evacuation from the j-th HS along the i-th STER, which began at the moment $t_j^{\left(0\right)}$ in a time not exceeding ${\widehat{\tau }}_j-t_j^{\left(0\right)}$. If the condition ${\widehat{\tau }}_j>R^{\left(lim\right)}$ is satisfied for some j, then only $t_j^{\left(0\right)}$ are considered that ensure the inequality ${\widehat{\tau }}_j-t_j^{\left(0\right)}\le R^{\left(lim\right)}$. Since the condition $T_{ij}^{\left(1\right)}\le R^{\left(lim\right)}$ is satisfied, the number of vehicles used at the 1st stage of evacuation at each EM ensures evacuation at the 2nd stage as well.

The first part of the algorithm implements the sequential distribution of a possible reserve of vehicles, formed at the end of the 1st stage of evacuation from the j-th HS along the i-th STER at the moment ${\widehat{\tau }}_j$, to other ERs. An improved admissible traffic schedule for the corresponding ER is compiled for each potential reserve assignment with the largest possible reduction in the number of cars in $\parallel m_{itj}^{\left(0\right)}\parallel$. If the maximum of this reduction is reached on a single ER, then the corresponding assignment is accepted. If the maximum is reached on several ERs (equal to $J_k$), then the improved admissible schedule is calculated for each of them and the moment of appearance of its reserve ${\vartheta }_{ikj}$ is determined, the value of which is equal to $\mathsf{\Delta }{m}_{ikj}$$(j=1,2,\dots ,J_k)$. The assignment is made to the route with the maximum reserve index ${\displaystyle I{R}_{ikj}=\frac{\mathsf{\Delta }{m}_{ikj}}{{\vartheta }_{ikj}}}$. If the maximum of this index is reached on several routes, then the reserve is assigned to the route with the minimum value of ${\vartheta }_{ikj}$. If these values are the same for several routes, then the assignment is made to any of them. The reserve of cars is distributed over several routes if it is redundant for one route. The algorithm assumes that the routes of movement of assigned vehicles to new routes during the time $\mathsf{\Delta }t$ are free from flooding. This is performed in particular when passing all paths through a single constant EP. If there are several such EPs, then this problem requires investigation. The complete release of evacuation transport for the population from the j-th HS occurs after the completion of the 2nd stage of evacuation. This transport is also distributed in the manner described above. The algorithm for constructing a sequence of improved schedules ($\parallel m_{itj}^{\left(0\right)}\parallel$, $\parallel m_{itj}^{\left(1\right)}\parallel$, ...) ends after the release of all vehicles.

The second part of the algorithm is a multi-step process of heuristic construction of admissible schedules with a further decrease in the total number of vehicles by one at each step. Since the number of vehicles is fixed, the admissible schedule is constructed by the sequential distribution of vehicles between evacuation routes at each moment of time in the direction from the end of the evacuation to its beginning. Routes with the longest allowable evacuation time have higher priority. Thus, the algorithm maximizes the reserve of machines at each time step. Such a reserve can be used on routes with minimum allowable evacuation times. The algorithm stops when it is impossible to construct an admissible schedule in this way.

3.5. Construction of the Function of Minimum Resource Support for Safe Evacuation and Selection of Evacuation Routes

The solution to the problem (2) is based on the construction of the hierarchical graph of ERs $\tilde{G}$. Directed edges of $\tilde{G}$ form paths from each HS to one permanent EP. The evacuation route graph $\tilde{G}$ is a subgraph of the graph G. Its edges define escape routes, which are the last to be flooded. We replace the maps $K_{fl}(t,Q)$ with the data array $\Lambda =\langle {\overrightarrow{\tau }}_1,\dots ,{\overrightarrow{\tau }}_n,{\overrightarrow{L}}_1,\dots ,{\overrightarrow{L}}_n\rangle$ in the tuple E using the ERs and emergency evolution models. The array $\Lambda$ contains the coordinates of each pair of vectors $({\overrightarrow{\tau }}_j,{\overrightarrow{L}}_j)$, which define the moments of flooding of road sections on the j-th evacuation route and the distance from them to jth HS.

Possible locations for temporary EPs are marked on each ER. The path from these EPs to at least one permanent EP is not flooded according to the flood development forecast. Temporary EPs are intended to be used when there is a lack of vehicles to safely evacuate the population to permanent EPs. The evacuation of the population to temporary EPs is carried out in two stages: (1) evacuation of all residents from HSs to temporary EPs, (2) evacuation from temporary EPs to permanent EPs. If the temporary EP is located close to the HS (no more than 5 km), then the 1st stage of evacuation can be carried out both by vehicle and on foot. We reduce the values of ${\vartheta }_{ij}$ by a factor of k when calculating the parameter $R_e$ in the case of placing temporary EPs in buildings.

Each collection of n temporary ERs forms a set of temporary evacuation routes (STER). Thus, the total number of STER is $K={\displaystyle \prod _{j=1}^n}{K}_j$, where $K_j$ is the number of temporary EPs on j-th ER. The maximum duration of the 1st stage of evacuation on the j-th ER is determined by the moment of flooding of its first section ${\widehat{\tau }}_j=\underset{i}{min}{\tau }_{ij}$$(i=1,...,I_j)$, where $I_j$ is the number of flooded areas on the j-th ER.

4. Results

The Minimum Resource Support Function for the Safe Evacuation of the Northern Part of the Volga-Akhtuba Floodplain Population

The Volga-Akhtuba floodplain (VAP) is located in the lower reaches of the Volga river between the Volga and Akhtuba rivers (see Figure 3). The discharge of the Volga river is regulated by the flow of water through the dam of the Volga Hydroelectric Power Plant (VHPP). The development of the federal road network in the floodplain has become a powerful stimulus for the socio-economic development and population growth of this region. The creation of a system for managing the hydrological safety of this territory is an urgent task in the context of the increasing urbanization of the NPVAF.

The map of the NPVAF road network with the numbers of its local sections (edges of the graph G) and the location of settlements are shown in Figure 8. The triangles mark bridges across the Volga and Akhtuba rivers, where permanent EPs are also located.

Table 2. Compositions of permanent evacuation routes.

Number of HS	Number of Permanent ER	List of Permanent ER Road Numbers	Number of HS	Number of Permanent ER	List of Permanent ER Road Numbers
1	1	66, 75, 26	27	25	38, 8, 9, 10
2, 49, 50	2	75, 26	28	26	19, 8, 9, 10
3	3	26	29, 51	27	65, 75, 26
4	4	10	30	28	67, 78, 65, 75, 26
5	5	55, 74, 77, 54, 68, 1, 0	31	29	63, 20, 9, 10
6	6	59, 17, 76, 16, 4, 3, 2, 1, 0	32	30	18, 7, 8, 9, 10
7	7	17, 76, 16, 4, 3, 2, 1, 0	33	31	74, 77, 54, 68, 1, 0
8	8	58, 57, 14, 0	34, 35	32	77, 54, 68, 1, 0
9	9	76, 16, 4, 3, 2, 1, 0	36	33	61, 55, 74, 77, 54, 68, 1, 0
10	10	80, 2, 1, 0	37	34	70, 60, 17, 76, 16, 4, 3, 2, 1, 0
11	11	32, 5, 4, 3, 2, 1, 0	38	35	60, 17, 76, 16, 4, 3, 2, 1, 0
12, 23	12	51, 49, 47, 46, 36, 7, 8, 9, 10	39	36	62, 70, 60, 17, 76, 16, 4, 3, 2, 1, 0
13	13	50, 69, 48, 46, 36, 7, 8, 9, 10	40	37	72, 62, 70, 60, 17, 76, 16, 4, 3, 2, 1, 0
14	14	49, 47, 46, 36, 7, 8, 9, 10	41, 44	38	68, 1, 0
15, 47	15	14, 0	42	39	53, 14, 0
16	16	44, 73, 43, 0	43	40	23, 10
17	17	45, 27, 0	45	41	2, 1, 0
18	18	69, 48, 46, 36, 7, 8, 9, 10	46	42	52, 69, 48, 46, 36, 7, 8, 9, 10
19	19	27, 0	48	43	43, 0
20	20	46, 36, 7, 8, 9, 10	52	44	78, 65, 75, 26
21, 22	21	42, 10	53	45	29, 3, 2, 1, 0
24	22	41, 9, 10	54	46	25
25	23	22, 10	55	47	35, 6, 5, 4, 3, 2, 1, 0
26	24	40, 9, 10	56	48	33, 5, 4, 3, 2, 1, 0

describes the composition of ERs. The algorithm for solving the problem (7)–(9) involves moving the released vehicles to new routes during the evacuation. The road connecting permanent EPs is not flooded by hydrographs at $Q\le$ 60,000 m${}^3$ s${}^{-1}$. Therefore, the paths of vehicles between ERs are free from flooding.

The calculation of the parameters of the evacuation situation was based on the dynamics of HSs and ERs flooding of the NPVAF area for $Q^{(min)}=5000$ m${}^3/$c, $Q^{(max)}\in$ [28,000 ÷ 60,000] m${}^3/$s with the construction of flood maps $K_{fl}(t,Q^{(max)})$ every hour. The modeling results show that the flooding of settlements starts at $Q^{(max)}\ge$ 29,000 m${}^3/$s and occurs during the first four days. Figure 9 shows examples of flood maps for several emergency parameter values.

Figure 10 shows the dynamics of flooding of ERs for various settlements (j) flooded in the first 27 h of the emergency. The time after the start of flooding is plotted on the horizontal axes. The distance from HS to flooded areas of ERs is shown on the vertical axes. Each breakpoint on the graphs corresponds to the appearance of a new flooding section of the evacuation route with the coordinates of the corresponding vectors of the array $\Lambda$. The duration of the 1st stages of evacuation is determined by the first coordinates of these vectors and their length is given by the last coordinates. Figure 11 shows the features of the process of flooding temporary evacuation point on the example of one road to understand the dependencies in Figure 10.

The search for an approximate solution of the problem (7)–(9) was carried out for the NPVAF, using a series of computational experiments with parameters $Q^{(max)}\in$ [30,000–60,000] m${}^3/$s; $t^{\left(0\right)}\in$ [0–8] h; $R^{\left(lim\right)}\in$ [6–30] hs for three STERs: 1-stage evacuation to permanent EPs (STER1); 2-stage evacuation with the shortest ERs of the 1st stage (STER2), the temporary EPs of which are located on the highway; 2-stage evacuation to temporary EPs located in the nearest settlements (STER3). In the latter case (STER3), the value of $R^{\left(lim\right)}$ changes by a factor of k, where the coefficient k is given by an expert.

All settlements are divided into five classes by types (transport, pedestrian) and time groups ($t_j^{\left(0\right)}=t^{\left(0\right)}$, $t_j^{\left(0\right)}=\widehat{{\tau }_j}-R^{\left(lim\right)}>t^{\left(0\right)}$) of the 1st evacuation stage for each fixed set of parameters $Q^{(max)}$, $t^{\left(0\right)}$, $R^{\left(lim\right)}$, STER2, STER3.

Table 3. Distribution of settlements and their inhabitants into classes by types and time groups of the 1st stage of evacuation for two STER and three values of the river discharges. The arrows show the dynamics of the composition of the classes with the growth of $Q^{(max)}$. The cells contain: “the number of settlements / the total number of inhabitants in these settlements”.

Class Definition	Class Dynamics	${\mathit{Q}}^{\left(\mathbf{max}\right)}=$ 35,000 m${}^3$/s	${\mathit{Q}}^{\left(\mathbf{max}\right)}=$ 45,000 m${}^3$/s	${\mathit{Q}}^{\left(\mathbf{max}\right)}=$ 60,000 m${}^3$/s
$K_0$: No evacuation	$K_0$→ ($K_1$, $K_2$)	31/10,015	19/4825	5/815
$K_1$: Pedestrian  evacuation  at  the 1st stage at $t_j^{\left(0\right)}={\widehat{\tau }}_j-R_j^{\left(lim\right)}>t^{\left(0\right)}$	$K_1$→($K_2$, $K_3$)	STER2: 9/3690 STER3: 8/3490	STER2: 6/3130 STER3: 5/2930	STER2: 9/2415 STER3: 7/2015
$K_2$: Transport evacuation at the 1st stage at $t_j^{\left(0\right)}=t^{\left(0\right)}$	$K_2$→$K_4$	STER2: 6/3375 STER3: 10/21,505	STER2: 14/7070 STER3: 18/25,200	STER2: 27/10,855 STER3: 30/28,885
$K_3$: Pedestrian  evacuation  at  the 1st stage at $t_j^{\left(0\right)}=t^{\left(0\right)}$	$K_3$→$K_4$	STER2: 8/21,170 STER3: 7/3040	STER2: 12/21,660 STER3: 8/3530	STER2: 12/24,270 STER3: 9/6240
$K_4$: Transport evacuation at the 1st stage at $t_j^{\left(0\right)}={\widehat{\tau }}_j-R_j^{\left(lim\right)}>t^{\left(0\right)}$	$K_4$	STER2: 2/555 STER3: 3/755	STER2: 5/2120 STER3: 6/2320	STER2: 3/450 STER3: 5/850

contains the distribution of the number of settlements and their inhabitants by these classes for three hydrographs and $R^{\left(lim\right)}=24$ h. The composition of the 2nd classes changes monotonically with the increase in water discharge, and the composition of the 3rd classes shows non-monotonic behavior (see Table 3). Most of the population moves from classes $K_1$ and $K_2$ to classes $K_3$ and $K_4$ as $Q^{\left(max\right)}$ grows, which is due to the earlier start of evacuation.

The heuristic algorithm described above made it possible to construct the series of schedules $S{h}_i^{sec}(t^{\left(0\right)},Q^{(max)},R^{\left(lim\right)})$, $i=1,2,3$ for each of the three STERs and the corresponding approximate functions of the minimum transport provision for the safe evacuation of the population:

(1)

$M_1^{min}({\overrightarrow{t}}^{(0)},Q^{(max)})$ for 1-stage safe vehicle-pedestrian evacuation in NPVAF;

(2)

$M_2^{min}({\overrightarrow{t}}^{(0)},Q^{(max)},24)$ for 2-stage safe pedestrian evacuation of the population in the NPVAF;

(3)

$M_3^{min}({\overrightarrow{t}}^{(0)},Q^{(max)},k24)$ for 2-stage safe vehicle-pedestrian evacuation in NPVAF with accommodation of the temporary evacuation points in settlements, ${\displaystyle \left(k=3/2\right)}$ (see Figure 12).

The constructed schedules are conditionally safe, since the inhabitants of HS no.16 are excluded at $t^{\left(0\right)}>5$ and the inhabitants of HS no.17 at $t^{\left(0\right)}>6$ (see Table 2).

The value of $R^{\left(lim\right)}$ can decrease under adverse weather conditions and increase under favorable ones. This affects the minimum level of transport support for the 2nd stage of evacuation. Figure 13 shows the dependencies $M_2^{min}(t^{\left(0\right)}=0,Q^{(max)},R^{\left(lim\right)})$ and $M_3^{min}(t^{\left(0\right)},Q^{(max)},R^{\left(lim\right)})$.

Figure 12 shows that the minimum number of machines corresponds to the function $M_3^{min}({\overrightarrow{t}}^{(0)},Q^{(max)},24)$, although its difference from the function $M_2^{min}({\overrightarrow{t}}^{(0)},Q^{(max)},24)$ is insignificant. Within our evacuation model, this property is preserved for $k\ge 1.1$ over the entire range of parameters $t^{(0)}$, $Q^{(max)}$, $R^{(lim)}$ indicated above. Thus, the solution to the problem (7)–(9) is the function $M_3^{min}({\overrightarrow{t}}^{(0)},Q^{(max)},R^{(lim)})$ in the considered range of parameters $t^{(0)}$, $Q^{(max)}$, $R^{(lim)}$. A slight difference from the function $M_2^{min}({\overrightarrow{t}}^{(0)},Q^{(max)},R^{(lim)})$ can be explained by the fact that STER3 is provided at the 1st stage by a larger number of vehicles than STER2, and at the 2nd stage by a smaller number.

5. Discussion

Comparison of the results of the algorithm’s implementation of 1-stage vehicle evacuation [79] and 2-stage vehicle-pedestrian evacuation in this work showed that the use of 2-stage vehicle-pedestrian evacuation reduces the limit of transport support for safe evacuation by 2–3.5 times in the Northern part of the Volga-Akhtuba floodplain. We emphasize that the placement of temporary EPs in places closest to flooded settlements is appropriate only if the evacuation starts late. The placement of temporary EPs in settlements in other cases increases the values of the minimum transport resources of the 1st stage of evacuation to a lesser extent than it reduces them in the 2nd stage. However, we did not analyze the possibility of real placement of the evacuated population in social infrastructure facilities in the floodplain.

The critical shortage of time for organizing evacuation in the VAF occurs quite early. Therefore, it cannot be compensated for by the redundancy of the transport resource. At the same time, the number of settlements with rapid flooding is small. The best recommendation for such settlements is to embank or raise the flooded sections of the evacuation route to the required height.

Comparison of the results of hydrodynamic modeling with various roughness models shows differences in the forecasts of the flooding moments of some settlements and evacuation routes within 2 h (see Table 1). Moreover, flooding can be both earlier and later. A later start of an evacuation can have a significant impact on the evacuation organization. This indicates the relevance of creating heterogeneous roughness models for each specific area of potential flooding.

We did not consider the influence of social factors on the solution of our optimization problem. First, a significant effect is the dependence of the pedestrian evacuation time on the age and gender characteristics of the population, the season, the presence of panic, etc. [89,90,91,92,93,94,95]. Another significant effect is associated both with the refusal of part of the population to evacuate, and with the uncertainty in the behavior of people evacuating on their own, for example, when choosing evacuation routes and evacuation points [89,96,97]. The third factor is the problem of the intensity of the traffic flow of self-evacuation of the population by private transport (see transport models of self-evacuation [31,32]). The study of the influence of these factors on the evacuation organization and its resource provision is the subject of future research. The created algorithms and the results obtained can be the basis for future studies on the problem of evacuation of the population of floodplain territories during floods. The first goal of such research is to find the optimal organization of land evacuation of a part of the population in the event of a shortage of transport resources for the safe evacuation of the entire population. Then, the optimization aim is to minimize the life risk of the population remaining in the flooded area. The second goal is to search for and optimize hydrotechnical projects for protecting the area from floods in order to reduce the life risk of the population (embankment of settlements, protective dams, raising the roadbed sections). The third goal is to develop a plan for the safest urbanization of the floodplain.

6. Conclusions

We propose a method for estimating the minimum resource support for the safe evacuation of the population of floodplain areas under flood threat based on the mathematical model of surface water dynamics and minimizing the number of evacuation vehicles on a set of safe evacuation schedules. The method involves evacuation in two stages with a shortage of cars with the optimal placement of temporary evacuation points on evacuation routes. The used flood development model contains the multilayer geoinformation model of the territory with layers of relief, channel structure, road network, and settlements, as well as the numerical hydrodynamic model. We describe the two-level problem of finding optimal locations for temporary evacuation points on evacuation routes and finding the minimum number of vehicles that ensure the safe evacuation of the population. The heuristic algorithm for solving it has been developed and studied.

A feature of the flood dynamics in floodplain areas is the territorial and temporal distribution of flood threats, which makes it possible to replace the risk minimization problem with a more deterministic and simple task of finding an acceptable (safe) evacuation schedule with the available reserve of transport evacuation means. The uncertainty of this reserve is overcome by solving the problem of minimizing the number of vehicles on a set of safe evacuation schedules. We consider the flooding of the Northern part of the Volga-Akhtuba floodplain as a result of a strong water discharge through the dam as a testing ground for applying the developed method. Our results show that the construction of the function of the minimum resource support for safe evacuation for a specific floodplain area makes it possible to carry out such an evacuation with a shortage of transport resources using different moments at the start of flooding.

The possibility of practical implementation of admissible schedules of organized evacuation with minimal resource support requires their advanced preparation in a wide range of parameters of the emergency water discharge.

It is necessary to single out another danger for floodplains associated with heavy precipitation. Then water can enter the floodplain not only through the riverbed, but also through spatially distributed sources. Such situations require special modeling of the hydrological regime, but our approach to organizing the evacuation of the population remains unchanged.