Optimal placement of multiple types of detectors under a small vessel attack threat to port security

Transportation Research Part E 93 (2016) 71–94

Contents lists available at ScienceDirect

Transportation Research Part E journal homepage: www.elsevier.com/locate/tre

Optimal placement of multiple types of detectors under a small vessel attack threat to port security Xihong Yan a, Xiaofeng Nie b,⇑ a b

Higher Education Key Laboratory of Engineering and Scientific Computing, Taiyuan Normal University, Taiyuan 030012, China School of Mechanical and Aerospace Engineering, Nanyang Technological University, Singapore 639798, Singapore

a r t i c l e

i n f o

Article history: Received 14 January 2015 Received in revised form 5 May 2016 Accepted 7 May 2016

Keywords: Port security Small vessel attack Detector deployment Exact algorithms

a b s t r a c t We focus on a threat scenario where a terrorist would utilize a small vessel to attack a maritime target. We consider how to place multiple types of detectors to protect maritime targets from such an attack. Detectors are not perfectly reliable. The resulting detector placement problem is formulated as a nonlinear binary integer program such that the expected damage cost caused by the small vessel attack is minimized. Two exact algorithms and a greedy adding heuristic are proposed. Moreover, we conduct a detailed computational study and provide a case study in New York Harbor. Ó 2016 Elsevier Ltd. All rights reserved.

1. Introduction and literature review Driven by the economic boom in emerging markets, global maritime trade has grown significantly over the years. According to the Review of Maritime Transport 2011 of the United Nations, more than 80% of international trade in goods is being carried by sea and the 2010 international sea-borne trade was estimated to be 8.41 billion tons of goods (United Nations, 2011). As gateways for the movement of such a huge amount of goods, ports play critical roles to facilitate the flow. Since many manufacturing companies use just-in-time distribution systems to lower inventory holding costs, port shutdowns may have catastrophic economic impacts. For example, the closure of 29 U.S. West Coast ports in 2002 was estimated to cost the U.S. economy $1 billion to $2 billion per day (Gaudette, 2002). It was estimated that a nine-month labor dispute at the same 29 ports in 2014–2015 costed Honda 25,000 vehicles in February 2015 (Trott, 2015). After the tragic terrorist attacks of September 11, 2001, homeland security became a top priority in the U.S. Ports are extremely vulnerable to terrorist attacks due to their characteristics: size, open accessibility by water and land, location in metropolitan areas, the amount of materials being transported, and the ready transportation links to many locations (U.S. GAO, 2002). Because of their strategic importance and criticalness to economy, ports may become potential terrorist targets. Providing an effective and efficient port security system has become one of the primary tasks to safeguard the nation. One critical threat to port security comes from small vessels. Terrorists have shown a great interest in conducting attacks using small boats. Recently, the small boat threat has regained a lot of media attention, see, for example, USA Today (Frank, 2007), CBS News (CBS News, 2008), and GlobalPost (Neild, 2011). In April 2008, the U.S. Department of Homeland Security (DHS) announced the Small Vessel Security Strategy. In the report, the DHS identified four serious concerns associated with

⇑ Corresponding author. E-mail address: [email protected] (X. Nie). http://dx.doi.org/10.1016/j.tre.2016.05.003 1366-5545/Ó 2016 Elsevier Ltd. All rights reserved.

72

X. Yan, X. Nie / Transportation Research Part E 93 (2016) 71–94

using small vessels in terrorist-related activities: using small vessels as water-borne improvised explosive devices, to smuggle weapons, to smuggle terrorists, and as water-borne platforms to carry out stand-off attacks (U.S. DHS, 2008). In this paper, we concentrate on the threat scenario where a terrorist would utilize a small vessel as a water-borne improvised explosive device to attack a maritime target (for example, a military vessel, a cruise ship, a passenger ferry, and an oil tanker). Many up-to-date technologies are available to detect the presence of small amounts of radiological, chemical, and biological materials from hundreds of meters away. These devices include fixed-position detectors and mobile sensors. A detailed description of sensors to detect the presence of illicit materials on small boats is provided by Hill (2009). These devices are small and technologically advanced enough to be utilized at ports. For example, the DHS has announced a West Coast Maritime pilot program that involves developing a detection architecture to reduce the risk of threats transported on small vessels (U.S. DHS, 2007). The program is currently underway in Washington’s Puget Sound and California’s San Diego areas. We consider the situation where the government has collected intelligence information regarding the potential small vessel attack and plans to locate detectors to interdict the small vessel attack. As for the terrorist, he/she subsequently conducts the small vessel attack without knowing detector placement. Since various technologies (for example, electro-optical, electro-magnetic, radar, and seismic technologies) are involved, different types of detectors have their own costs and detection characteristics (for example, detection rate and effective detection radius). Our paper focuses on how to cost-effectively place multiple types of fixed-position detectors to help protect maritime targets in a port area from a terrorist attack using a small vessel. We assume that detectors are not perfectly reliable and the probability of detection depends not only on the type of the detector but also on how long the small vessel would stay in the detector’s effective detection area. The expected damage cost caused by the small vessel attack is minimized while keeping the cost of detectors within a budget limit. The resulting optimization problem is formulated as a nonlinear binary integer program. In the homeland security literature, many quantitative studies have focused on improving security from different domains, for example, aviation security and border security. In the aviation security domain, Nie et al. (2009a) incorporate joint responses of devices into passenger grouping strategies based on classification theory. Nie et al. (2009b) propose a mixed integer linear program to study how to group passengers with different risk levels. McLay et al. (2010) formulate a sequential stochastic multi-level passenger screening problem as a Markov decision process. Lee and Jacobson (2011) propose three probability-based performance measures for assessing sequential passenger assignment policies. Nie (2011) investigate risk-based grouping for a checked baggage screening system through a cost-effectiveness model. Nikolaev et al. (2012) consider dynamic passenger risk updates in a multi-stage sequential passenger screening problem. Lee and Jacobson (2012) propose three types of estimators for exploring passenger risk uncertainty. Nie et al. (2012) construct a simulation-based selectee lane queueing design framework to determine how to assign passengers with different risk characteristics to a selectee lane. In the border security domain, Morton et al. (2007) propose two types of stochastic network interdiction models for placing radiation sensors at border crossings. For the type under asymmetric information, Sullivan et al. (2014) introduce a smaller and tighter formulation on a bipartite network. Wein et al. (2009) introduce a mathematical optimization model for resource allocation such that a terrorist’s likelihood of successfully crossing the U.S.–Mexico border is minimized. Zhang et al. (2011) model a security-check queue at a U.S.–Canada border crossing as a two-stage queueing model. Most port security literature concentrates on container inspection. A container inspection system consists of a sequence of inspection sensors with different capabilities. Elsayed et al. (2009) focus on determining the threshold levels of inspection sensors and their corresponding sequence simultaneously. The objective is to minimize the overall system cost which includes both inspection and misclassification costs. Young et al. (2010) extend the model to incorporate the time required to complete inspection as another objective. Boros et al. (2009) consider finding container inspection strategies (represented as decision trees) such that the expected per-container-inspection cost is minimized. They develop a large-scale linear programming model based on the polyhedral description of the decision trees. Several evolutionary algorithms are proposed by Ramirez-Marquez (2008), Concho and Ramirez-Marquez (2010), and Van Weele and Ramirez-Marquez (2011). Merrick and McLay (2010) use multiple-objective decision analysis to examine whether screening cargo containers for smuggled nuclear threats is worthwhile. They conclude that the result relies on some key inputs. McLay et al. (2011) propose a risk-based framework to investigate how to define a system alarm. They develop a container reliability knapsack problem which decides what fractions of high-risk and low-risk containers undergo secondary screening. Gaukler et al. (2012) investigate two container inspection policies (a hardness control policy and a hybrid inspection policy) for detecting nuclear materials smuggling. Both policies are shown to perform better than the existing Automated Targeting System based system in most realistic situations. McLay and Dreiding (2012) introduce two linear programming models (a multi-level knapsack screening model and a multi-level threshold knapsack screening model) for detecting nuclear materials in cargo containers. Some structural properties for both models are provided and it is shown through a computational example that enforcing a threshold policy may not result in significant differences. Some related literature deals with the deployment of detectors in monitoring environments. Chakrabarty et al. (2002) consider determining the placement and type of sensors in a surveillance region such that the cost is minimized while achieving the desired coverage. An integer linear programming model is formulated and a divide-and-conquer approach is proposed for solving large problem instances. Kaplan and Kress (2005) study the operational effectiveness of suicidebomber detector schemes in two urban environments: a grid model and a plaza model. They conclude that even though suicide bomber sensors could play an important role in the defense of known targets, they may not be effective in random

X. Yan, X. Nie / Transportation Research Part E 93 (2016) 71–94

73

pedestrian attacks. Nie et al. (2007) consider the optimal placement of suicide-bomber detectors in a threat area where the locations of the potential targets are known. The objective is to minimize the expected number of casualties. Wilhelm and Gokce (2010) introduce an integer linear program to prescribe the types of sensors, the number of each type, and the location of each sensor in a surveillance system for port and waterway security. A branch-and-price decomposition method is formulated to solve the problem. Brown et al. (2011) propose a bi-level (a defender and an attacker) integer linear program to position patrol vessels for port radar surveillance. To solve large-scale problems, a decomposition scheme is developed. Danczyk and Liu (2011) introduce a linear binary integer programming model to optimally locate sensors along freeway corridors with the objective of minimizing performance measurement errors. The model turns out to be a resource-constrained shortest path problem which is solved through a customized branch-and-bound algorithm. In the setting of network traffic surveillance, Li and Ouyang (2011) build a linear binary integer program to decide sensors’ optimal locations such that the expected flow and path coverage benefits are maximized. A greedy heuristic and a Lagrangian relaxation based algorithm are proposed for large-scale problems. All previously mentioned literature on detector deployment concentrates on stationary sensors. Zou and Chakrabarty (2004) study how to reposition sensors such that coverage is improved from an initial random placement. Moreover, a practical virtual force algorithm is proposed. Focusing on a border surveillance system, Szechtman et al. (2008) investigate the patrol policies of moving sensors with the objective of minimizing the infiltration rate. They consider two types of sensor trajectories: leap-to-origin trajectories and back-and-forth trajectories. Hochbaum and Fishbain (2011) address the issue of deploying mobile sensors to identify and locate nuclear threats in an urban environment. They propose a concentrated alert problem and show that the problem is solvable in polynomial time. Zhu et al. (2014) investigate a mobile traffic sensor routing problem with the aim of enhancing traffic information acquisition benefits. The problem is formulated as a variant of the vehicle routing problem and a hybrid two-stage heuristic is proposed. In a series of papers, Wein et al. consider designing radiation detection–interdiction systems to protect cities from nuclear terrorist attacks. They focus on a threat scenario where a terrorist attempts to drive a nuclear or radiological weapon toward a target in a city center. They formulate three models: a sensor model (Wein et al., 2006), an optimal stopping problem (Atkinson et al., 2008), and a spatial interdiction model (Atkinson and Wein, 2008). Specifically, the sensor model determines the detection probability and the false positive probability in terms of sensor threshold levels. In the optimal stopping model, the terrorist determines whether to proceed or to detonate the bomb at each node along his route. The interdiction model proposes a spatial queueing model to determine the mean damage inflicted by the terrorist in terms of the number of interdiction vehicles. Based on these three models, Wein and Atkinson (2007) formulate a Stackelberg game where the leader is the government and the follower is the terrorist. The government chooses the number of sensors, the sensor thresholds, and the number of interdiction vehicles to minimize the expected damage caused by the terrorist. The terrorist then decides whether to move forward or to detonate the bomb with the objective of maximizing the expected damage. Our paper is closely related to Nie et al. (2007), Wilhelm and Gokce (2010), and Brown et al. (2011). In Nie et al. (2007), they consider how to place one type of suicide-bomber detector such that the expected number of casualties is minimized. We adopt their framework and extend it to incorporate multiple types of detectors in a maritime domain. Moreover, a much more efficient exact algorithm is proposed to solve the corresponding nonlinear binary integer program. Wilhelm and Gokce (2010) consider minimizing the cost of multiple types of sensors while keeping the probabilities of failing to detect an intrusion within preset bounds. Brown et al. (2011) consider minimizing the maximum probability that attackers can evade detection. Since Wilhelm and Gokce (2010) and Brown et al. (2011) involve only the probabilities of evasion, both models can be formulated as linear integer programs (a linear integer program and a bi-level linear integer program, respectively). But, the objective of our model is to minimize the expected damage cost, which is a nonlinear function and cannot be linearized through simple transformations. Our paper has three contributions. First, incorporating different characteristics of detectors (for example, detection rate, effective detection radius, and cost), we consider how to place multiple types of detectors to protect maritime targets. Second, we introduce two dominance definitions (strong dominance and weak dominance), investigate their relationship, and derive a property which can reduce feasible regions. Third, we propose another problem formulation which is amenable to a linear outer approximation, and design a linearized branch and bound algorithm which computationally proves to be more efficient than a standard branch and bound algorithm. The paper is organized as follows: The modeling is presented in Section 2. Section 3 concentrates on some properties. Two branch and bound based exact algorithms and one greedy adding heuristic are proposed in Section 4. A computational study is provided in Section 5 and a case study is provided in Section 6. Finally, conclusions and several future research directions are stated in Section 7. All proofs are in the appendices.

2. Modeling In a port area, there is a collection of valuable assets (for example, military vessels, cruise ships, and oil tankers) anchoring in the sea. These assets could be potential targets for a terrorist attack utilizing a small vessel as a water-borne improvised explosive device. Possible entrances for a small vessel are docks and marinas. A detector can sense whether or not a passing small vessel carries explosive materials. To protect maritime targets from such an attack, we consider how to deploy detectors. We assume that detectors are stationary and perfectly concealed (for example, attached to buoys). Different types of

74

X. Yan, X. Nie / Transportation Research Part E 93 (2016) 71–94

detectors may own different characteristics (for example, detection rate, effective detection radius, and cost). We investigate how to place multiple types of detectors in a port area such that the expected damage cost is minimized. Without loss of generality, we assume that the port area is rectangular. To simplify the modeling, the rectangle is modeled as a grid with equal-size squares. To model obstructions (for example, islands and infrastructures), we categorize squares into two types: blocked squares (a small vessel cannot travel through) and unblocked squares. Let D be the set of all blocked squares. There is a set of entrances, denoted by E. Within the rectangle, there exists a set of potential targets, denoted by S. All entrances and targets are located in the centers of some unblocked squares. An illustrative example is provided in Fig. 1 where a shaded region represents a blocked square, a triangle represents an entrance, and a star represents a target. In this example, D ¼ f13; 14; 18; 39; 40; 52g; E ¼ f32; 33; 58; 62g, and S ¼ f12; 29g. In this paper, we use square and square center interchangeably. A terrorist selects one entrance and steers a small vessel toward a chosen target. Traveling from a square center to another square center is allowed as long as the small vessel does not intersect any blocked square. We assume that the small vessel tries to travel on the straight line path for the selected entrance and target pair. If this straight line path happens to cross one or more blocked squares, the small vessel chooses a shortest path that is composed of a sequence of straight line segments not intersecting any blocked square. The shortest path can be constructed through some shortest path algorithm, for example, Dijkstra’s algorithm (Dijkstra, 1959). Continuing the example shown in Fig. 1, the shortest path from entrance 32 to target 12 consists of two straight line segments: 32 ! 20 and 20 ! 12. For the simplicity of exposition, we assume that the shortest path from entrance k to target j is unique and is labeled as P kj . Let ckj be the probability that the terrorist uses P P path P kj . Hence, ckj P 0; 8k 2 E and j 2 S, and k2E j2S ckj ¼ 1. These probabilities could be estimated based on collected intelligence information. We assume that the centers of squares (blocked and unblocked squares) are potential locations for detectors. Let H be the set of all detector types, indexed by h. For a type h detector, the corresponding effective detection radius, instantaneous detection rate, and unit purchasing cost are sh ; gh , and qh , respectively. In order to allow the interdiction team to have enough time to take action after detection, we define timely detection as detection such that there are at least k seconds remaining before the small vessel reaches the target. Let the driving speed of the small vessel be w m/s. A time period of k seconds converts to a distance of kw meters away from the target.

1

2

3

4

5

6

7

8

9

10

11

12

13 13

14

15

16

17

18

19

20

21

22

23

24

25

26

27

28

29

30

31

32

33

34

35

36

37

38

39

40

41

42

43

44

45

46

47

48

49

50

51

52

53

54

55

56

57

58

59

60

61

62

63

64

Fig. 1. An example and the associated shortest paths.

75

X. Yan, X. Nie / Transportation Research Part E 93 (2016) 71–94

Define N hkj as the collection of all squares from which a detector of type h can detect the small vessel traveling on path Pkj in a timely manner. The set N hkj can be obtained geometrically. First, we construct a timely detection zone which includes all points in the rectangle from which a type h detector can detect the small vessel on path P kj in a timely fashion (we need to check if the distance from one point to its nearest point on path P kj is less than sh and if the first point of detection is at least kw meters away from target j). Then, we check if a square center is within the timely detection zone. If yes, we include it into N hkj . Fig. 2 provides an example where the square size is 200 m  200 m, s2 ¼ 200 m, k ¼ 10 s, and w ¼ 20 m/s. In this example, N 232;12 ¼ f20; 21; 22; 23; 24; 28; 29; 30; 31; 32g. We assume that all types of detectors are not perfectly reliable and the probability of detecting the small vessel in a timely fashion depends on how long the small vessel stays within a specific type of detector’s effective detection radius. For every i 2 N hkj , we define phikj as the probability that a type h detector placed in square i detects the small vessel traveling on path P kj in a timely manner. The procedure of calculating phikj is as follows: First, we construct a sub-path of path Pkj , labeled as P kj , which excludes the portion within kw meters from target j. Then, we obtain the length of the part of Pkj that is within sh meters from square i and define this length as lhikj . If entrance k happens to be inside of the circle with radius sh centered at square i, we stretch the first segment of Pkj to intersect the circle. Finally, based on Przemieniecki (2000), we calculate phikj by using the following formula:

phikj ¼ 1  egh lhikj : The rationale of using the formula is that the probability of timely detection is monotonically increasing with respect to a detector’s instantaneous detection rate and the length of stay within its detection radius. Moreover, if the length of stay approaches zero (infinity), the timely detection probability will approach zero (one). Following the example shown in Fig. 2, l2;22;32;12 equals 349.790 m. Given that g2 ¼ 0:005, we obtain that p2;22;32;12 ¼ 1  e0:005349:790 ¼ 0:826. Define T h ¼ [kj N hkj , which includes all possible locations from which a type h detector can detect the small vessel in a timely manner. Moreover, defining T ¼ [h T h , we can confine potential detector locations to T. We assume that at most

1

2

3

4

5

6

7

8

9

10

11

12

13

14

15

16

17

18

19

20

21

22

23

24

l2,22,32,12 25

26

27

28

29

30

31

32

33

34

35

36

37

38

39

40

41

42

43

44

45

46

47

48

49

50

51

52

53

54

55

56

57

58

59

60

61

62

63

64

Fig. 2. An example for N hkj and detection length lhikj .

76

X. Yan, X. Nie / Transportation Research Part E 93 (2016) 71–94

one detector can be placed in one square and there is a budget M for procuring detectors. For each h and each i 2 T h , we define a binary decision variable xhi , which equals 1 if a detector of type h is placed in square i and equals 0 otherwise. We assume that if the small vessel is detected in a timely fashion, the intervention team foils the attack with probability h. There are two possible scenarios in which the attack via path Pkj is successful. One scenario is that the small vessel is not detected with probability

YY

ð1  phikj Þxhi :

h2Hi2N h

kj

The other one is that the small vessel is detected but the team fails to stop the attack with probability

"

1

YY

#

ð1  phikj Þ

xhi

ð1  hÞ:

h2Hi2Nh

kj

For each target j 2 S, let C j be the resulting damage cost if the attack is successful. Since the small vessel uses path P kj with probability ckj , we obtain the expected damage cost as follows:

8 k2E j2S

> :h2Hi2Nh

" ð1  phikj Þxhi ckj C j þ 1 

YY

ð1  phikj Þxhi

h2Hi2Nh

kj

9 > = ð1  hÞckj C j : > ;

#

kj

ð1Þ

Rearranging (1), we obtain the following expression:

" # XX YY ð1  hÞckj C j þ hckj C j ð1  phikj Þxhi : k2E j2S

ð2Þ

h2Hi2N h

kj

Inside the square brackets of (2), the first term ð1  hÞckj C j is constant and hence it can be removed without affecting the optimal solution. Moreover, for convenience, we define W kj ¼ hckj C j . The resulting detector placement problem is formulated as the following nonlinear binary integer program:

½DP

min

XX YY W kj ð1  phikj Þxhi k2E j2S

s:t:

X

ð3Þ

h2Hi2Nh

kj

xhi 6 1;

8i 2 T;

ð4Þ

h:i2T h

XX

qh xhi 6 M;

ð5Þ

h2H i2T h

xhi 2 f0; 1g;

8h 2 H; i 2 T h :

ð6Þ

The objective is to minimize the expected damage cost of the small vessel attack. Constraints (4) ensure that there is at most one detector in one square. Constraint (5) ensures that the cost of placing detectors is within the budget limit. Constraints (6) are binary constraints. 3. Properties Dominance has been frequently exploited to reduce the feasible region of an optimization problem in the network and discrete location literature. In this section, we first provide two definitions for dominance and then explore some properties based on the definitions. Definition 1 (Strong dominance). Consider two squares u and v. If phukj P phv kj for every path Pkj and every detector type h, and moreover there exists at least one path Pkj and one detector type h such that phukj > phv kj , we say that square u strongly dominates square v. Definition 2 (Weak dominance for detector type a). Consider two squares u and v and one specific detector type a. If paukj P pav kj for every path P kj , and moreover there is at least one path P kj such that paukj > pav kj , we say that square u weakly dominates square v for detector type a. From the two dominance definitions, we observe that if square u strongly dominates square v, square u weakly dominates square v for all types. But weak dominance for a specific type does not imply strong dominance. Moreover, even though square u weakly dominates square v for the type with the largest radius, square u may not weakly dominate v for other types. Hence, square u may not strongly dominate square v. A counter example is provided in Fig. 3.

77

X. Yan, X. Nie / Transportation Research Part E 93 (2016) 71–94

Example 1. Consider two detector types. A type b detector has detection radius sb ¼ 200 m and instantaneous detection rate gb ¼ 0:005, while for a type c detector, sc ¼ 180 m and gc ¼ 0:004. As shown in Fig. 3, we consider two squares u and v and one unique path Pk0 j0 which has one line segment within a circle centered at square u with radius 200 m and two line segments within a circle centered at square v with the same radius. We obtain that lbuk0 j0 ¼ 360:000; lcuk0 j0 ¼ 314:960; lbv k0 j0 ¼ 2  179:003 ¼ 358:006, and lcv k0 j0 ¼ 2  158:643 ¼ 317:268. Hence, we have that

pbuk0 j0 ¼ 0:835; pbv k0 j0 ¼ 0:833; pcuk0 j0 ¼ 0:716; and pcv k0 j0 ¼ 0:719: Based on these detection probabilities, we observe that square u weakly dominates square v for type b, but square v weakly dominates square u for type c. Let type h1 be the detector type with the largest radius. Even though in general weak dominance for type h1 does not imply strong dominance, the following theorem provides sufficient conditions under which the implication holds. Theorem 1. Consider the situations where for each ðk; jÞ pair the portions of path P kj within the circles centered at squares u and v with radius sh1 are at most two line segments. Moreover, for one line segment and two line segments with their intersection point not on the circle boundary cases, there are either two or four intersections with the circle. If the intersection point of two segments is on the circle boundary, there are three intersections with the circle (including the intersection point of the two segments). Suppose that square u weakly dominates square v for type h1 . For each ðk; jÞ pair, if one of the following two conditions is satisfied: (i) There is no portion of path Pkj within the circle centered at square v with radius sh1; (ii) The longer distance from square u to the lines passing through the two segments (if there is only one segment, the longer distance is the distance itself) is less than or equal to the shorter distance from square v to the corresponding lines, we have that square u strongly dominates square v. Based on the definition of strong dominance, we have the following theorem. Theorem 2. If square u strongly dominates square

v;

P

xH h:u2T h hu

P

P

xH h:v 2T h hv

in an optimal solution.

According to Theorem 2, we obtain the following corollaries without proof. Corollary 1. Suppose that the budget can afford to buy at most m detectors of the type with the cheapest unit price. If square v is ¼ 0 for every h 2 fh : v 2 T h g in an optimal solution. strongly dominated by at least m squares, xH hv Corollary 2. If there is a square which strongly dominates all other squares, and moreover the budget allows us to place at least one detector, we will place one detector in that square in an optimal solution. Let us consider a special case of Corollary 2. Suppose that the port area has only one entrance and several maritime targets, and moreover the entrance square strongly dominates other squares. Then, we will place one detector in the entrance square. H Corollary 3. Suppose that there is only one detector type b, that is, jHj ¼ 1. If square u strongly dominates square v ; xH bu P xbv in an optimal solution.

360.000

179.003

314.960

158.643

u

v

Fig. 3. A counter example.

78

X. Yan, X. Nie / Transportation Research Part E 93 (2016) 71–94

4. Algorithms Since the proposed [DP] in Section 2 is a nonlinear binary integer program, we provide two branch and bound based algorithms (a standard branch and bound algorithm and a linearized branch and bound algorithm) to obtain an optimal solution for a small- or middle-scale problem. For large-scale problems, it is implausible to obtain optimal solutions. So, we propose a greedy adding heuristic which can obtain good solutions in a short period of time. 4.1. Standard branch and bound algorithm When we relax the binary constraints for all xhi in [DP], the resulting relaxation problem is as follows:

½DP-R

min

XX k2E j2S

s:t:

X

W kj

YY

ð1  phikj Þxhi

ð7Þ

h2Hi2Nh

kj

xhi 6 1;

8i 2 T;

ð8Þ

h:i2T h

XX

qh xhi 6 M;

ð9Þ

h2H i2T h

0 6 xhi 6 1;

8h 2 H; i 2 T h :

ð10Þ

It is clear that [DP-R] is a convex nonlinear program. Hence, we can use a standard branch and bound algorithm to obtain an exact optimal solution for [DP]. For the ease of exposition, we use a vector X to represent a feasible solution, that is, X ¼ ðx1;1 ; . . . ; x1;jT 1 j ; . . . ; xjHj;1 ; . . . ; xjHj;jT jHj j Þ. The procedure of the standard branch and bound algorithm is as follows:

Algorithm 1. Standard branch and bound algorithm Step 0. (Initialization) Set stage s ¼ 0, active node set Os ¼ f0g, current node nðsÞ ¼ 0, and n ¼ 0. Let upper bound UB be the objective function value of the incumbent solution where all decision variables are equal to zero. Let e be some optimality tolerance. Solve [DP-R] and obtain a solution X0 . If all values of X0 happen to be integer, we stop the procedure and X0 is an optimal solution. Otherwise, go to Step 1. Step 1. (Division) Split node nðsÞ into two child nodes (labeled by n þ 1 and n þ 2, respectively) through the follow0 0 ing branching rule. Let XnðsÞ be the solution of node nðsÞ. Find index i such that XnðsÞ ði Þ is the first non-integer ele0 0 ment. Add a constraint Xði Þ ¼ 0 to the problem for node n þ 1 and Xði Þ ¼ 1 to the problem for node n þ 2. Then, we update the active node set by Os ¼ Os [ fn þ 1; n þ 2g n fnðsÞg. Step 2. (Bound updating and fathoming) Solve the relaxation programs of nodes n þ 1 and n þ 2. Update the incumbent solution and its objective function value UB. If the relaxation program is infeasible, or if the objective function value is larger than UBð1  eÞ, we fathom the node. Update Osþ1 by removing all fathomed nodes from Os . Let n ¼ n þ 2 and s ¼ s þ 1. Step 3. (Termination and node selection) If Os is an empty set, stop the procedure and return the incumbent solution. Otherwise, select a node nðsÞ from Os with the smallest objective function value and go to Step 1.

4.2. Linearized branch and bound algorithm For the [DP] formulation, the relaxation problems repeatedly solved during the procedure are nonlinear programs. In this subsection, we propose another detector placement problem formulation which is amenable to a linear outer approximation. The idea is borrowed from Sherali et al. (2004) where they deal with a nonconvex continuous program and the corresponding relaxation problems are linear programs. Here, we apply the idea to a nonlinear binary integer program and the associated relaxation problems are mixed integer linear programs. For each ðk; jÞ pair, we introduce two decision variables ykj and zkj as follows:

ykj ¼

YY

ð1  phikj Þxhi and zkj ¼ lnðykj Þ ¼

h2Hi2N h

XX

lnð1  phikj Þxhi :

h2H i2Nh

kj

kj

For the purpose of algorithm implementation, we need both lower and upper bounds for each ykj , that is, Lkj 6 ykj 6 U kj . These bounds will be modified during algorithm execution. Define L and U as vectors that include all lower bounds Lkj and upper bounds U kj , respectively. According to the expression of ykj , the initial bounds are given by:

L0kj ¼

YY

h2Hi2N h

kj

ð1  phikj Þ and U 0kj ¼ 1;

8k 2 E; j 2 S:

ð11Þ

79

X. Yan, X. Nie / Transportation Research Part E 93 (2016) 71–94

Tangential supports

kj kj

kj

kj

kj

Convex envelope

Fig. 4. A polyhedral outer approximation for zkj ¼ lnðykj Þ over ½Lkj ; U kj  (Sherali et al., 2004).

Based on the introduced new decision variables ykj and zkj , an equivalent formulation for [DP] is stated as follows:

½DPðL; UÞ

min

XX

ð12Þ

W kj ykj

k2E j2S

s:t:

X

8i 2 T;

xhi 6 1;

ð13Þ

h:i2T h

XX

qh xhi 6 M;

ð14Þ

h2H i2T h

zkj ¼

XX

8k 2 E; j 2 S;

lnð1  phikj Þxhi ;

ð15Þ

h2H i2N h

kj

8k 2 E; j 2 S; 8k 2 E; j 2 S;

ð17Þ

8h 2 H; i 2 T h :

ð18Þ

zkj ¼ lnðykj Þ; Lkj 6 ykj 6 U kj ; xhi 2 f0; 1g;

ð16Þ

We notice that [DP(L; U)] is a mixed integer linear program except that constraints (16) are nonlinear. To find an effective method to solve [DP(L; U)], we relax constraints (16) through constructing a polyhedral outer approximation for each ðk; jÞ pair which consists of a convex envelope and three tangential supports. As suggested by Sherali et al. (2004), we adopt three U L

kj , and U kj . The middle point ykj is the one that minimizes the maximum approximation support points: Lkj ; ykj ¼ lnðU kjÞlnðL Þ kj

kj

error. Fig. 4 shows such a polyhedral outer approximation. The corresponding relaxation problem of [DP(L; U)] is as follows:

½DPðL; UÞ-R

min

XX

ð19Þ

W kj ykj

k2E j2S

s:t:

X

8i 2 T;

xhi 6 1;

ð20Þ

h:i2T h

XX

qh xhi 6 M;

ð21Þ

h2H i2T h

zkj ¼

XX

lnð1  phikj Þxhi ;

8k 2 E; j 2 S;

ð22Þ

h2H i2N h

kj

lnðU kj Þ  lnðLkj Þ ðykj  Lkj Þ; U kj  Lkj ykj  Lkj ; 8k 2 E; j 2 S; zkj 6 lnðLkj Þ þ Lkj ykj  ykj zkj 6 lnðykj Þ þ ; 8k 2 E; j 2 S; ykj ykj  U kj zkj 6 lnðU kj Þ þ ; 8k 2 E; j 2 S; U kj

zkj P lnðLkj Þ þ

Lkj 6 ykj 6 U kj ; xhi 2 f0; 1g;

8k 2 E; j 2 S;

8h 2 H; i 2 T : h

8k 2 E; j 2 S;

ð23Þ ð24Þ ð25Þ ð26Þ ð27Þ ð28Þ

80

X. Yan, X. Nie / Transportation Research Part E 93 (2016) 71–94

1

2

3

4

5

6

7

8

9

10

11

12

13

14

15

16

17

18

19

20

21

22

23

24

25

26

27

28

29

30

31

32

33

34

35

36

37

38

39

40

41

42

43

44

45

46

47

48

49

50

51

52

53

54

55

56

57

58

59

60

61

62

63

64

65

66

67

68

69

70

71

72

73

74

75

76

77

78

79

80

81

82

83

84

85

86

87

88

89

90

91

92

93

94

95

96

97

98

99

100

101

102

103

104

105

106

107

108

109

110

111

112

113

114

115

116

117

118

119

120

121

122

123

124

125

126

127

128

129

130

131

132

133

134

135

136

137

138

139

140

141

142

143

144

Fig. 5. An illustrative example.

Table 1 Parameter values. Parameter

Description

Value

E S D

Set of entrances Set of targets Set of blocked squares Probability of choosing path P kj Damage cost at target j Probability of successful intervention Timely detection time Speed of the small vessel Total budget Unit purchasing price for type h Effective detection radius for type h Instantaneous detection rate for type h

f1; 4; 9; 48; 49; 84; 85; 135; 141g f77; 105g f16; 21; 50; 90; 95; 112g

ckj Cj h k w M

qh sh gh

1 18

35 million 0.6 10 s 20 m/s 200 thousand q1 ¼ 50; q2 ¼ 40; q3 ¼ 30 thousand s1 ¼ 220; s2 ¼ 200; s3 ¼ 180 m g1 ¼ 0:006; g2 ¼ 0:005; g3 ¼ 0:004

81

X. Yan, X. Nie / Transportation Research Part E 93 (2016) 71–94

We observe that [DP(L; U)-R] is a mixed integer linear program. Now, we present a linearized branch and bound algorithm to solve [DP] through repeatedly solving relaxation problems [DP(L; U)-R]. While doing so, L and U will be updated correspondingly. We define Y and Z as vectors that include all elements ykj and zkj , respectively. The procedure of the linearized branch and bound algorithm is as follows:

Table 2 Computational results for small-size problems: 12  12 and 13  13. Dim

E

S

D

SBB

12  12 72, 132, 139 8, 48, 135 3, 108, 141 8, 84, 144 1, 60, 120 2, 60, 96 60, 97, 132 7, 48, 140 5, 120, 137 3, 37, 96, 108, 120, 139 24, 48, 84, 133, 134, 137 2, 3, 11, 85, 132, 134 5, 72, 84, 136, 137, 141 5, 72, 134, 136, 139, 140 3, 12, 72, 133, 139, 144 8, 25, 96, 132, 136, 141 5, 11, 97, 135, 137, 143 7, 72, 96, 109, 135, 138 3, 12, 13, 25, 49, 61, 85, 138, 143 11, 12, 13, 24, 36, 49, 72, 137, 138 8, 9, 10, 13, 24, 109, 134, 135, 140 1, 4, 6, 11, 37, 109, 132, 133, 141 1, 9, 10, 11, 49, 109, 121, 132, 142 1, 5, 85, 108, 120, 121, 134, 135, 139 5, 6, 48, 61, 72, 97, 120, 138, 139 6, 48, 60, 72, 97, 133, 135, 138, 140 3, 6, 11, 36, 72, 132, 134, 139, 144

80 85 79 36, 86 31, 123 57, 128 5, 36, 139 57, 94, 131 8, 63, 125 125 101 50 93, 96 24, 73 45, 130 59, 61, 125 68, 91, 107 1, 31, 44 114 70 52 49, 81 19, 105 61, 69 47, 50, 106 70, 136, 141 59, 84, 137

59, 111 47, 102, 128, 137 26, 31, 38, 99, 128, 139 58, 115 33, 35, 57, 101 26, 29, 30, 34, 74, 121 54, 63 26, 61, 62, 129 53, 66, 99, 104, 108, 115 33, 85 35, 55, 96, 105 1, 27, 40, 110, 121, 128 132, 140 16, 18, 26, 124 26, 80, 99, 119, 141, 142 15, 100 12, 42, 124, 129 19, 22, 23, 58, 68, 87 39, 100 18, 41, 127, 143 26, 48, 50, 73, 94, 102 64, 129 86, 118, 128, 137 16, 60, 86, 99, 101, 114 89, 104 29, 61, 96, 120 70, 89, 109, 122, 129, 131

16.023 775.152 16.023 36.417 16.853 102.035 16.853 16.505 18.389 100.322 18.389 9.182 15.704 269.532 15.704 100.877 16.251 1520.866 16.251 203.273 15.957 129.769 15.957 90.306 15.995 251.096 15.995 170.265 19.024 35.870 19.024 20.291 16.868 320.453 16.868 140.680 15.524 66.584 15.524 44.880 19.388 13.533 19.388 35.359 17.078 39.498 17.078 27.133 18.208 328.471 18.208 79.712 21.606 476.305 21.606 275.510 23.985 1616.113 23.985 300.808 20.331 5580.966 20.331 1947.719 23.980 590.020 23.980 587.703 22.980 155.074 22.980 311.824 18.199 1135.613 18.199 69.831 20.868 33.467 20.868 30.684 19.842 1269.117 19.842 110.971 24.677 5444.539 24.677 817.478 20.958 2864.878 20.958 510.012 23.632 4227.774 23.632 820.818 21.727 16023.899 21.727 2541.104 20.022 766.308 20.022 438.222 23.470 1926.104 23.470 1535.137

13  13 157, 160, 163 11, 117, 167 6, 105, 118 91, 166, 168 9, 78, 161 12, 39, 131 40, 165, 168 3, 92, 164 6, 9, 162 4, 14, 27, 79, 165, 169 39, 130, 144, 158, 164, 165 3, 11, 66, 143, 157, 169 9, 26, 92, 105, 117, 162 2, 6, 65, 161, 163, 166 2, 14, 40, 78, 117, 131 7, 9, 92, 117, 118, 144 4, 11, 26, 53, 165, 169 1, 92, 118, 156, 162, 164 3, 4, 5, 9, 10, 11, 13, 65, 160 14, 26, 27, 78, 79, 118, 130, 162, 165 1, 6, 39, 52, 53, 65, 92, 144, 164 9, 13, 26, 65, 91, 131, 156, 161, 165 3, 8, 12, 39, 130, 159, 161, 164, 166 1, 4, 7, 13, 14, 53, 143, 159, 163 3, 7, 8, 9, 79, 92, 156, 159, 163 2, 8, 26, 66, 91, 117, 144, 159, 160 5, 6, 12, 78, 159, 163, 166, 168, 169

124 67 90 61, 142 10, 124 62, 93 9, 67, 79 50, 77, 123 15, 30, 161 101 36 106 45, 109 7, 156 21, 43 32, 148, 151 28, 135, 147 33, 81, 148 104 71 15 20, 61 123, 144 125, 133 33, 142, 167 73, 93, 165 40, 47, 142

52, 87 33, 134, 144, 151 2, 19, 57, 77, 84, 112 72, 132 31, 36, 98, 128 7, 54, 55, 91, 98, 138 11, 66 34, 95, 101, 159 4, 26, 27, 52, 114, 129 39, 121 3, 51, 123, 137 7, 57, 82, 92, 117, 135 7, 113 52, 54, 77, 96 16, 30, 41, 82, 114, 164 11, 25 77, 98, 116, 157 21, 37, 69, 112, 113, 124 20, 86 20, 47, 63, 65 32, 85, 96, 122, 155, 165 8, 97 16, 81, 136, 153 74, 85, 86, 93, 126, 157 73, 130 19, 87, 140, 163 31, 62, 77, 101, 105, 156

Value

LBB CPU (s)

Value

GAH CPU (s)

Value

CPU (s) GR (%)

16.123 16.853 18.525 15.836 16.274 15.957 15.995 19.024 16.908 15.700 19.388 17.078 18.421 21.841 24.129 21.497 24.481 22.980 18.808 20.868 19.876 24.702 21.466 23.753 22.053 20.329 23.470

0.000 0.031 0.016 0.031 0.047 0.031 0.078 0.031 0.062 0.031 0.016 0.031 0.062 0.125 0.109 0.203 0.172 0.187 0.062 0.047 0.062 0.234 0.234 0.172 0.296 0.218 0.296

0.625 0.000 0.739 0.844 0.140 0.000 0.000 0.000 0.235 1.129 0.000 0.000 1.172 1.090 0.604 5.738 2.089 0.000 3.346 0.000 0.172 0.103 2.424 0.513 1.499 1.531 0.000

14.740 13.829 14.740 25.882 14.753 14.683 59.915 14.683 28.139 14.882 14.335 55.702 14.335 19.505 14.335 15.858 3.413 15.858 41.398 15.858 17.042 446.370 17.042 97.041 17.083 15.864 517.295 15.864 79.706 15.915 15.816 83.411 15.816 93.296 15.816 16.117 655.498 16.117 195.964 16.272 15.893 10.170 15.893 69.074 15.977 15.987 134.824 15.987 33.596 15.987 15.907 363.410 15.907 89.948 15.907 16.936 201.011 16.936 47.852 16.936 18.871 1518.490 18.871 412.507 19.351 18.253 130.943 18.253 76.237 18.664 22.668 19.182 22.668 42.181 22.668 19.594 6673.649 19.594 934.747 19.634 18.925 2057.123 18.925 1364.440 19.160 18.777 1660.410 18.777 325.547 18.777 15.595 34.235 15.595 102.719 16.030 19.067 109.976 19.067 39.182 20.548 17.686 600.288 17.686 33.264 17.686 23.012 661.465 23.012 441.235 23.298 22.293 628.885 22.293 492.069 22.293 20.844 11844.955 20.844 1472.953 21.878 20.628 5336.797 20.628 1861.998 20.628 26.349 7081.304 26.349 1794.156 26.349 23.233 3106.415 23.233 975.715 23.233

0.000 0.016 0.016 0.031 0.031 0.016 0.078 0.094 0.047 0.047 0.047 0.062 0.094 0.140 0.109 0.203 0.265 0.203 0.047 0.094 0.078 0.218 0.296 0.250 0.374 0.499 0.515

0.091 1.354 0.000 0.000 0.241 0.322 0.000 0.966 0.529 0.000 0.000 0.000 2.541 2.252 0.000 0.206 1.242 0.000 2.789 7.765 0.000 1.243 0.000 4.963 0.000 0.000 0.000

82

X. Yan, X. Nie / Transportation Research Part E 93 (2016) 71–94

Table 3 Computational results for small-size problems: 14  14 and 15  15. Dim

E

S

D

SBB

LBB

GAH

Value

CPU (s)

Value

CPU (s)

Value

CPU (s)

GR (%)

14  14 7, 11, 99 4, 9, 85 3, 9, 188 70, 71, 183 184, 186, 195 6, 56, 188 28, 182, 190 3, 56, 113 15, 57, 169 6, 13, 14, 85, 183, 188 43, 56, 85, 127, 188, 189 1, 3, 6, 112, 187, 195 3, 42, 99, 168, 186, 187 2, 8, 14, 98, 126, 169 9, 14, 56, 57, 140, 182 5, 71, 84, 112, 127, 193 57, 127, 184, 185, 192, 196 4, 15, 113, 154, 169, 190 1, 4, 5, 28, 70, 154, 168, 185, 193 6, 12, 15, 71, 98, 99, 185, 188, 190 6, 9, 28, 57, 84, 154, 155, 187, 192 2, 5, 14, 28, 155, 182, 189, 190, 193 2, 11, 12, 57, 113, 169, 182, 189, 196 5, 6, 7, 84, 99, 169, 191, 192, 194 11, 12, 141, 182, 188, 189, 190, 195, 196 11, 14, 15, 56, 99, 113, 127, 185, 193 4, 5, 6, 11, 28, 141, 188, 194, 196

53 147 154 9, 155 117, 154 20, 151 4, 36, 169 100, 141, 147 76, 95, 184 70 104 162 11, 189 3, 171 105, 112 62, 95, 145 2, 18, 27 33, 48, 193 178 48 172 4, 85 97, 179 114, 164 62, 69, 170 31, 181, 191 99, 173, 190

17, 37 16, 111, 130, 182 7, 21, 53, 58, 114, 193 36, 167 36, 50, 64, 106 7, 22, 44, 132, 169, 170 115, 122 70, 91, 105, 175 38, 86, 151, 168, 187, 192 60, 114 13, 67, 143, 165 21, 56, 70, 138, 150, 151 65, 66 65, 104, 136, 161 68, 102, 138, 146, 165, 176 108, 159 23, 33, 56, 65 23, 38, 64, 100, 104, 150 69, 191 53, 69, 90, 128 23, 24, 29, 87, 93, 137 16, 119 14, 54, 109, 151 71, 74, 111, 123, 176, 177 45, 184 43, 71, 87, 154 13, 104, 133, 135, 153, 175

15.052 14.720 14.697 16.306 15.236 17.171 15.914 16.205 15.751 16.440 17.312 16.944 19.182 19.305 22.971 20.341 17.723 24.068 18.870 19.554 18.138 22.864 19.863 24.154 18.354 20.990 20.577

108.939 2.800 252.673 833.852 35.631 467.636 1575.392 225.307 277.033 208.823 73.537 561.907 4408.313 1023.131 1804.479 1119.694 51154.304 17720.082 3800.152 445.808 767.713 4884.007 1827.174 1547.033 4894.198 42577.062 18695.649

15.052 14.720 14.697 16.306 15.236 17.171 15.914 16.205 15.751 16.440 17.312 16.944 19.182 19.305 22.971 20.341 17.723 24.068 18.870 19.554 18.138 22.864 19.863 24.154 18.354 20.990 20.577

33.263 21.841 30.195 117.436 81.370 110.161 316.186 132.365 116.761 88.428 45.251 45.364 483.802 302.783 334.874 994.566 4111.577 2361.321 164.316 122.154 71.157 1367.643 331.710 597.161 576.444 3677.363 3254.772

15.052 14.720 14.697 16.306 15.443 17.171 16.036 16.205 16.094 16.735 17.312 16.944 19.932 19.814 22.998 21.542 17.816 24.206 19.281 19.867 18.735 22.919 19.863 24.154 18.598 21.805 21.074

0.000 0.016 0.016 0.047 0.016 0.047 0.156 0.062 0.062 0.062 0.031 0.047 0.218 0.234 0.125 0.250 0.406 0.250 0.109 0.109 0.140 0.359 0.421 0.250 0.452 0.655 0.562

0.000 0.000 0.000 0.000 1.361 0.000 0.765 0.000 2.172 1.790 0.000 0.000 3.910 2.640 0.118 5.906 0.524 0.572 2.175 1.601 3.289 0.241 0.000 0.000 1.328 3.885 2.418

15  15 166, 210, 215 15, 30, 195 45, 75, 181 3, 150, 214 13, 16, 136 3, 165, 219 10, 165, 180 5, 212, 216 4, 12, 76 60, 106, 165, 195, 215, 221 6, 8, 10, 60, 91, 225 90, 105, 136, 216, 218, 225 4, 6, 75, 180, 195, 215 3, 14, 15, 30, 135, 219 7, 210, 213, 214, 215, 223 11, 14, 150, 195, 214, 221 1, 3, 105, 211, 216, 219 1, 5, 106, 165, 210, 221 6, 7, 9, 75, 105, 212, 214, 215, 220 7, 12, 91, 150, 151, 165, 166, 218, 225 3, 10, 14, 45, 46, 218, 220, 221, 222 1, 3, 31, 45, 46, 91, 150, 214, 216 6, 150, 151, 211, 213, 216, 219, 222, 223 1, 4, 5, 8, 10, 166, 216, 221, 224 46, 60, 75, 105, 180, 215, 216, 221, 222 8, 14, 45, 105, 181, 213, 215, 221, 222 14, 16, 135, 150, 213, 214, 219, 222, 224

71 147 200 78, 135 25, 207 136, 147 25, 132, 148 4, 42, 124 9, 92, 181 159 18 41 54, 208 29, 109 41, 183 56, 142, 220 138, 149, 205 114, 116, 148 3 99 208 48, 110 198, 225 45, 223 88, 110, 196 61, 141, 220 168, 176, 211

218, 221 13, 79, 86, 221 11, 40, 60, 70, 138, 170 105, 158 28, 38, 110, 119 10, 14, 19, 67, 83, 216 142, 162 70, 156, 162, 172 56, 58, 93, 103, 106, 214 101, 169 25, 57, 67, 193 39, 54, 70, 83, 176, 187 16, 140 9, 43, 133, 207 13, 58, 93, 117, 154, 186 48, 139 37, 124, 151, 208 101, 120, 156, 168, 193, 214 148, 157 26, 71, 190, 204 6, 18, 41, 66, 166, 187 130, 225 45, 61, 71, 127 14, 31, 63, 83, 87, 126 106, 178 96, 127, 214, 223 20, 27, 69, 73, 83, 205

14.894 14.481 14.576 17.004 18.131 19.112 15.411 16.665 16.002 18.984 15.665 19.302 18.135 23.968 21.202 20.629 19.441 18.263 15.400 23.062 16.363 21.402 19.292 18.849 20.784 24.981 24.168

22.499 10.932 13.040 186.147 752.086 101.608 12.727 1454.017 21.309 313.635 13.897 30.436 1083.082 664.753 150.351 6968.771 2103.693 800.127 681.917 248.298 223.093 8634.507 4860.170 13323.602 51167.478 43098.449 18191.604

14.894 14.481 14.576 17.004 18.131 19.112 15.411 16.665 16.002 18.984 15.665 19.302 18.135 23.968 21.202 20.629 19.441 18.263 15.400 23.062 16.363 21.402 19.292 18.849 20.784 24.981 24.168

24.824 13.434 22.564 66.310 61.959 30.075 63.379 233.259 52.756 57.718 54.234 25.407 166.914 166.803 183.164 2226.363 434.739 443.893 158.080 119.566 98.096 472.322 621.953 1313.895 8438.560 4019.433 1752.683

14.894 14.481 14.576 17.004 18.131 19.112 15.411 16.829 16.002 20.317 16.660 19.302 18.135 24.681 21.202 22.197 19.612 18.737 15.400 23.062 16.363 21.647 19.740 19.563 21.187 24.981 24.168

0.031 0.000 0.016 0.062 0.094 0.062 0.016 0.109 0.062 0.047 0.047 0.078 0.140 0.156 0.172 0.203 0.343 0.172 0.125 0.109 0.140 0.234 0.203 0.484 0.593 0.796 0.562

0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.985 0.000 7.024 6.349 0.000 0.000 2.973 0.000 7.598 0.878 2.596 0.000 0.000 0.000 1.146 2.323 3.788 1.938 0.000 0.000

Algorithm 2. Linearized branch and bound algorithm Step 0. (Initialization) Set stage s ¼ 0, active node set Os ¼ f0g, current node nðsÞ ¼ 0, and n ¼ 0. Let the initial bounds be L0 and U0 with elements derived through (11). Solve [DP(L0 ; U0 )-R] and obtain an optimal solution C0 ¼ ðX0 ; Y0 ; Z0 Þ and the corresponding objective function value V½DPðL0 ; U0 Þ-R. We observe that X0 is feasible to [DP]. So, we let upper bound UB be the objective function value of the incumbent X0 evaluated through (3). Let e be some optimality tolerance. If V½DPðL0 ; U0 Þ-R P UBð1  eÞ, stop the procedure. Otherwise, go to Step 1.

83

X. Yan, X. Nie / Transportation Research Part E 93 (2016) 71–94 Table 4 Computational results for middle-size problems: 20  20 and 25  25. Dim

E

S

D

LBB

GAH

Value

CPU (s)

Value

CPU (s)

GR (%)

20  20

81, 180, 385 41, 80, 391 80, 200, 382 360, 388, 395 61, 260, 383 12, 61, 101 12, 40, 80 120, 383, 392 12, 15, 280 2, 3, 20, 381, 384, 395 17, 80, 81, 200, 321, 381 9, 18, 81, 100, 386, 392 18, 21, 180, 201, 393, 395 2, 14, 121, 201, 320, 399 16, 20, 241, 320, 381, 389 1, 3, 221, 321, 386, 388 17, 18, 81, 380, 385, 388 13, 221, 361, 385, 389, 397 9, 120, 160, 320, 321, 380, 384, 388, 393 5, 12, 320, 341, 361, 383, 385, 394, 395 140, 180, 181, 321, 386, 388, 389, 397, 399 5, 20, 40, 80, 201, 261, 341, 381, 388 9, 140, 180, 260, 382, 385, 388, 390, 393 8, 15, 20, 100, 121, 261, 281, 381, 390 17, 80, 141, 260, 341, 360, 387, 393, 394 8, 15, 17, 20, 180, 240, 281, 387, 397 6, 7, 12, 18, 20, 81, 120, 390, 400

328 246 137 350, 398 129, 198 67, 100 17, 29, 373 110, 118, 295 58, 192, 240 255 78 193 15, 32 177, 190 106, 143 16, 60, 359 56, 95, 238 167, 172, 183 22 227 354 25, 365 169, 297 81, 119 54, 178, 284 258, 302, 325 55, 182, 197

185, 323 59, 240, 355, 371 79, 120, 190, 257, 349, 386 4, 46 21, 52, 60, 127 192, 208, 246, 255, 294, 315 158, 229 11, 13, 200, 242 134, 142, 163, 164, 203, 295 22, 374 175, 190, 267, 298 93, 147, 174, 220, 298, 356 86, 189 24, 228, 291, 354 4, 7, 27, 80, 108, 319 140, 152 62, 129, 327, 377 34, 111, 145, 188, 255, 299 122, 232 32, 279, 286, 384 49, 162, 209, 219, 325, 336 220, 367 29, 244, 294, 376 111, 118, 120, 211, 343, 384 324, 340 33, 85, 158, 264 21, 38, 92, 93, 152, 275

16.014 16.014 16.014 15.682 16.030 15.257 15.323 16.157 16.020 17.178 18.300 18.090 19.787 24.184 17.791 18.045 19.290 20.540 15.643 17.152 16.872 24.141 23.023 19.049 25.711 22.011 25.901

42.158 51.781 38.490 51.969 126.238 49.255 130.322 133.502 254.363 37.538 24.066 51.158 657.357 296.224 385.649 1204.516 1273.457 3657.791 70.556 91.700 57.141 139.347 393.418 203.123 2559.565 3670.307 1881.195

16.258 16.258 16.187 15.752 16.259 15.503 15.323 16.310 16.020 18.537 19.168 18.090 20.670 24.184 18.000 19.033 19.739 21.016 15.764 17.152 16.872 24.547 23.360 19.049 26.061 22.590 26.057

0.031 0.031 0.016 0.016 0.078 0.047 0.125 0.218 0.078 0.109 0.109 0.094 0.328 0.374 0.328 1.248 0.686 0.546 0.468 0.172 0.109 0.686 0.640 1.014 1.825 1.888 1.700

1.521 1.521 1.081 0.447 1.426 1.614 0.000 0.948 0.000 7.911 4.742 0.000 4.463 0.000 1.174 5.479 2.329 2.321 0.772 0.000 0.000 1.682 1.463 0.000 1.359 2.631 0.601

25  25

15, 201, 551 151, 275, 604 21, 225, 605 6, 151, 610 19, 450, 605 51, 615, 618 10, 375, 500 3, 550, 613 3, 4, 300 25, 151, 201, 301, 602, 621 176, 225, 251, 451, 606, 624 15, 126, 300, 400, 614, 622 4, 175, 401, 476, 610, 611 50, 175, 376, 401, 426, 613 8, 151, 250, 251, 275, 601 11, 51, 76, 225, 476, 525 2, 13, 26, 300, 475, 602 50, 326, 351, 426, 526, 625 11, 12, 23, 151, 251, 275, 400, 605, 614 19, 526, 551, 604, 607, 608, 609, 619, 621 25, 50, 226, 475, 476, 550, 605, 608, 613 4, 7, 15, 18, 19, 25, 51, 450, 621 9, 19, 21, 226, 251, 300, 450, 605, 606 5, 200, 375, 501, 612, 614, 621, 622, 623 15, 18, 125, 450, 525, 613, 614, 621, 625 16, 17, 19, 51, 326, 425, 450, 526, 604 1, 7, 9, 10, 12, 526, 551, 619, 624

312 190 265 18, 389 132, 298 249, 254 344, 355, 441 156, 241, 247 299, 537, 553 265 163 309 182, 573 13, 524 128, 287 176, 195, 439 124, 141, 402 10, 439, 556 467 166 552 122, 460 320, 440 36, 205 153, 199, 570 412, 587, 618 121, 351, 483

79, 158 219, 337, 347, 438 63, 101, 305, 501, 526, 620 435, 514 210, 231, 571, 578 5, 63, 70, 256, 307, 452 450, 469 17, 341, 399, 562 177, 291, 436, 463, 543, 589 510, 581 320, 510, 556, 601 57, 216, 302, 429, 511, 518 160, 473 61, 128, 380, 549 14, 24, 59, 313, 317, 356 517, 518 227, 414, 515, 544 21, 130, 307, 371, 455, 571 102, 403 80, 227, 469, 591 73, 76, 112, 228, 297, 306 23, 26 81, 93, 169, 254 7, 62, 102, 193, 475, 490 198, 447 3, 100, 227, 462 70, 118, 153, 168, 253, 517

16.015 16.015 15.787 15.965 16.282 16.030 16.083 19.766 15.715 17.389 17.241 17.414 18.734 17.442 23.668 19.165 19.449 19.034 18.569 15.582 16.525 19.829 25.599 19.041 22.501 22.712 20.357

48.354 58.844 22.033 223.191 197.012 198.039 308.430 125.172 114.717 109.730 100.931 34.316 728.586 349.204 287.136 630.217 762.117 1008.547 136.671 176.828 83.982 286.342 880.312 1258.734 9043.404 9001.400 1642.575

16.119 16.247 15.823 16.444 16.282 16.261 16.263 19.766 15.715 18.063 17.462 18.852 19.264 18.112 23.668 20.502 20.130 19.678 18.833 15.582 17.959 20.747 25.683 19.400 22.791 22.970 20.909

0.031 0.031 0.047 0.156 0.187 0.172 0.250 0.452 0.390 0.172 0.187 0.172 0.780 0.671 0.468 1.529 1.888 1.654 0.406 0.328 0.437 1.732 1.092 1.888 3.276 2.792 3.307

0.651 1.449 0.226 3.001 0.000 1.438 1.121 0.000 0.000 3.874 1.282 8.259 2.828 3.837 0.000 6.975 3.497 3.380 1.418 0.000 8.676 4.631 0.326 1.889 1.287 1.132 2.711

Step 1. (Division) Split node nðsÞ into two child nodes n þ 1 and n þ 2 through the following branching rule. From an 0 0 optimal solution CnðsÞ ¼ ðXnðsÞ ; YnðsÞ ; ZnðsÞ Þ of [DP(LnðsÞ ; UnðsÞ )-R], we find indices k and j such that

     nðsÞ  nðsÞ nðsÞ  nðsÞ  zk0 j0  lnðyk0 j0 Þ ¼ max zkj  lnðykj Þ: k2E;j2S

h i h i h i nðsÞ nðsÞ nðsÞ nðsÞ nðsÞ nðsÞ Then, we split the bound interval Lk0 j0 ; U k0 j0 into two subintervals: Lk0 j0 ; yk0 j0 and yk0 j0 ; U k0 j0 , and add a subinterval for each child node. Update the active node set by Os ¼ Os [ fn þ 1; n þ 2g n fnðsÞg and denote the bounds for two child nodes as ðLnþ1 ; Unþ1 Þ and ðLnþ2 ; Unþ2 Þ, respectively. (continued on next page)

84

X. Yan, X. Nie / Transportation Research Part E 93 (2016) 71–94

Table 5 Computational results for middle-size problems: 30  30 and 35  35. Dim

E

S

D

LBB

GAH

Value

CPU (s)

Value

CPU (s) GR (%)

30  30 60, 180, 882 16, 541, 871 181, 721, 880 23, 570, 630 120, 890, 892 17, 61, 631 30, 877, 900 4, 26, 571 690, 691, 781 1, 10, 271, 690, 721, 876 1, 17, 30, 510, 660, 895 120, 330, 511, 870, 883, 895 24, 90, 120, 690, 840, 870 1, 24, 360, 721, 840, 890 2, 12, 211, 301, 780, 898 8, 210, 240, 480, 875, 891 5, 8, 450, 872, 873, 899 4, 17, 29, 90, 721, 889 15, 28, 361, 630, 780, 841, 880, 886, 891 18, 20, 121, 361, 481, 811, 874, 894, 895 4, 8, 9, 24, 61, 240, 331, 391, 691 90, 451, 481, 541, 571, 870, 881, 890, 893 4, 6, 8, 91, 391, 451, 570, 750, 881 10, 61, 210, 330, 420, 421, 690, 840, 897 11, 19, 21, 270, 271, 630, 780, 889, 890 3, 9, 19, 60, 151, 541, 879, 883, 896 12, 27, 28, 60, 841, 873, 884, 887, 891

495 345 280 357, 443 387, 689 175, 255 121, 567, 880 12, 50, 830 176, 438, 752 50 568 131 144, 298 353, 433 694, 737 56, 152, 254 154, 493, 500 423, 434, 519 532 160 284 287, 679 455, 682 294, 879 257, 276, 536 124, 561, 643 18, 501, 513

507, 572 159, 209, 260, 660 28, 97, 193, 320, 836, 869 349, 893 110, 433, 674, 836 64, 228, 431, 529, 727, 780 231, 737 300, 406, 418, 727 114, 192, 226, 266, 571, 696 444, 713 64, 338, 592, 805 154, 259, 566, 689, 705, 869 264, 445 287, 470, 658, 861 139, 461, 682, 777, 819, 896 24, 612 375, 432, 439, 677 253, 272, 338, 377, 514, 773 511, 760 528, 599, 810, 869 53, 123, 482, 589, 605, 708 383, 698 387, 423, 505, 661 245, 281, 289, 298, 784, 874 274, 404 73, 441, 475, 542 87, 141, 287, 370, 531, 651

14.578 14.973 15.226 14.945 15.455 16.037 16.360 16.312 16.429 16.212 18.180 16.151 16.526 20.616 17.626 19.540 19.233 20.564 18.867 17.780 21.358 20.706 24.012 21.089 23.201 27.085 24.910

35.281 35.014 36.483 83.675 145.378 137.113 330.327 297.993 149.664 69.919 41.652 16.666 214.836 837.639 246.488 1574.229 898.796 3027.228 65.043 127.739 76.275 992.294 968.134 475.532 7051.625 7733.953 833.841

14.578 14.973 15.226 14.945 15.597 16.258 16.360 16.312 16.429 16.212 18.946 16.803 16.526 21.002 18.239 20.876 20.782 21.317 19.363 17.913 21.358 21.172 24.161 21.089 23.472 27.652 24.910

0.047 0.047 0.047 0.047 0.094 0.187 0.983 0.577 0.406 0.328 0.172 0.374 0.156 1.170 1.061 3.104 2.044 2.168 0.577 0.421 0.343 1.810 1.888 3.619 4.165 6.848 5.320

0.000 0.000 0.000 0.000 0.921 1.376 0.000 0.000 0.000 0.000 4.213 4.036 0.000 1.873 3.479 6.841 8.057 3.658 2.629 0.747 0.000 2.247 0.620 0.000 1.170 2.092 0.000

35  35 13, 30, 34 350, 735, 805 4, 16, 210 14, 1121, 1200 34, 1015, 1202 23, 946, 1016 29, 315, 490 70, 1208, 1212 2, 6, 24 6, 456, 666, 1050, 1217, 1218 11, 26, 30, 280, 771, 1015 1, 140, 210, 770, 1195, 1198 15, 175, 315, 1120, 1212, 1221 23, 30, 386, 665, 910, 1216 32, 140, 666, 735, 1210, 1224 4, 16, 31, 71, 1219, 1220 28, 175, 350, 526, 911, 1212 7, 17, 70, 631, 700, 1219 16, 20, 175, 211, 490, 525, 1192, 1200, 1219 30, 32, 246, 806, 1086, 1193, 1197, 1215, 1223 2, 3, 9, 596, 665, 1015, 1085, 1198, 1210 105, 106, 175, 210, 385, 420, 911, 946, 1016 10, 16, 70, 245, 246, 1201, 1204, 1211, 1212 19, 24, 29, 32, 34, 141, 1086, 1212, 1223 8, 13, 70, 281, 456, 526, 1205, 1210, 1218 6, 22, 23, 36, 525, 595, 806, 1194, 1222 5, 14, 27, 28, 30, 490, 806, 1016, 1193

293 601 1139 692, 826 1115, 1193 104, 674 78, 450, 1157 836, 988, 1029 639, 644, 785 1197 696 1158 595, 803 176, 995 588, 622 439, 801, 812 39, 162, 200 557, 558, 739 693 221 571 38, 1075 27, 993 214, 569 373, 881, 1107 52, 634, 897 252, 381, 1111

285, 305 545, 599, 640, 1215 117, 404, 576, 606, 709, 836 396, 732 374, 408, 432, 980 86, 239, 506, 769, 984, 1010 110, 853 252, 322, 702, 800 401, 655, 672, 842, 855, 893 622, 692 208, 279, 822, 936 97, 450, 530, 626, 685, 728 177, 739 246, 258, 567, 1034 236, 345, 444, 812, 921, 940 85, 477 434, 450, 900, 1055 46, 53, 234, 355, 437, 516 775, 1209 87, 255, 810, 1156 409, 420, 505, 1138, 1157, 1225 727, 976 140, 353, 411, 443 172, 179, 561, 951, 1140, 1161 606, 761 699, 716, 830, 978 328, 411, 735, 1103, 1138, 1205

14.494 14.098 14.517 15.776 16.064 15.018 16.165 18.069 15.784 14.726 16.994 16.135 18.360 18.966 18.705 23.958 20.516 24.281 18.981 17.024 18.190 18.141 24.504 19.093 23.603 25.616 22.647

26.263 10.905 29.953 242.465 256.835 242.494 849.906 271.526 336.767 50.290 65.721 146.320 533.294 413.065 248.169 716.990 1127.391 1287.140 223.781 150.842 109.485 1421.988 1722.755 673.552 4110.776 2901.040 5592.256

14.494 0.031 14.098 0.078 14.517 0.140 15.792 0.218 16.925 0.374 15.018 0.234 16.237 0.983 18.069 0.562 16.056 0.406 14.726 0.234 16.994 0.296 16.461 0.374 19.561 0.546 20.033 2.387 19.364 0.562 24.486 3.713 22.018 3.526 24.538 3.276 19.032 0.749 17.181 0.998 18.849 0.671 18.517 3.120 25.028 3.697 19.728 3.292 24.277 9.251 25.731 10.421 23.579 11.591

0.000 0.000 0.000 0.102 5.360 0.000 0.447 0.000 1.725 0.000 0.000 2.022 6.539 5.629 3.522 2.202 7.319 1.059 0.269 0.920 3.625 2.072 2.139 3.328 2.857 0.450 4.117

Step 2. (Bound updating and fathoming) Solve the two relaxation programs: [DP(Lnþ1 ; Unþ1 )-R] and [DP(Lnþ2 ; Unþ2 )-R]. Evaluate the X parts of the solutions through (3) and update the incumbent solution and its objective function value UB. If the relaxation program is infeasible, or if the objective function value of the relaxation program is larger than UBð1  eÞ, we fathom the node. Update Osþ1 by removing all fathomed nodes from Os . Let n ¼ n þ 2 and s ¼ s þ 1. Step 3. (Termination and node selection) If Os is an empty set, stop the procedure and return the incumbent. Otherwise, select a node nðsÞ from Os with the smallest relaxation objective function value and go to Step 1.

X. Yan, X. Nie / Transportation Research Part E 93 (2016) 71–94

85

The linearized branch and bound algorithm with e ¼ 0 possesses a global convergence property. The convergence result is stated in the following theorem. Theorem 3. Consider any infinite branch whose nodes are labeled as fnðsÞg where s is in some index set N. Let sequence fXnðsÞ g be the X-variable parts of the corresponding sequential solutions generated by the linearized branch and bound algorithm with e ¼ 0. Then, fXnðsÞ g converges to a global minimizer of [DP]. 4.3. Greedy adding heuristic In this subsection, we introduce a greedy adding heuristic to quickly obtain solutions for large-scale problems. The idea is to sequentially choose pairs of ðh; iÞ such that the total expected damage cost is minimized given that the previously chosen pairs are fixed. The procedure of the greedy adding heuristic is as follows: Algorithm 3. Greedy adding heuristic Step 0. (Initialization) Let Q be the detector placement set. Initially, we set Q ¼ £. Step 1. (Termination checking) Check if M < minh2H qh . If yes, stop the procedure. Otherwise, go to Step 2. ~ # H be the set of detector types whose unit purchasing prices are less than or equal to M. Step 2. (Selection) Let H ~ ~ , we calculate the corresponding objective function value given ~ Define Q ¼ fðh; iÞ : h 2 H; i 2 T h g. For each ðh; iÞ 2 Q 0 0 that the placement set is Q [ fðh; iÞg. Then, we determine a pair ðh ; i Þ such that the calculated objective function value is minimum, that is, 0

0

ðh ; i Þ 2 arg min

n

o ~ : the objective function value of adding ðh; iÞ to Q : ðh; iÞ 2 Q 0

Step 3. (Updating) Let Q ¼ Q [ fðh ; i Þg; T h ¼ T h n fi g for every h such that i 2 T h , and M ¼ M  qh0 . Then, go to Step 1. 0

0

0

There are special situations where the greedy adding heuristic obtains optimal solutions. For example, if there is only one detector type and none of the sets N hkj intersect (that is, \kj N hkj ¼ £), as shown in Nie et al. (2007), the greedy adding heuristic yields optimal solutions.

Fig. 6. The satellite image of the area of study.

86

X. Yan, X. Nie / Transportation Research Part E 93 (2016) 71–94

325

326

327

328

329

330

331

332

333

334

335

336

337

338

339

340

341

342

307

308

309

310

311

312

313

314

315

316

317

318

319

320

321

322

323

324

289

290

291

292

293

294

295

296

297

298

299

300

301

302

303

304

305

306

271

272

273

274

275

276

277

278

279

280

281

282

283

284

285

286

287

288

253

254

255

256

257

258

259

260

261

262

263

264

265

266

267

268

269

270

235

236

237

238

239

240

241

242

243

244

245

246

247

248

249

250

251

252

217

218

219

220

221

222

223

224

225

226

227

228

229

230

231

232

233

234

199

200

201

202

203

204

205

206

207

208

209

210

211

212

213

214

215

216

181

182

183

184

185

186

187

188

189

190

191

192

193

194

195

196

197

198

163

164

165

166

167

168

169

170

171

172

173

174

175

176

177

178

179

180

145

146

147

148

149

150

151

152

153

154

155

156

157

158

159

160

161

162

127

128

129

130

131

132

133

134

135

136

137

138

139

140

141

142

143

144

109

110

111

112

113

114

115

116

117

118

119

120

121

122

123

124

125

126

91

92

93

94

95

96

97

98

99

100

101

102

103

104

105

106

107

108

73

74

75

76

77

78

79

80

81

82

83

84

85

86

87

88

89

90

55

56

57

58

59

60

61

62

63

64

65

66

67

68

69

70

71

72

37

38

39

40

41

42

43

44

45

46

47

48

49

50

51

52

53

54

19

20

21

22

23

24

25

26

27

28

29

30

31

32

33

34

35

36

1

2

3

4

5

6

7

8

9

10

11

12

13

14

15

16

17

18

2

2

2

1

Fig. 7. The grid of the case study.

5. Computational study In this section, we give a detailed computational study. During the implementation of both exact algorithms, we incorporate the strong dominance result in Theorem 2 to tighten our formulations. All algorithms are coded in Matlab and run on a PC with a 3.10 GHz Intel Core 2 Duo processor and 4 GB of RAM. To solve continuous nonlinear programs for the standard branch and bound algorithm, we use Matlab’s fmincon function with a sequential quadratic programming approach. To solve mixed integer linear programs for the linearized branch and bound algorithm, we run IBM ILOG CPLEX V12.4 in a Matlab environment. The optimality tolerances for both branch and bound based algorithms are set to 103 . The structure of this section is as follows: First, we provide an illustrative example. Second, computational results are provided to investigate the performance of the exact and heuristic algorithms proposed in Section 4.

X. Yan, X. Nie / Transportation Research Part E 93 (2016) 71–94

87

Table 6 Parameter values for the case study. Parameter

Description

Value

E S

Set of entrances Set of targets Probability of choosing path P kj

Cj h k w M

Damage cost at target j Probability of successful intervention Timely detection time Speed of the small vessel Total budget Unit purchasing price for type h Effective detection radius for type h Instantaneous detection rate for type h

f4; 12; 19; 270; 337; 338g f178; 208g c4;178 ¼ 0:005; c4;208 ¼ 0:095 c12;178 ¼ 0:108; c12;208 ¼ 0:018 c19;178 ¼ 0:096; c19;208 ¼ 0:086 c270;178 ¼ 0:100; c270;208 ¼ 0:116 c337;178 ¼ 0:083; c337;208 ¼ 0:122 c338;178 ¼ 0:053; c338;208 ¼ 0:118 C 178 ¼ 80; C 208 ¼ 100 million 0.6 10 s 20 m/s 350 thousand q1 ¼ 100; q2 ¼ 80 thousand s1 ¼ 500; s2 ¼ 460 m g1 ¼ 0:006; g2 ¼ 0:005

ckj

qh sh gh

5.1. Illustrative example In our illustrative example, we focus on a square threat area which is divided into 12  12 squares with each square of size 200 m  200 m. The sets of entrances, targets, and blocked squares are f1; 4; 9; 48; 49; 84; 85; 135; 141g; f77; 105g, and f16; 21; 50; 90; 95; 112g, respectively. The associated shortest paths are shown in Fig. 5. We assume that the terrorist chooses 1 every pair of entrance and target with the same probability, that is, ckj ¼ jEjjSj ; 8k; j. The damage cost for every target is the same with the value of 35 million, that is, C j ¼ 35 million; 8j. The probability h that the intervention team foils the attack is 0.6. The timely detection time k is 10 s and the speed of the small vessel w is 20 m/s. The total budget M for all types of detectors is 200 thousand. We consider three types of detectors where ðq1 ; s1 ; g1 Þ ¼ ð50; 220; 0:006Þ; ðq2 ; s2 ; g2 Þ ¼ ð40; 200; 0:005Þ, and ðq3 ; s3 ; g3 Þ ¼ ð30; 180; 0:004Þ. The unit of qh is in thousands and the unit of sh is in meters. All parameter values are summarized in Table 1. For the illustrative example, we find 131 pairs of squares such that one strongly dominates the other, for example, square 1 strongly dominates squares 13, 14, 15, 26, 27, 38, and 51. The optimal solution is x1;65 ¼ x1;68 ¼ x1;83 ¼ x1;88 ¼ 1 (that is, four type 1 detectors are placed in squares 65, 68, 83, and 88) and the corresponding objective function value is 25.980 million. The running CPU time for the standard branch and bound algorithm is 7867.074 s, while the CPU time for the linearized one is 1317.041 s. The example illustrates that the linearized branch and bound algorithm takes much less time to obtain the optimal solution. 5.2. Computational results Here, we continue to consider square threat areas with each square of size 200 m  200 m. We consider two categories of experiments: small-size test problems and middle-size test problems. For small-size test problems, the dimensions are 12  12; 13  13; 14  14, and 15  15 squares, while for middle-size test problems, the dimensions are 20  20; 25  25; 30  30, and 35  35 squares. For each dimension, the number of entrances may be 3, 6, or 9, the number of targets may be 1, 2, or 3, and the number of blocked squares may be 2, 4, or 6. The corresponding locations for entrances, targets, and blocked squares are randomly generated. In total, we obtain 27 instances for each dimension. Other parameter values are the same as those in the illustrative example. Computational results for small-size problems are provided in Tables 2 and 3. In both tables, Dim represents dimension, SBB represents standard branch and bound, LBB represents linearized branch and bound, and GAH represents greedy adding heuristic. GR represents the relative deviation of the greedy heuristic from an exact algorithm, which is computed by the following formula:

GR ¼

Heuristic Value  Exact Value  100%: Exact Value

From Tables 2 and 3, we observe that in most instances the CPU times of the linearized branch and bound algorithm are much shorter than those of the standard one. We find two instances for dimension 12  12, six for 13  13, two for 14  14, and seven for 15  15 where the CPU times of the linearized algorithm are longer than those of the standard one. But, for all these instances, both algorithms obtain optimal solutions within six minutes. As for the relative performance of the greedy adding heuristic, we observe that in 48 out of 108 instances it obtains optimal solutions. The average GR is 1.137% and the maximum GR is 7.765%. Computational results for middle-size problems are provided in Tables 4 and 5. Due to the preference of the linearized algorithm over the standard one for small-size problems, only the results for the linearized algorithm are provided. As for

88

X. Yan, X. Nie / Transportation Research Part E 93 (2016) 71–94

Table 7 Sensitivity analysis on the total budget. Budget

Optimal solution

Value

100 200 300 400 500 600 700 800

(1, 319) (1, 176) (1, 78) (1, 78) (1, 78) (1, 78) (1, 176) (1, 78)

75.296 60.480 50.702 46.687 43.526 41.231 41.100 39.513

(1, 227) (1, 176) (1, 176) (1, 158) (1, 158) (1, 227) (1, 158)

(1, 227) (1, 227) (1, 227) (1, 172) (1, 230) (1, 230)

(1, 319) (1, 230) (1, 227) (2, 78) (1, 319)

(1, 319) (1, 230) (2, 97) (2, 97)

(1, 319) (2, 190) (2, 190)

(2, 245) (2, 227)

(2, 319) (2, 270)

(2, 282)

Table 8 Sensitivity analysis on the type 2 detector characteristics. Characteristics

Optimal composition

Characteristics

Optimal composition

ðq2 ; s2 ; g2 Þ

Type 1

Type 2

ðq2 ; s2 ; g2 Þ

Type 1

Type 2

(140, 700, 0.008) (140, 700, 0.007) (140, 700, 0.006) (140, 700, 0.005) (140, 700, 0.004) (140, 600, 0.008) (140, 600, 0.007) (140, 600, 0.006) (140, 600, 0.005) (140, 600, 0.004) (140, 500, 0.008) (140, 500, 0.007) (140, 400, 0.008) (140, 400, 0.007) (140, 300, 0.008) (140, 300, 0.007) (120, 700, 0.008) (120, 700, 0.007) (120, 700, 0.006) (120, 700, 0.005) (120, 700, 0.004) (120, 600, 0.008) (120, 600, 0.007) (120, 600, 0.006) (120, 600, 0.005) (120, 600, 0.004) (120, 500, 0.008) (120, 500, 0.007) (120, 400, 0.008) (120, 400, 0.007) (120, 300, 0.008) (120, 300, 0.007) (100, 700, 0.005) (100, 700, 0.004) (100, 600, 0.005) (100, 600, 0.004)

2 2 2 2 2 2 2 2 2 3 2 2 2 3 3 3 1 1 1 1 1 1 1 1 1 3 1 1 1 3 3 3 0 1 1 3

1 1 1 1 1 1 1 1 1 0 1 1 1 0 0 0 2 2 2 2 2 2 2 2 2 0 2 2 2 0 0 0 3 2 2 0

(100, 400, 0.008) (100, 400, 0.007) (100, 300, 0.008) (100, 300, 0.007) (80, 700, 0.005) (80, 700, 0.004) (80, 600, 0.005) (80, 600, 0.004) (80, 500, 0.005) (80, 500, 0.004) (80, 400, 0.008) (80, 400, 0.007) (80, 400, 0.006) (80, 400, 0.005) (80, 400, 0.004) (80, 300, 0.008) (80, 300, 0.007) (80, 300, 0.006) (80, 300, 0.005) (80, 300, 0.004) (60, 700, 0.005) (60, 700, 0.004) (60, 600, 0.005) (60, 600, 0.004) (60, 500, 0.005) (60, 500, 0.004) (60, 400, 0.008) (60, 400, 0.007) (60, 400, 0.006) (60, 400, 0.005) (60, 400, 0.004) (60, 300, 0.008) (60, 300, 0.007) (60, 300, 0.006) (60, 300, 0.005) (60, 300, 0.004)

1 3 3 3 0 1 1 1 1 3 1 1 1 3 3 1 3 3 3 3 0 1 1 1 1 1 1 1 1 1 2 1 1 2 3 3

2 0 0 0 4 3 3 3 3 0 3 3 3 0 0 3 0 0 0 0 5 4 4 4 4 4 4 4 4 4 2 4 4 2 0 0

the relative performance of the greedy adding heuristic, we observe that it obtains optimal solutions in 31 out of 108 instances and the maximum running time among all instances is less than 12 s. The average GR is 1.936% and the maximum GR is 8.676%. 6. Case study We apply our model to a case study in New York Harbor, the traditional heart of the Port of New York and New Jersey. The satellite image of the area of study is shown in Fig. 6. Due to their strategic importance, let us choose the Statue of Liberty and the Brooklyn Cruise Terminal as the two potential targets. As shown in Fig. 7, the area of study is divided into 19  18 squares with each square of size 400 m  400 m. All shaded squares represent blocked squares and the sets of entrances and targets are E ¼ f4; 12; 19; 270; 337; 338g and S ¼ f178; 208g, respectively. The damage costs associated with the two targets are C 178 ¼ 80 million and C 208 ¼ 100 million, respectively. The shortest paths for each pair of entrance and target are

89

X. Yan, X. Nie / Transportation Research Part E 93 (2016) 71–94 Table 9 Sensitivity analysis on the probabilities ckj .

ckj

Situ 1

Situ 2

Situ 3

Situ 4

Situ 5

Situ 6

Situ 7

Situ 8

Situ 9

Situ 10

Situ 11

c4;178 c12;178 c19;178 c270;178 c337;178 c338;178 c4;208 c12;208 c19;208 c270;208 c337;208 c338;208

0.083 0.083 0.083 0.083 0.083 0.083 0.083 0.083 0.083 0.083 0.083 0.083

0.400 0.055 0.055 0.055 0.055 0.055 0.055 0.055 0.055 0.055 0.055 0.055

0.800 0.018 0.018 0.018 0.018 0.018 0.018 0.018 0.018 0.018 0.018 0.018

0.055 0.055 0.055 0.055 0.055 0.055 0.400 0.055 0.055 0.055 0.055 0.055

0.018 0.018 0.018 0.018 0.018 0.018 0.800 0.018 0.018 0.018 0.018 0.018

0.033 0.100 0.033 0.100 0.067 0.067 0.067 0.100 0.033 0.133 0.133 0.133

0.011 0.033 0.011 0.033 0.022 0.022 0.022 0.033 0.011 0.267 0.267 0.267

0.033 0.067 0.067 0.133 0.100 0.133 0.033 0.033 0.067 0.100 0.133 0.100

0.011 0.022 0.022 0.267 0.033 0.267 0.011 0.011 0.022 0.033 0.267 0.033

0.133 0.100 0.100 0.033 0.067 0.033 0.133 0.133 0.100 0.067 0.033 0.067

0.267 0.033 0.033 0.011 0.022 0.011 0.267 0.267 0.033 0.022 0.011 0.022

Opt Sol

(1, 176) (2, 78) (2, 119) (2, 227)

(1, 158) (2, 78) (2, 101) (2, 227)

(1, 158) (2, 4) (2, 101) (2, 227)

(1, 176) (2, 4) (2, 78) (2, 227)

(1, 60) (2, 4) (2, 78) (2, 319)

(1, 176) (2, 172) (2, 227) (2, 319)

(1, 227) (2, 208) (2, 230) (2, 319)

(1, 176) (2, 78) (2, 227) (2, 319)

(1, 176) (2, 230) (2, 245) (2, 319)

(1, 158) (2, 78) (2, 172) (2, 227)

(1, 172) (2, 4) (2, 119) (2, 158)

Value

51.101

46.882

38.255

50.952

45.923

51.913

48.881

48.222

42.455

50.641

46.500

obtained and the corresponding probabilities of choosing each path are estimated to be c4;178 ¼ 0:005; c4;208 ¼ 0:095; c12;178 ¼ 0:108; c12;208 ¼ 0:018; c19;178 ¼ 0:096; c19;208 ¼ 0:086; c270;178 ¼ 0:100; c270;208 ¼ 0:116; c337;178 ¼ 0:083; c337;208 ¼ 0:122; c338;178 ¼ 0:053; and c338;208 ¼ 0:118. There are two types of detectors where ðq1 ; s1 ; g1 Þ ¼ ð100; 500; 0:006Þ and ðq2 ; s2 ; g2 Þ ¼ ð80; 460;0:005Þ. The unit of qh is in thousands and the unit of sh is in meters. The probability h that the intervention team foils the attack is 0.6. The timely detection time k is 10 s and the speed of the small vessel w is 20 m/s. The total budget M for both types of detectors is 350 thousand. All parameter values are summarized in Table 6. The optimal solution is x1;176 ¼ x2;78 ¼ x2;227 ¼ x2;319 ¼ 1, that is, one type 1 detector is placed in square 176 and three type 2 detectors are placed in squares 78, 227, and 319. The optimal solution is shown in Fig. 7 where a circle with 1 in the center means a type 1 detector and a circle with 2 in the center means a type 2 detector. The corresponding optimal objective function value is 49.436 million. When we use the greedy adding heuristic, the obtained solution is x1;78 ¼ x1;230 ¼ x1;319 ¼ 1, that is, three type 1 detectors are placed in squares 78, 230, and 319. The corresponding objective function value is 53.106 mil 100% ¼ 7:424%. lion. The relative deviation of the greedy adding heuristic is 53:10649:436 49:436 Based on the parameter values provided in Table 6, we conduct sensitivity analysis with respect to the total budget M, the type 2 detector characteristics ðq2 ; s2 ; g2 Þ, and the probabilities of choosing paths ckj . We vary the total budget from 100 to 800 with step size 100. The corresponding optimal solutions and objective values are shown in Table 7. In the table, the first number within a pair of parentheses represents a type and the second number represents a square number. We observe that when the total budget increases, the optimal objective value decreases. Moreover, in general the marginal benefit of each additional 100 thousand budget is decreasing. We focus on sensitivity analysis on the type 2 detector characteristics. We keep ðq1 ; s1 ; g1 Þ ¼ ð100; 500; 0:006Þ and vary q2 from 140 to 60 with step size 20, s2 from 700 to 300 with step size 100, and g2 from 0.008 to 0.004 with step size 0.001. We remove all dominant combinations where one detector type dominates the other type in terms of the three characteristics. In total, we obtain 72 combinations. The corresponding optimal compositions of both types of detectors are shown in Table 8 where the numbers below Type 1 and Type 2 represent the numbers of Type 1 and Type 2 detectors, respectively. We observe that if we keep the values of s2 and g2 unchanged and reduce the value of q2 , the number of type 2 detectors in the corresponding optimal composition will increase. But if we reduce the value of s2 or if we reduce the value of g2 , the number of type 2 detectors will decrease. Which type is more cost effective depends on detector characteristics. For example, if ðq2 ; s2 ; g2 Þ ¼ ð100; 700; 0:005Þ, we will utilize only type 2 detectors, that is, type 2 detectors are more cost effective. Moreover, if we reduce the value of q2 to 80 and 60, we will stick to only type 2 detectors. If ðq2 ; s2 ; g2 Þ ¼ ð100; 400; 0:007Þ, we will utilize only type 1 detectors, that is, type 1 detectors are more cost effective. Moreover, if we increase the value of q2 to 120 and 140 or if we reduce the value of s2 to 300, we will stick to only type 1 detectors. Now, we turn to sensitivity analysis on the probabilities of choosing paths. We construct 11 situations as shown in Table 9. Situation 1 represents the situation where all probabilities are equal. Situations 2 and 3 represent the situations where the probabilities of choosing path P 4;178 (that is, attacking the target with the lower damage cost) are 0.4 and 0.8, respectively, while situations 4 and 5 represent the situations where the probabilities of choosing path P4;208 (that is, attacking the target with the higher damage cost) are 0.4 and 0.8, respectively. Situations 6 and 7 represent the situations where the terrorist is more likely to choose paths with short distances. The total probabilities of choosing the three shortest paths for situations 6 and 7 are 0.4 and 0.8, respectively. Situations 8 and 9 represent the situations where the terrorist is more likely to choose paths close to the shore, while situations 10 and 11 represent the situations where the terrorist is more likely to choose paths far from the shore. For situations 8 and 9 (10 and 11), the total probabilities for choosing the three paths closest to (farthest

90

X. Yan, X. Nie / Transportation Research Part E 93 (2016) 71–94

from) the shore are 0.4 and 0.8, respectively. From the optimal solutions, we observe that when the probabilities of choosing some specific paths are higher, detectors will be located closer to the corresponding paths. 7. Conclusions and future work In this paper, we investigate how to effectively deploy multiple types of detectors under a small vessel attack threat to port security. We formulate the detector placement problem as a nonlinear binary integer program. Moreover, we propose two branch and bound based exact algorithms (standard and linearized branch and bound algorithms) and a greedy adding heuristic. Based on a detailed computational study, we observe that the linearized branch and bound algorithm is preferred to the standard one and the greedy adding heuristic performs well. Moreover, we provide a case study in New York Harbor. There are several model limitations that lead to future research directions. First, one model limitation is that we assume that detectors work independently, that is, each detector separately declares a passing small vessel as either carrying explosive materials or not. This assumption may not be reasonable if the effective detection areas of the placed detectors have significant intersection. Under these circumstances, a joint detection probability should be incorporated. Second, in current work, we only consider stationary detectors, that is, the locations of detectors are fixed. Recent technology advances have enabled the adoption of mobile detectors for military and civil applications. How to effectively deploy mobile detectors (attached to patrol boats, leisure boats, or even unmanned vehicles) in port areas is another direction. Third, an implicit assumption in our model is that the terrorist is blind with respect to the locations of detectors. This assumption may not be realistic because the terrorist is intelligent. In order to incorporate the terrorist’s behavior, a natural approach is through game-theoretic modeling (for example, a leader–follower game). Last, we assume that the terrorist’s probabilities of choosing paths are known through intelligence analysis. In real situations, these probabilities are uncertain data and a robust optimization modeling technique would be beneficial. Appendix A. Proof of Theorem 1 Before we prove Theorem 1, we introduce some notations. For every ðk; jÞ pair, the portion of path Pkj within the circle centered at square u with radius sh1 has at most two line segments. Moreover, for one line segment and two line segments with their intersection point not on the circle boundary cases, there are either two or four intersections with the circle. If the intersection point of two segments is on the circle boundary, there are three intersections with the circle (including the intersection point of the two segments). There are three possible patterns as shown in Fig. A.1. Pattern 1 indicates one line segment, Pattern 2 indicates two line segments and their intersection point is either outside the circle or on the circle boundary, and Pattern 3 indicates two line segments and the intersection is within the circle. For Patterns 2 and 3, we use cukj1 and cukj2 to represent the two distances from square u to the lines passing through the two segments. For convenience, we label the shorter one with subscript 1, that is, cukj1 6 cukj2 . For Pattern 1, we have that cukj1 ¼ cukj2 . Following these notations, the condition in (ii) of Theorem 1 becomes cukj2 6 cv kj1 . If there are two line segments within the circle, there are two perpendicular lines from the square center to the lines passing through the two segments. If the vertical point is on the line segment, we claim that the distance from the vertical point to the intersection point is positive. Otherwise, if the vertical point is not on the line segment (that is, on the extension of the segment), we claim that the distance is negative. For square i, define the algebraic sum of these two distances as di . Three possible structures for di are shown in Fig. A.2. From the figure, du ¼ du1 þ du2 is the sum of two positive distances, dv ¼ ðdv 1 þ dv 2 Þ is the sum of two negative distances, and dw ¼ dw1  dw2 is the sum of a positive and a negative distances. The proof of Theorem 1 is as follows: Proof. From Definition 2, we have that ph1 ukj P ph1 v kj for all ðk; jÞ pairs. Based on the formula ph1 ikj ¼ 1  egh1 lh1 ikj , we obtain that lh1 ukj P lh1 v kj for all ðk; jÞ pairs. Our target is to show that lhukj P lhv kj for each ðk; jÞ pair and each type h whose radius is smaller than sh1 . We note that if sh 6 cv kj1 ; lhv kj ¼ 0 but lhukj P 0. It is obvious that lhukj P lhv kj . Therefore, we focus on all

c cukuj2kj1

cu

kj1

cukj2

u

u

Pattern 1

Pattern 2

Fig. A.1. Possible patterns and the associated distances.

cuk u j1 j 2 c uk

Pattern 3

91

X. Yan, X. Nie / Transportation Research Part E 93 (2016) 71–94

du

du

1

2

w

v

u

dw

2

dv

1

dw

dv

1

2

Fig. A.2. Possible structures for di .

types h such that sh > cv kj1 . Furthermore, if the condition (i) of Theorem 1 is true, it is obvious that lhukj P lhv kj . Hence, we will concentrate on the situations where there are either one or two line segments within the circles centered at squares u and v with radius sh1 . From Fig. A.1, length lh1 ikj can be calculated by the following formula:

lh1 ikj ¼ where

qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi s2h1  c2ikj1 þ s2h1  c2ikj2 þ 1i ;

1i ¼ 0 for Pattern 1, 1i ¼

8ðk; jÞ and i ¼ u; v ;

qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi s2h1  c2ikj1 þ s2h1  c2ikj2 for Pattern 2, and

ðA:1Þ

1i ¼ di for Pattern 3.

Based on the length of detection radius sh , we consider the following two cases: Case 1. sh P cv kj2 . Case 2. cv kj1 < cv kj2 and cv kj1 < sh < cv kj2 . Case 1. sh P cv kj2 . For this case, based on different patterns, we have the following four subcases: Case 1.1. 1u ¼ 0. qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi Case 1.2. 1u ¼ s2h1  c2ukj1 þ s2h1  c2ukj2 or 1u ¼ du and du P s2h  c2ukj1 þ s2h  c2ukj2 . qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi Case 1.3. 1u ¼ du and jdu j < s2h  c2ukj1 þ s2h  c2ukj2 . qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi s2h  c2ukj1 þ s2h  c2ukj2 . Case 1.4. 1u ¼ du and du 6  Case 1.1. 1u ¼ 0. For this subcase, as in Pattern 1 in Fig. A.1, the portion of Pkj within the circle centered at square u with radius sh has one line segment. Lengths lhukj and lhv kj satisfy the following two expressions:

lhukj ¼

qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi s2h  c2ukj1 þ s2h  c2ukj2 and lhv kj 6 s2h  c2v kj1 þ s2h  c2v kj2 þ 1v :

From (A.1) and the fact that lh1 ukj P lh1 v kj , we obtain that

gðsh1 Þ  1v P 0; where

c2v kj1  c2ukj1 c2v kj2  c2ukj2 gðsÞ ¼ qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi þ qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi s2  c2ukj1 þ s2  c2v kj1 s2  c2ukj2 þ s2  c2v kj2

ðA:2Þ

is a monotonically decreasing function. Hence,

lhukj  lhv kj P gðsh Þ  1v P gðsh1 Þ  1v P 0: Case 1.2.

1u ¼

qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi s2h1  c2ukj1 þ s2h1  c2ukj2 or

1u ¼ du and du P

qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi s2h  c2ukj1 þ s2h  c2ukj2 .

For this subcase, as in Pattern 2 in Fig. A.1, the portion of Pkj within the circle centered at square u with radius sh has two line segments and the intersection point is either outside the circle or on the circle boundary. Lengths lhukj and lhv kj satisfy the following two expressions:

lhukj ¼ 2

qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi s2h  c2ukj1 þ 2 s2h  c2ukj2 and lhv kj 6 2 s2h  c2v kj1 þ 2 s2h  c2v kj2 :

Since cukj1 6 cv kj1 and cukj2 6 cv kj2 , we have that lhukj P lhv kj . qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi Case 1.3. 1u ¼ du and jdu j < s2h  c2ukj1 þ s2h  c2ukj2 . For this subcase, as in Pattern 3 in Fig. A.1, the portion of Pkj within the circle centered at square u with radius sh has two line segments and the intersection point is within the circle. Lengths lhukj and lhv kj satisfy the following two expressions:

92

X. Yan, X. Nie / Transportation Research Part E 93 (2016) 71–94

lhukj ¼

qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi s2h  c2ukj1 þ s2h  c2ukj2 þ 1u and lhv kj 6 s2h  c2v kj1 þ s2h  c2v kj2 þ 1v :

From (A.1) and the fact that lh1 ukj P lh1 v kj , we have that

gðsh1 Þ þ 1u  1v P 0; where gðsÞ is defined by (A.2). Therefore,

lhukj  lhv kj P gðsh Þ þ 1u  1v P gðsh1 Þ þ 1u  1v P 0: Case 1.4.

qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi s2h  c2ukj1 þ s2h  c2ukj2 .

1u ¼ du and du 6 

For this subcase, there is no portion of path Pkj within the circle centered at square u with radius sh . For this subcase, lhukj ¼ 0. At the same time, we must have that lhv kj ¼ 0. Suppose that it is not true, that is, lhv kj > 0. lhukj ¼ 0 means that qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi s2h  c2ukj1 þ s2h  c2ukj2 þ 1u 6 0, while lhv kj > 0 means that s2h  c2v kj1 þ s2h  c2v kj2 þ 1v > 0. It follows that

qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi s2h  c2ukj1 þ s2h  c2ukj2 þ 1u  s2h  c2v kj1  s2h  c2v kj2  1v < 0;

that is,

gðsh Þ þ 1u  1v < 0: Furthermore,

lh1 ukj  lh1 v kj ¼ gðsh1 Þ þ 1u  1v 6 gðsh Þ þ 1u  1v < 0; which contradicts with the fact that lh1 ukj P lh1 v kj . Case 2. cv kj1 < cv kj2 and cv kj1 < sh < cv kj2 . From the proof process for Case 1, there is a type e with se ¼ cv kj2 and leukj P lev kj for all ðk; jÞ pairs. Hence, in this case, our proof can be converted to show that lhukj P lhv kj for all ðk; jÞ pairs and all types h such that cv kj1 < sh < cv kj2 under the condition that leukj P lev kj , which is a special case of Case 1. Combining Case 1 and Case 2, we obtain that lhukj P lhv kj , which implies that phukj P phv kj for all ðk; jÞ pairs. Hence, square u strongly dominates square v. Appendix B. Proof of Theorem 2 Before we prove Theorem 2, we provide the following lemma. Lemma 1. Let N kj ¼ [h N hkj be the potential detector location set for path Pkj . If square u strongly dominates square v, we have that (i) fh : u 2 T h g  fh : v 2 T h g; (ii) u R N kj implies that v R N kj .

Proof. (i) (By contradiction) If fh : u 2 T h g + fh : v 2 T h g, there is a type b with

v 2 T b and u

R T b . From the construction of T b ,

there is a path Pk0 j0 such that pbv k0 j0 > 0 and pbuk0 j0 ¼ 0. It contradicts with the fact that square u strongly dominates square v. (ii) (By contradiction) If v 2 N kj , by the construction of N kj , we have that pbv kj > 0 for some type b. On the other hand, u R N kj implies that phukj ¼ 0 for every type h. It contradicts with the fact that square u strongly dominates square v. The proof of Theorem 2 is as follows: P P H Proof. (By contradiction) Suppose that xH < h:v 2T h xH hv in an optimal solution X . Since both sides of the strict h:u2T h hu P P H H inequality are binary, we have that h:u2T h xhu ¼ 0 and h:v 2T h xhv ¼ 1. Moreover, because all decision variables are binary, we h h H H have that xH hu ¼ 0 for h 2 fh : u 2 T g and there exists a type b 2 fh : v 2 T g such that xbv ¼ 1 and xhv ¼ 0 for h h 2 fh : v 2 T g n fbg. From the conclusion (i) of Lemma 1, we have that b 2 fh : u 2 T h g. Hence, we can find another feasible solution X HH with H HH xHH bu ¼ 1; xbv ¼ 0, and other values being the same as those of X . From the conclusion (ii) of Lemma 1, we observe that there is no path P kj such that u R N kj but v 2 N kj . Hence, the objective function can be rewritten in the following form:

93

X. Yan, X. Nie / Transportation Research Part E 93 (2016) 71–94

X

f ðXÞ ¼

W kj

kj:uRNkj

W kj

kj:u2N kj

YY

v 2Nkj

Y

ð1  phikj Þxhi

YY

ð1  phikj Þxhi

h2Hi2Nh

kj i–u;v

Y h:u2N hkj

ð1  phukj Þxhu

Y h:v

ð1  phv kj Þxhv

2N hkj

ð1  phukj Þxhu :

h:u2Nhkj

h2Hi2Nh

kj i–u

v RNkj

W kj

kj:u2Nkj

kj

X

X

ð1  phikj Þxhi þ

h2Hi2Nh

v RNkj

þ

YY

For convenience, we let

I¼

X

W kj

kj:uRNkj

ð1  phikj Þxhi ;

h2Hi2Nh kj

v RNkj

J kj ¼ W kj

YY

YY

ð1  phikj Þxhi ;

h2Hi2Nh

kj i–u;v

and

K kj ¼ W kj

YY

ð1  phikj Þxhi :

h2Hi2Nh

kj i–u

Then, the objective function values of X H and X HH are as follows:

f ðX H Þ ¼ I þ

X

kj:u2Nkj

f ðX HH Þ ¼ I þ

X

K kj

ðB:1Þ

kj:u2Nkj

v 2Nkj

and

X

J kj ð1  pbv kj Þ þ

v RNkj

J kj ð1  pbukj Þ þ

kj:u2N kj

v 2Nkj

X

K kj ð1  pbukj Þ:

ðB:2Þ

kj:u2N kj

v RNkj

Based on (B.1) and (B.2), we observe that the objective function value of X HH is smaller than that of X H . It contradicts with the fact that X H is an optimal solution.

Appendix C. Proof of Theorem 3 Proof. The proof procedure is similar to that of Theorem 17 in Sherali et al. (2004) except that the X-part decision variables are integer in our algorithm while they are continuous in Sherali et al. (2004). Therefore, we need to prove that there is a convergent subsequence for the solution sequence generated by the linearized branch and bound algorithm in infinite iterations. For each node nðsÞ, let CnðsÞ ¼ ðXnðsÞ ; YnðsÞ ; ZnðsÞ Þ be an optimal solution of [DPðLnðsÞ ; UnðsÞ Þ-R]. Because the X-variable parts are finite in any infinite branch, there is a set N 1 # N such that

lim XnðsÞ ¼ XH : s!1

ðC:1Þ

s2N1

Furthermore, because the feasible regions for fYnðsÞ ; ZnðsÞ ; LnðsÞ ; UnðsÞ gs2N1 are compact sets, there is a set N 2 # N 1 such that

lim ðYnðsÞ ; ZnðsÞ ; LnðsÞ ; UnðsÞ Þ ¼ ðYH ; ZH ; LH ; UH Þ: s!1

ðC:2Þ

s2N2

Based on (C.1) and (C.2), we obtain a convergent subsequence. The remaining process is identical to that in Sherali et al. (2004). References Atkinson, M.P., Cao, Z., Wein, L.M., 2008. Optimal stopping analysis of a radiation detection system to protect cities from a nuclear terrorist attack. Risk Anal. 28 (2), 353–371. Atkinson, M.P., Wein, L.M., 2008. Spatial queueing analysis of an interdiction system to protect cities from a nuclear terrorist attack. Oper. Res. 56 (1), 247– 254. Boros, E., Fedzhora, L., Kantor, P.B., Saeger, K., Stroud, P., 2009. A large-scale linear programming model for finding optimal container inspection strategies. Nav. Res. Logist. 56 (5), 404–420. Brown, G., Carlyle, M., Abdul-Ghaffar, A., Kline, J., 2011. A defender-attacker optimization of port radar surveillance. Nav. Res. Logist. 58 (3), 223–235. CBS News, 2008. U.S. Wary of Small Boat Terrorism. Washington, DC, 28 April 2008. Chakrabarty, K., Iyengar, S.S., Qi, H., Cho, E., 2002. Grid coverage for surveillance and target location in distributed sensor networks. IEEE Trans. Comput. 51 (12), 1448–1453.

94

X. Yan, X. Nie / Transportation Research Part E 93 (2016) 71–94

Concho, A.L., Ramirez-Marquez, J.E., 2010. An evolutionary algorithm for port-of-entry security optimization considering sensor thresholds. Reliab. Eng. Syst. Safety 95 (3), 255–266. Danczyk, A., Liu, H.X., 2011. A mixed-integer linear program for optimizing sensor locations along freeway corridors. Transp. Res. Part B: Methodol. 45 (1), 208–217. Dijkstra, E.W., 1959. A note on two problems in connexion with graphs. Numer. Math. 1 (1), 269–271. Elsayed, E.A., Young, C.M., Xie, M., Zhang, H., Zhu, Y., 2009. Port-of-entry inspection: sensor deployment policy optimization. IEEE Trans. Autom. Sci. Eng. 6 (2), 265–276. Frank, T., 2007. Small Boats Seen as a Terror Threat. USA Today, 31 October 2007. Gaudette, K., 2002. Court Approves Bush Request to Open West Coast Ports. The Seattle Times, 8 October 2002. Gaukler, G.M., Li, C., Ding, Y., Chirayath, S.S., 2012. Detecting nuclear materials smuggling: performance evaluation of container inspection policies. Risk Anal. 32 (3), 531–554. Hill, B.P., 2009. Maritime Terrorism and the Small Boat Attack Threat to the United States: A Proposed Response. Master Thesis. Naval Postgraduate School, Monterey, CA. Hochbaum, D.S., Fishbain, B., 2011. Nuclear threat detection with mobile distributed sensor networks. Ann. Oper. Res. 187 (1), 45–63. Kaplan, E.H., Kress, M., 2005. Operational effectiveness of suicide-bomber-detector schemes: a best-case analysis. Proc. Natl. Acad. Sci. USA 102 (29), 10399– 10404. Lee, A.J., Jacobson, S.H., 2011. Evaluating the effectiveness of sequential aviation security screening policies. IIE Trans. 43 (8), 547–565. Lee, A.J., Jacobson, S.H., 2012. Addressing passenger risk uncertainty for aviation security screening. Transp. Sci. 46 (2), 189–203. Li, X., Ouyang, Y., 2011. Reliable sensor deployment for network traffic surveillance. Transp. Res. Part B: Methodol. 45 (1), 218–231. McLay, L.A., Dreiding, R., 2012. Multilevel, threshold-based policies for cargo container security screening systems. Eur. J. Oper. Res. 220 (2), 522–529. McLay, L.A., Lee, A.J., Jacobson, S.H., 2010. Risk-based policies for airport security checkpoint screening. Transp. Sci. 44 (3), 333–349. McLay, L.A., Lloyd, J.D., Niman, E., 2011. Interdicting nuclear material on cargo containers using knapsack problem models. Ann. Oper. Res. 187 (1), 185–205. Merrick, J.R.W., McLay, L.A., 2010. Is screening cargo containers for smuggled nuclear threats worthwhile? Decis. Anal. 7 (2), 155–171. Morton, D.P., Pan, F., Saeger, K.J., 2007. Models for nuclear smuggling interdiction. IIE Trans. 39 (1), 3–14. Neild, B., 2011. The Small Vessel Threat. GlobalPost, 9 February 2011. Nie, X., 2011. Risk-based grouping for checked baggage screening systems. Reliab. Eng. Syst. Safety 96 (11), 1499–1506. Nie, X., Batta, R., Drury, C.G., Lin, L., 2007. Optimal placement of suicide bomber detectors. Mil. Oper. Res. 12 (2), 65–78. Nie, X., Batta, R., Drury, C.G., Lin, L., 2009a. The impact of joint responses of devices in an airport security system. Risk Anal. 29 (2), 298–311. Nie, X., Batta, R., Drury, C.G., Lin, L., 2009b. Passenger grouping with risk levels in an airport security system. Eur. J. Oper. Res. 194 (2), 574–584. Nie, X., Parab, G., Batta, R., Lin, L., 2012. Simulation-based selectee lane queueing design for passenger checkpoint screening. Eur. J. Oper. Res. 219 (1), 146– 155. Nikolaev, A.G., Lee, A.J., Jacobson, S.H., 2012. Optimal aviation security screening strategies with dynamic passenger risk updates. IEEE Trans. Intell. Transp. Syst. 13 (1), 203–212. Przemieniecki, J.S., 2000. Mathematical Methods in Defense Analyses. American Institute of Aeronautics and Astronautics, Reston, VA. Ramirez-Marquez, J.E., 2008. Port-of-entry safety via the reliability optimization of container inspection strategy through an evolutionary approach. Reliab. Eng. Syst. Safety 93 (11), 1698–1709. Sherali, H.D., Desai, J., Glickman, T.S., 2004. Allocating emergency response resources to minimize risk with equity considerations. Am. J. Math. Manage. Sci. 24 (3–4), 367–410. Sullivan, K.M., Morton, D.P., Pan, F., Smith, J.C., 2014. Securing a border under asymmetric information. Nav. Res. Logist. 61 (2), 91–100. Szechtman, R., Kress, M., Lin, K., Cfir, D., 2008. Models of sensor operations for border surveillance. Nav. Res. Logist. 55 (1), 27–41. Trott, B., 2015. West Coast Port Dispute Costing Honda 25,000 Vehicles. Reuters, 21 February 2015. United Nations Conference on Trade and Development, 2011. Review of Maritime Transport 2011. New York and Geneva, November 2011. United States Department of Homeland Security, 2007. DHS Announces West Coast Maritime Radiation Detection Project. Washington, DC, September 2007. United States Department of Homeland Security, 2008. Small Vessel Security Strategy. Washington, DC, April 2008. United States Government Accounting Office, 2002. Nation Faces Formidable Challenges in Making New Initiatives Successful. Washington, DC, August 2002. Van Weele, S.F., Ramirez-Marquez, J.E., 2011. Optimization of container inspection strategy via a genetic algorithm. Ann. Oper. Res. 187 (1), 229–247. Wein, L.M., Atkinson, M.P., 2007. The last line of defense: designing radiation detection-interdiction systems to protect cities from a nuclear terrorist attack. IEEE Trans. Nucl. Sci. 54 (3), 654–669. Wein, L.M., Liu, Y., Motskin, A., 2009. Analyzing the homeland security of the U.S.-Mexico border. Risk Anal. 29 (5), 699–713. Wein, L.M., Wilkins, A.H., Baveja, M., Flynn, S.E., 2006. Preventing the importation of illicit nuclear materials in shipping containers. Risk Anal. 26 (5), 1377– 1393. Wilhelm, W.E., Gokce, E.I., 2010. Branch-and-price decomposition to design a surveillance system for port and waterway security. IEEE Trans. Autom. Sci. Eng. 7 (2), 316–325. Young, C.M., Li, M., Zhu, Y., Xie, M., Elsayed, E.A., Asamov, T., 2010. Multiobjective optimization of a port-of-entry inspection policy. IEEE Trans. Autom. Sci. Eng. 7 (2), 392–400. Zhang, Z.G., Luh, H.P., Wang, C.H., 2011. Modeling security-check queues. Manage. Sci. 57 (11), 1979–1995. Zhu, N., Liu, Y., Ma, S., He, Z., 2014. Mobile traffic sensor routing in dynamic transportation systems. IEEE Trans. Intell. Transp. Syst. 15 (5), 2273–2285. Zou, Y., Chakrabarty, K., 2004. Sensor deployment and target localization in distributed sensor networks. ACM Trans. Embed. Comput. Syst. 3 (1), 61–91.