User balance flow distribution-based damaged road network recovery order optimization method with maximum toughness as target

Abstract

The invention discloses a damaged road network recovery sequence optimization method with maximum toughness as a target based on user balance flow distribution, and aims to solve the problems of low recovery efficiency of a post-disaster traffic road network and lack of scientific basis for decision making. The method comprises the following implementation steps: firstly, constructing a damaged road network toughness evaluation framework, and defining a toughness index as a reciprocal of a product of recovery time and total transit time in combination with a traffic network topological structure and a flow demand; secondly, establishing damaged road network toughness optimization scene parameters, abstracting a road network into a connected graph, and defining a recovery scheme set; then, a mixed integer linear programming model with the maximum toughness as the target is constructed, a nonlinear problem is converted into a linear problem through secant approximation, and constraint conditions such as flow balance and travel requirements are set; and finally, performing global optimization solution on the model by using a high-performance solver to generate an optimal recovery scheme. According to the method, the road network toughness can be scientifically quantified, the recovery sequence is optimized, blindness of experience decision making is avoided, the passing efficiency of a post-disaster traffic system is rapidly improved, efficient and scientific decision support is provided for post-disaster traffic recovery, and the method has important practical application value.

Description

Technical Field

[0001] The present invention provides a MILP method for optimizing the restoration order of a damaged road network based on user balanced flow distribution with the goal of maximum resilience. The present invention relates to the field of traffic network restoration technology, and specifically to a method for optimizing the restoration order of a damaged road network based on user balanced flow distribution with the goal of maximum resilience. This method is suitable for rapid recovery planning of traffic networks after a disaster. Through mathematical modeling and optimization algorithms, the optimal restoration order of the damaged road network is determined to improve the overall resilience and traffic efficiency of the road network. The present invention combines traffic engineering, operations research and computer science, and uses a mixed integer linear programming model to quantitatively analyze road network restoration problems. It is applicable to fields such as urban traffic management, emergency response and infrastructure repair. Its technical core lies in evaluating the resilience of the road network through a user balanced flow distribution model, and optimizing the restoration order, thereby maximizing the traffic capacity of the road network in the shortest time, and providing scientific decision-making support for post-disaster traffic recovery. Background Art

[0002] After natural or man-made disasters, transportation networks often suffer severe damage, resulting in reduced traffic capacity and hindering rescue efforts and normal urban operations. Traditional road network restoration methods often rely on empirical evidence or simple priority rules, lacking systematic optimization and making it difficult to achieve efficient restoration in complex scenarios. Existing technologies often rely on static analysis for resilience assessment and restoration order optimization, failing to fully consider the dynamics of traffic flow and user behavior, resulting in a disconnect between restoration plans and actual needs. Furthermore, nonlinear road resistance functions and flow balance constraints increase the complexity of model solutions, making it difficult for existing methods to achieve both computational efficiency and global optimality. While recent research has attempted to combine user balance allocation models with optimization algorithms, the following issues remain: First, the definition of resilience is unclear, failing to quantify the impact of restoration order on road network efficiency; second, the nonlinear nature of the model makes it difficult to solve, hindering its application to large-scale road networks; and third, the lack of comprehensive optimization of restoration time and efficiency. Therefore, a scientific and efficient road network restoration order optimization method is urgently needed to enhance the rapid recovery capabilities of post-disaster transportation systems. In response to the above problems, the present invention proposes an optimization method based on mixed integer linear programming, which provides a new technical path for the restoration of damaged road networks. Summary of the Invention

[0003] (1) Objects of the present invention:

[0004] The purpose of the present invention is to provide a method for optimizing the restoration order of a damaged road network based on user balanced flow distribution, which solves the problems of low efficiency and poor adaptability in traditional methods by quantifying the resilience of the road network and optimizing the restoration order. First, the present invention defines a scientific resilience index based on the topological structure of the traffic network and the flow demand, quantifies the impact of the restoration order on the traffic efficiency of the road network, and provides a basis for decision-making. Then, the present invention establishes an efficient optimization model, transforms the nonlinear road resistance function and the flow balance constraint into a mixed integer linear programming problem, reduces the solution complexity, and ensures the applicability of the model in large-scale road networks. Finally, relying on a high-performance solver, the present invention can efficiently obtain the best restoration plan, ensure that the road sections that have the most significant impact on the traffic capacity of the road network are repaired first within a limited time after the disaster, minimize the traffic interruption time, and achieve rapid recovery decisions. Through the above goals, the present invention aims to provide a scientific and efficient technical means for post-disaster traffic recovery, and enhance the overall resilience and emergency response capabilities of the road network.

[0005] (2) Technical Solution: Based on the above theories and ideas, the present invention provides a method for optimizing the restoration order of a damaged road network based on user-balanced traffic distribution and aiming at maximum resilience, namely, a MILP-based method. The specific implementation steps are as follows:

[0006] Step 1: Construction of a damaged road network resilience assessment framework; Specifically, a comprehensive assessment of the damage to the traffic network topology and traffic demand after a disaster is conducted, and a method for assessing the resilience of the damaged road network is determined based on a user balance distribution model. First, a common road resistance function is given, which is the relationship between the travel time t of a road and the traffic flow f on that road. A common road resistance function is as follows, and this road resistance function is used as an example in this invention:

[0007] t=t0[1+α(f / c) β ] (1)

[0008] Next, we give the objective function of user balance flow allocation, which is the integral of travel time and flow:

[0009]

[0010] Assume that there are k damaged road sections, and the required recovery time for each damaged road section is τ k , each time a damaged road section is restored, the total travel time for all traffic demands on the entire road network should be reduced. Therefore, the inverse of the total travel time obtained by restoring each damaged road section multiplied by the cumulative recovery time can be used as the recovery resilience R of the entire road network:

[0011]

[0012] Step 2: Define the parameter representation of the damaged road network resilience optimization scenario;

[0013] Specifically, the traffic network is mathematically abstracted into a connected graph G = (V, E), where i, j∈V is the serial number of the road node, i, j = 1, 2, ..., |V| (|*| represents the total number of elements in the set), and the nodes correspond to intersections and key locations with traffic demand in the traffic network, such as schools, factories, hospitals, banks, office buildings, and residential areas. Set E' represents the set of damaged roads. The set K represents the set of repair sequence numbers for the damaged road sections, k∈K, k=1,2,…,|K|, and |K|=|E'|. The set P represents the set of repair plans, p∈P, and each element of the set corresponds to an arrangement of the repair order for the damaged road sections. Parameter h ijpk The value is 0 or 1, indicating whether the edge (i, j) is restored when the kth edge is restored for the recovery plan p. If it is restored, h ijpk =1, otherwise h ijpk =0; correspondingly, τ pk It represents the recovery time required to restore the kth edge for the recovery solution p. Each edge (i, j)∈E in the graph represents a road in the traffic network, t 0ij represents the free flow travel time of this edge, that is, the shortest time to pass this edge; c ij represents the capacity of edge (i, j), that is, the maximum number of vehicles that can pass through this road within a period of time (usually one day). For the damaged road section, t 0ij and c ij They represent the free flow time and capacity of the road section after repair. α is a parameter in the road resistance function, which is generally set to 0.15; β is a parameter in the road resistance function, which is generally set to 4; γ ij The parameters introduced when performing secant approximation on the impedance function are composed of α, β and c ij Calculated together. Parameters and They represent the slope and intercept of the secant line introduced when the impedance function is approximated by the secant line. Specifically, represents the slope of the dth secant line when the secant approximation is performed on the edge (i, j); represents the intercept of the dth secant line when the secant line approximates the edge (i, j). W represents the set of travel demands, w∈W, w=1,2,…,|W|. (s w ,t w ) represent the starting point and end point of the w-th travel demand, s w ,t w ∈V. Parameter q w represents the flow of the wth travel demand. M is a large number.

[0014] Step 3: Establish a mixed integer linear programming model for the restoration order of the damaged road network with maximum resilience as the goal. The decision variables are the flow rate on each road segment, the resilience corresponding to each restoration scheme, and the integral time obtained at each restoration step. The optimization goal is to maximize the restoration resilience (which is treated as minimum in the mathematical expression).

[0015] Specifically: In this problem, each time a damaged road section is restored, the total travel time will decrease. Our goal is to use the shortest possible time to restore those sections that will have the greatest improvement in road network efficiency after restoration. Since we define road network efficiency as total travel time, our goal is to restore those sections that will reduce the total travel time the most after restoration. Therefore, the objective function is defined as restoration time × total travel time. The goal is to make this objective function as small as possible. The corresponding resilience is defined as the inverse of this objective function.

[0016] Step 3.1: Define the decision variables for optimizing the restoration sequence of the damaged road network with the goal of maximizing resilience;

[0017] A non-negative continuous variable representing the flow of the wth pair of travel demands on the road segment (i, j) when the kth edge is restored for the restoration plan p;

[0018] f ijpk : a non-negative continuous variable, representing the total traffic flow on the road segment (i, j) when the kth edge is restored for the restoration plan p;

[0019] t ijpk : Non-negative continuous variable, indicating that when the total flow on the road section (i, j) is f ijpk The travel time on the road section;

[0020] y ijpk : non-negative continuous variable, an intermediate variable introduced when linearizing the nonlinear impedance function;

[0021] g ijpk : non-negative continuous variable, an intermediate variable introduced when linearizing the nonlinear impedance function;

[0022] T pk : a non-negative continuous variable, representing the total travel time of all traffic demands on the road network when restoring to the kth edge for the restoration plan p;

[0023] B p : a non-negative continuous variable, representing the inverse of the recovery resilience for the recovery scheme p;

[0024] Step 3.2: Establish an objective function for optimizing the restoration order of the damaged road network with the goal of maximizing resilience;

[0025]

[0026] Among them, B p is the inverse of the resilience of the recovery solution p, T pk It represents the total travel time of all traffic demands on the road network when restoring to the kth edge for the restoration plan p, which reflects the travel efficiency of the road network. pk It represents the time required to restore the kth edge for the restoration plan p. Minimizing their product reflects the goal of quickly improving the traffic efficiency of the road network.

[0027] Step 4: Coordinate factors such as flow demand, flow balance, and travel impedance to constrain the model. Steps 2 and 3 determined the model parameters, variables, and optimization objectives, respectively. This step constrains the model in multiple dimensions.

[0028] Step 4.1: Establish travel demand constraints;

[0029]

[0030]

[0031] in, For the wth pair of travel demands from starting point s w Outflow flow, The wth pair of travel demand flows into the starting point s w The first equation in constraint (5) indicates the flow rate from the starting point s w The sum of the outflows should be equal to the required flow q for the wth trip demand w , the second formula represents the flow into the starting point s w The total flow should be 0. Similarly, the first equation in constraint (6) represents the flow into the terminal t w The sum of the total flows should be equal to the flow q required for the wth travel demand w , the second formula represents the outflow end point t w The total flow should be 0.

[0032] Step 4.2: Establish traffic balance constraints for each node in the network;

[0033]

[0034] This constraint means that for each node in the network and for each travel demand, the flow into the node should be equal to the flow out of the node. This constraint ensures that the flow in the network will not disappear suddenly.

[0035] Step 4.3: Establish the total flow constraint for each edge in the network;

[0036]

[0037] Among them, f ijpk is the total traffic flow on the road segment (i, j) when the kth edge is restored for the restoration plan p, is the flow of the wth pair of travel demands on the road segment (i, j) when the kth edge is restored for the restoration plan p.

[0038] Step 4.4: Establish a flow upper bound constraint for unrepaired edges;

[0039]

[0040] Step 4.5: Establish the relationship constraint between the travel time and flow of each edge in the network;

[0041]

[0042] Among them, t ijpk Indicates that when the total flow on the road section (i, j) is f ijpk The travel time on this road section.

[0043] Step 4.6: Establish the total travel time constraint for all traffic demands on the road network. This constraint is also the objective function in the user balance allocation model.

[0044]

[0045] Step 4.7: For T pk With flow f ijpk The nonlinear relationship is approximated by secant;

[0046] Specifically: First define a new parameter γ ij and two new variables y ijpk and g ijpk :

[0047]

[0048]

[0049]

[0050] Then formula (11) can be written as formula (15):

[0051]

[0052] Finally we use a set of secant lines to y ijpk With g ijpk By approximating the relationship between , we can obtain the following linear constraints:

[0053]

[0054] Step 5: Model Solving—Solving the proposed mixed-integer linear programming model aims to efficiently obtain the optimal recovery solution. This step relies on high-performance commercial solvers (such as CPLEX and Gurobi) to perform global optimization on the model.

[0055] (3) Advantages and effects:

[0056] This invention can construct a mixed integer linear programming (MILP) model based on the actual traffic distribution characteristics, achieving accurate quantification and optimal recovery decision-making for post-disaster road network capacity recovery. Compared with existing technologies, this invention has the following beneficial effects:

[0057] ① Scientificity and accuracy: This invention uses a user balance distribution model and resilience quantitative indicators to accurately evaluate the impact of different restoration orders on road network traffic efficiency, avoiding the blindness of traditional empirical methods and providing a scientific basis for decision-making.

[0058] ② Efficiency and practicality: This paper adopts a mixed-integer linear programming model to transform complex nonlinear problems into efficiently solvable linear problems. Combined with high-performance solvers (such as CPLEX and Gurobi), it can generate optimal restoration solutions in a short time and is suitable for practical applications in large-scale road networks.

[0059] ③ Comprehensive optimization and globality: This invention comprehensively considers recovery time, traffic demand and traffic efficiency, takes maximizing resilience as the goal, and ensures that the sections that will most significantly improve the overall traffic capacity of the road network are repaired first, thereby quickly restoring traffic functions after the disaster and reducing social and economic losses. BRIEF DESCRIPTION OF THE DRAWINGS

[0060] Figure 1 This is a flow chart of a method for optimizing the restoration order of a damaged road network based on user balanced distribution and aiming at maximum resilience, as described in the present invention. DETAILED DESCRIPTION

[0061] (1) Project background

[0062] This implementation simulates travel demand on a small-scale road network in urban Beijing, obtaining network structural and travel demand parameters. Based on this, the model and method described above are tested, assuming that several road sections are damaged due to various reasons (such as heavy rain). The goal is to determine the optimal restoration sequence for the damaged sections. This small-scale road network consists of 22 nodes and 68 edges, including 6 damaged sections, resulting in a total of 182 pairs of traffic demands.

[0063] (2) Basic information

[0064] The specific data of the damaged road sections are shown in Table 1 below: Table 1 Specific data of damaged road sections Serial number starting point end Repair costs Repair time 1 10 16 1000 12 2 16 10 1000 12 3 16 17 2000 15 4 17 16 2000 15 5 17 19 1000 18 6 19 17 1000 18

[0065] (3) Optimize calculation results Table 2 The effect of restoration order obtained by optimization calculation Restore order starting point end Recovery effect (total travel time) 1 16 10 4179271 2 10 16 4120190 3 17 19 4088597 4 19 17 4054765 5 17 16 4022758 6 16 17 4003429

[0066] Table 2 shows the results of optimizing the restoration order using the proposed optimization method and model. The restoration effect represents the total travel time of the road network when the road segment is restored. This demonstrates that the proposed method fully adheres to constraints such as flow conservation, user balance, and restoration resilience, achieving the best restoration effect (minimizing total travel time) for the damaged road network and achieving successful restoration results.

[0067] (4) Comparative effect analysis

[0068] To demonstrate the superiority of this invention, we conducted a comparative analysis with the "single-step optimal" approach. This approach is a logical and easy-to-implement restoration strategy. Its core concept is to find the edge that reduces travel time the most for the entire network when restoring the first edge. Based on this, we then search for the edge that reduces travel time the most for the entire network among the remaining unrestored edges, continuing until all damaged sections are restored. This resulting restoration order is also a good result, but it is not optimal.

[0069] Table 3 lists the results of the "single-step optimal" method that conforms to intuitive logic: Table 3 Results of the “single-step optimal” method Restore order starting point end Recovery effect (total travel time) 1 17 19 4176363 2 19 17 4130256 3 17 16 4127834 4 16 17 4091081 5 16 10 4053750 6 10 16 4003429

[0070] Tables 2 and 3 show the difference in restoration order between the "optimized calculation" and "single-step optimal" methods, respectively. The restoration order in Table 2 reflects the interdependencies among multiple edges, thus exhibiting a stable and significant downward trend. In contrast, the "single-step optimal" method in Table 3, while selecting the currently optimal edge for restoration at each step, achieves a smaller reduction in total travel time. The final result is the same as in Table 2, but the intermediate performance is weaker. Table 2 achieves better results because it considers the interdependencies among road segments in the network. In real traffic networks, each segment is not independent, and restoring a single edge may have a cascading impact on traffic conditions on other segments. The method in Table 2, through global optimization and simultaneous consideration of the restoration effects of multiple edges, can more effectively reduce total travel time. The "single-step optimal" method, on the other hand, focuses solely on the currently optimal solution, ignoring these interdependencies, resulting in local optimization that fails to achieve a global optimum. Therefore, the method in Table 2 performs better in optimizing overall travel time, demonstrating the importance of comprehensively considering network interdependencies.

[0071] In summary, the present invention proposes a method for optimizing the restoration order of damaged road networks based on user-balanced traffic distribution. By establishing a scientific resilience assessment framework and an efficient mixed-integer linear programming model, it solves the problems of insufficient decision-making basis and low efficiency in traditional post-disaster traffic recovery. This method innovatively quantifies the resilience of the road network as the inverse of the product of the recovery time and the total travel time, and uses the secant approximation technique to transform complex nonlinear problems into linear programming problems that can be efficiently solved, and combines high-performance optimization algorithms to quickly generate the optimal restoration plan. The implementation results show that this method can accurately evaluate the impact of different restoration orders on the traffic efficiency of the road network, and give priority to repairing the key sections that have the most significant impact on improving the overall traffic capacity, thereby significantly shortening the traffic interruption time and improving the overall resilience of the road network. The present invention provides a scientific and effective decision-making support tool for the rapid recovery of post-disaster transportation systems, and has important application value in the fields of emergency management and urban traffic planning.

[0072] The present invention has the following core innovations:

[0073] ① Resilience Assessment Framework: This framework uses the inverse of the product of restoration time and total travel time as a resilience metric for the first time, quantifying the dynamic impact of restoration order on network efficiency. A comprehensive resilience assessment method is constructed using a user balance flow model, combined with road resistance functions and traffic demand, providing a theoretical basis for optimizing restoration order.

[0074] ② Multi-constrained linear optimization: To address the difficulty of solving traditional nonlinear models, this paper transforms complex nonlinear constraints into linear constraints through secant approximation and variable reconstruction, establishing an efficient mixed-integer linear programming model. This model not only retains the practical physical meaning of the problem but also significantly improves solution efficiency and applicability.

[0075] ③ Dynamic Adaptability: This approach proposes an optimization strategy aimed at maximizing resilience. By defining decision variables and constraints, it coordinates factors such as traffic demand, node balance, and travel impedance to ensure the global optimality of the restoration plan. This approach transcends the limitations of traditional static priorities and enables dynamic and adaptive road network restoration planning, providing an innovative solution for post-disaster traffic management.

Claims

1. A method for optimizing the restoration order of a damaged road network based on user-balanced traffic allocation and aiming at maximum resilience, comprising the following steps: (1) Construction of a damaged road network resilience assessment framework: A comprehensive assessment is conducted based on the damage to the traffic network topology and traffic demand after a disaster, and a damaged road network resilience assessment method is determined based on the user balance distribution model. (1.1) First, a common road resistance function is given, that is, the relationship between the travel time t of a road and the flow rate f on that road. A common road resistance function is as follows, and this road resistance function is used as an example in this invention: t=t0[1+α(f / c) β ] (1.2) Next, we give the objective function of user balance flow assignment, which is the integral of travel time and flow: (1.3) Assume that there are k damaged road sections, and the required recovery time for each damaged road section is τ k , each time a damaged road section is restored, the total travel time for all traffic demands on the entire road network should be reduced. Therefore, the inverse of the total travel time obtained by restoring each damaged road section multiplied by the cumulative recovery time can be used as the recovery resilience R of the entire road network: (2) Parameter definition of damaged road network resilience optimization scenario The traffic network is mathematically abstracted as a connected graph G = (V, E), where i, j∈V are the serial numbers of the road nodes, i, j = 1, 2, ..., |V| (|*| represents the total number of elements in the set). The nodes correspond to intersections and key locations with traffic demand in the traffic network, such as schools, factories, hospitals, banks, office buildings, and residential areas. The set E' represents the set of damaged roads. The set K represents the set of repair sequence numbers for the damaged road sections, k∈K, k=1,2,…,|K|, and |K|=|E'|. The set P represents the set of repair plans, p∈P, and each element of the set corresponds to an arrangement of the repair order for the damaged road sections. Parameter h ijpk The value is 0 or 1, indicating whether the edge (i, j) is restored when the kth edge is restored for the recovery plan p. If it is restored, h ijpk =1, otherwise h ijpk =0; correspondingly, τ pk It represents the recovery time required to restore the kth edge for the recovery solution p. Each edge (i, j)∈E in the graph represents a road in the traffic network, t 0ij represents the free flow travel time of this edge, that is, the shortest time to pass this edge; c ij represents the capacity of edge (i, j), that is, the maximum number of vehicles that can pass through this road within a period of time (usually one day). For the damaged road section, t 0ij and c ij They represent the free flow time and capacity of the road section after repair. α is a parameter in the road resistance function, which is generally set to 0.15; β is a parameter in the road resistance function, which is generally set to 4; γ ij The parameters introduced when performing secant approximation on the impedance function are composed of α, β and c ij Calculated together. Parameters and They represent the slope and intercept of the secant line introduced when the impedance function is approximated by the secant line. Specifically, represents the slope of the dth secant line when the secant approximation is performed on the edge (i, j); represents the intercept of the dth secant line when the secant line approximates the edge (i, j). W represents the set of travel demands, w∈W, w=1,2,…,|W|. (s w ,t w ) represent the starting point and end point of the w-th travel demand, s w ,t w ∈V. Parameter q w represents the flow of the wth travel demand. M is a large number. (3) A mixed integer linear programming model for the restoration order of a damaged road network with maximum resilience as the goal is established. The flow rate on each road section, the resilience corresponding to each restoration scheme, and the integral time obtained at each restoration step are used as decision variables, and the maximization of restoration resilience (which is treated as minimum in mathematical expressions) is the optimization goal. In this problem, the total travel time will decrease with each damaged road section restored. Our goal is to restore the road sections that will have the greatest improvement in road network efficiency as soon as possible. Since we define road network efficiency as total travel time, our goal is to restore the road sections that will reduce the total travel time the most after restoration. Therefore, the objective function is defined as restoration time × total travel time. The goal is to make this objective function as small as possible. The corresponding resilience is defined as the inverse of this objective function. Define decision variables for optimizing the restoration sequence of damaged road networks with the goal of maximizing resilience; A non-negative continuous variable representing the flow of the wth pair of travel demands on the road segment (i, j) when the kth edge is restored for the restoration plan p; f ijpk : a non-negative continuous variable, representing the total traffic flow on the road segment (i, j) when the kth edge is restored for the restoration plan p; t ijpk : Non-negative continuous variable, indicating that when the total flow on the road section (i, j) is f ijpk The travel time on the road section; y ijpk : non-negative continuous variable, an intermediate variable introduced when linearizing the nonlinear impedance function; g ijpk : non-negative continuous variable, an intermediate variable introduced when linearizing the nonlinear impedance function; T pk : a non-negative continuous variable, representing the total travel time of all traffic demands on the road network when restoring to the kth edge for the restoration plan p; B p : a non-negative continuous variable, representing the inverse of the recovery resilience for the recovery scheme p; (3.1) Establish an objective function for optimizing the restoration order of damaged road networks with the goal of maximizing resilience; in, B p is the inverse of the resilience of the recovery solution p, T pk It represents the total travel time of all traffic demands on the road network when restoring to the kth edge for the restoration plan p, which reflects the travel efficiency of the road network. pk It represents the time required to restore the kth edge for the restoration plan p. Minimizing their product reflects the goal of quickly improving the traffic efficiency of the road network. (4) Coordinate factors such as flow demand, flow balance, and travel impedance to constrain the model; (4.1) Establish travel demand constraints; (4.2) Establish traffic balance constraints for each node in the network; (4.3) Establish the total flow constraint for each edge in the network; (4.4) Establish the upper limit constraint on the flow of unrepaired edges; (4.5) Establish the relationship constraints between the travel time and flow of each edge in the network; (4.6) Establish the total travel time constraint for all traffic demands on the road network, which is also the objective function in the user balance allocation model; (4.7) For T pk With flow f ijpk The nonlinear relationship is approximated by secant, and a new parameter γ is defined first. ij and two new variables y ijpk and g ijpk ; (4.8) Thus, T pk Written as; (4.9) Using a set of secant lines for y ijpk With g ijpk By approximating the relationship between , we can get the following linear constraints; (5) Model solving: Solve the proposed mixed integer linear programming model to efficiently obtain the optimal recovery solution; this step relies on high-performance commercial solvers (such as CPLEX, Gurobi, etc.) to perform global optimization solutions to the model.

Images