- Email: [email protected]
Planning safe navigation routes through mined waters
Accepted Manuscript Planning safe navigation routes through mined waters Luitpold Babel, Thomas Zimmermann PII: DOI: Reference:
S0377-2217(14)00618-3 http://dx.doi.org/10.1016/j.ejor.2014.08.002 EOR 12444
To appear in:
European Journal of Operational Research
Received Date: Accepted Date:
17 March 2014 2 August 2014
Please cite this article as: Babel, L., Zimmermann, T., Planning safe navigation routes through mined waters, European Journal of Operational Research (2014), doi: http://dx.doi.org/10.1016/j.ejor.2014.08.002
This is a PDF file of an unedited manuscript that has been accepted for publication. As a service to our customers we are providing this early version of the manuscript. The manuscript will undergo copyediting, typesetting, and review of the resulting proof before it is published in its final form. Please note that during the production process errors may be discovered which could affect the content, and all legal disclaimers that apply to the journal pertain.
Planning safe navigation routes through mined waters Luitpold Babel∗ , Thomas Zimmermann Institut f¨ ur Mathematik und Informatik, Fakult¨ at Betriebswirtschaft, Universit¨ at der Bundeswehr M¨ unchen, D-85579 Neubiberg, Germany
Abstract In this paper we investigate the problem of finding a safe transit of a ship through areas threatened by sea mines. The aim is to provide decision-making support by a tool that can be integrated into a naval command and control system. We present a route finding algorithm which avoids regions of risk higher than a given threshold. The algorithm takes into account the technical and operational restrictions of the ship’s movement. It allows to minimize the route length, the traveling time, the number of maneuvers, or other objectives. The basic idea is to embed a network in the operational area and compute a least-cost path. Instead of using a regular grid graph which strongly restricts the types of maneuvers and necessitates a path smoothing after optimization, we design a network which is especially tailored to the maneuverability of the vessel. Each path in this network represents a continuous-curvature track based on combinations of clothoids and straight line segments. The approach allows a large variety of maneuvers, hence high-quality solutions are achievable provided a sufficiently dense network. Keywords: Graph theory, Path planning, Navigation routes, Sea mines, Clothoid curves
∗ Corresponding
author. Tel.: +49 89 6004 3267. Email address: [email protected]
Preprint submitted to European Journal of Operational Research
August 6, 2014
1. Introduction The rapid increase in the proliferation of sea mines presents a unique challenge to maritime security. Sea mines are self-contained explosive devices placed in water to destroy surface ships or submarines. Since mines are inexpensive 5
they form an effective means to deny the use of large areas. For over a century, they have been used as a strategic weapon in military conflicts to blockade ports, waterways or other maritime zones. The naval forces must be prepared to counter the danger posed by mines. This, in particular, includes to find safe passages in order to sustain commercial and military shipping operations.
10
Planning a navigation route through an area threatened by sea mines is closely related to the classical problem of finding a collision-free path of a vehicle in an environment with obstacles. The latter, also known as path planning or trajectory planning problem, has long been the subject of intensive research. We refer to the reviews of Latombe [13] and LaValle [14].
15
Numerous ideas have been investigated to resolve the path planning problem. The most promising include the potential field method, cell decomposition, the roadmap method, the mass-spring-damper method, mixed integer linear programs, and several network-based methods with discretizations of the space by visibility graphs, Voronoi diagrams or regular grids. Newer approaches include
20
the rapidly-exploring random tree method, model predictive control methods, and mathematical programming methods. Comprehensive reviews are given by Goerzen et al. [10] and Ferguson et al. [9]. While previous research mainly focused on vehicles moving on ground (mobile robots, autonomous vehicles) or in air (unmanned aerial vehicles, drones,
25
autonomous helicopters), less work has been published for sea vehicles. Although the problems are closely related the proposed algorithms are not directly applicable since, from a technological point of view, ships behave quite differently than ground vehicles or aircraft. Moreover, shipping may imply specific operational requirements.
30
An important fact is that continuous methods can only find locally optimal
2
solutions which may be far from globally optimal. On the other hand, networkbased methods give an approximate problem, but means that one can find a globally optimal solution (see e.g. Zabarankin et al. [23]). Furthermore, networkbased methods appear more promising for realizing specific technological and 35
operational constraints. Fagerholt et al. [8] present an algorithm to find a shortest route for a ship located at sea to a destination port, in the presence of polygonal obstacles. It relies on the partial computation of a visibility graph and solves the problem using two variants of the algorithm of Dijkstra.
40
In a series of papers, Piatko et al. [18, 19, 21] present path planning techniques to find low-risk paths of a ship through a minefield. They discretize the operational area by evenly distributing vertices, with eight edges emerging from each vertex (along the axes and the diagonals), hence generating a regular grid graph. The authors develop probabilistic risk models and identify different path
45
quality attributes. The problem is then solved by defining appropriate costs on the edges of the graph and calculating a least-cost path. A similar method has been suggested by Bekker and Schmid [5]. They also use a regular grid graph, and remove all edges which lie within some safety radius around the mines. The remaining edges get as a cost either the length,
50
or a risk value depending on the distance to the mines, or a mixture of both. Again a least-cost path has to be determined. Their approach has been extended by Li [15] who added longer edges to the graph in order to improve the solution. A drawback of the above methods is that important factors including the characteristics of the ship (maneuverability, turning circle, etc.) or restrictions
55
which are implied by operational reasons (e.g. no sharp maneuvers due to passenger or cargo safety) are not considered. Each method provides a polygonal line, possibly with sharp turns or extreme maneuvers, which the ship may not be able to follow. Turn constraints have been successfully implemented in a recent research of
60
Ari et al. [4] where the goal is to find a shortest path through an area with obstacles. Their approach is based on a regular grid as above, however, they 3
use a vertex replication technique where each vertex is split into copies labeled by the direction the ship is coming from. This incorporates navigation history into the present location and allows to impose realistic ship turn constraints. 65
Moreover, the operational area can be discretized at any desired resolution. Indeed, there is still room for improvement since ship maneuvers are limited to 45◦ turns with full rudder deflection (which might imply that the solution is far from optimal) and a path smoothing is required after optimization (which might cause collision with an obstacle if not provided with a buffer zone).
70
The aim of our work is to develop an algorithm which finds a safe route through a mined area. Particular attention will be given to the technical capabilities and operational restrictions of the ship’s movement. Moreover, the maritime environment including shallows, islands, obstacles, etc. will be considered. The algorithm is intended to be integrated as a decision-making tool into
75
a naval command and control system. We propose a network optimization approach, focusing on the generation of a special randomized network which is adapted to the maneuverability of the ship. The paths in the network correspond to smooth (continuous-curvature) tracks consisting of clothoids and straight line segments. The approach allows to set
80
restrictions on the track including (non-symmetric) turn radius constraints, thus avoiding too sharp maneuvers, but also other restrictions such as the exclusion of long-lasting weak maneuvers (which might be undesirable due to operational reasons). The change of velocity during turning maneuvers is taken into account. The random structure of the network allows a large variety of maneuvers of
85
different strengths and without turn degree limitations. By applying standard shortest-path methods, we can determine routes which avoid regions of risk higher than a given threshold. With appropriate costs defined on the edges, different objectives can be considered such as minimizing route length, traveling time, the number of maneuvers, or some compromise of
90
these. Like all network-based methods, our approach is of a heuristic nature. On the other hand, the quality of the solutions can be continually improved by increasing the size of the network, as illustrated in Section 5. 4
In the next section, a mathematical model of the ship’s movement is presented, along with a brief description of the threat model and the optimization 95
problem. Section 3 contains the basic idea of the network generation and the route planning algorithm. In Section 4 we discuss how to speed up the algorithm by using efficient intersection procedures, and by separating the planning into preprocessing and online planning. Computational results are presented in Section 5. We conclude with final remarks in Section 6.
100
2. Mathematical modeling 2.1. Ship track model Our first task is to establish a realistic model of the movement of a ship. From an operational point of view, a planned ship track should fulfill certain requirements. First, the track should be easy to navigate, i.e. it should allow easy
105
handling and control of the ship by the crew. To the extent possible, long-lasting maneuvers and extreme course alterations should be avoided. Finally, the track’s curvature should be continuous in order to guarantee a smooth navigation. This is because following a planned track should not unduly burden the ship and the crew. These aspects become crucial in transporting people or dangerous goods.
110
Comfort and safety can only be ensured for continuous-curvature trajectories. A clothoid, also known as Euler or Cornu spiral, is a curve whose curvature changes linearly with its length [1]. It is well known that a vehicle following the curve at constant speed will have a constant change of centripetal acceleration. Hence, the clothoid is a perfect representation of a vehicle traveling at constant
115
speed with a constant rate of change of the steering angle [12]. The linear steering motion makes such a track comfortable for vehicle and crew. It has been shown in [17] that a clothoid very well approximates the curve of a turning ship up to a turning angle of about 100◦ to 120◦ , although the speed is not constant. In the first phase of the maneuver, when the rudder
120
is deflected, the track is close to a straight line since the ship starts turning very slowly. Then the curvature increases and the track gradually turns into a
5
Figure 1: Phases of ship turning, turning circle and clothoid
circle. During this part of the turn, the speed of the ship decreases. If the angle exceeds 120◦ the ship moves with constant speed in a circle with constant radius. The characteristic phases of a turning ship and the approximating clothoid are 125
sketched in Figure 1. Our goal is to realize the ship track as a combination of straight line segments and clothoids. The transitions between the track segments have to be smooth. A clothoid can be expressed as a parametrized curve ⎞ ⎛ x(t) ⎠ , 0≤t≤T γ(t) = ⎝ y(t) with domain [ 0, T ] and t √ πτ 2 x(t) = a π dτ , cos 2 0
130
√ y(t) = a π
t
sin 0
πτ 2 dτ . 2
It is well known and straightforward to verify that the arc length and the curvature of the curve are √ s = a πt ,
√ κ =
π t a
respectively. Hence κ =
1 s a2
i.e. the curvature is proportional to the arc length. The parameter a represents a scale factor which defines the rate of change of the curvature and hence the 6
135
size of the clothoid. The Fresnel integrals t t πτ 2 πτ 2 dτ , S(t) = dτ C(t) = cos sin 2 2 0 0 do not have an analytical solution. According to the power series expansions of cosine and sine cos(t) = 1 −
t2 t4 t6 + − + ... 2! 4! 6!
sin(t) = t −
t5 t7 t3 + − + ... 3! 5! 7!
they can be solved numerically using 1 π 4 t9 1 π 6 t13 1 π 2 t5 + − + ... 2! 2 5 4! 2 9 6! 2 13 1 π 3 t7 1 π 5 t11 1 π 7 t15 π t3 − + − + ... S(t) = 2 3 3! 2 7 5! 2 11 7! 2 15
C(t) = t − 140
For more details, particularly with regard to convergence, numerical stability and efficient computation, we refer to the standard text book [20]. We define a configuration as a combination of a position and a direction in the plane. For the route planning algorithm that will be presented in Section 3, 145
we need to solve a huge number of subproblems (SP) of the following type: Given two configurations, find a ship track from the first (initial) (SP)
to the second (final) configuration. The velocity at both endpoints should equal the cruising speed v0 .
To construct such a track from a configuration A to a configuration B, we first connect A and B by two line segments and a circle segment. We confine ourselves to circle segments with angle smaller than 180◦ . Otherwise 150
(SP) is assumed to have no solution. Hence the ray originating from A intersects with the ray ending in B. A circle with some predefined radius r is placed tangential to the line segments between A and B (see Figure 2). Note that the transition from a straight line path to a tangential circular path implies a sudden centripetal acceleration (which a ship is not able to follow). This is
155
avoided by replacing the circle segment by two symmetric clothoid curves. They
7
Figure 2: Connecting two configurations with two lines and a circle segment
are adjusted in such a way that the acceleration is continuous along the whole track and equal to zero at both ends. Consider a circle segment with radius r and angle 2α as depicted in Figure 3. The first clothoid starts in the origin of the coordinate system tangential to 160
the x-axis and ends in point S, perpendicular to the angle bisector. The second clothoid is mirrored at the bisector line. We denote the first clothoid by γ(t). The parameter a and the range T of γ(t) is then calculated as follows. The tangent vector of γ(t) is
⎞
⎞ ⎛ πt2 cos √ 2 ⎠ ⎠ = a π⎝ γ (t) = ⎝ . 2 y (t) sin πt2 ⎛
x (t)
Since the slope of γ(t) in S equals α we get ⎛ ⎞ ⎞ ⎛ 2 cos πT2 cos α ⎝ ⎠ . ⎠ = ⎝ 2 sin α sin πT2 165
This implies
T =
2α . π
The point S lies on the bisector line y = − cot α · x + r. Hence cot α · x(T ) + y(T ) = r and, by substituting x(T ) and y(T ),
T T √ πτ 2 πτ 2 a π cot α dτ + dτ = r. cos sin 2 2 0 0 8
Figure 3: Replacing a circle segment by two symmetric clothoids
The expression in brackets, abbreviated with I, is evaluated numerically. Then r a = √ πI which determines the parameter of the clothoid. 170
2.2. Ship characteristics The curvature of the curve is maximal at connection point S of the clothoids. At that point
√ κmax =
π 1√ T = 2π α I . a r
Since the ship must be able to follow the track, κmax must not exceed the curvature of the minimum turning circle (which is related to the maximal rudder 175
deflection). If Rturn denotes the minimum turning radius then 1√ 1 2π α I ≤ r Rturn must hold. We define Rmin =
√
2π α I Rturn .
For the construction from above, we restrict the radius r of the circle segment by Rmin ≤ r ≤ Rmax . The upper bound Rmax is the radius of the largest circle that can be placed 180
tangential to the line segments. Note that this bound might be used to forbid 9
long-lasting weak maneuvers, if those are not desired due to operational reasons. The radius r is selected randomly within the bounds in order to enable a large variety of maneuvers. Alternatively, to control the type of the occuring maneuvers, one might also define the radius deterministically, e.g. by always selecting 185
the smallest possible value which provides short but strong maneuvers, or the largest possible value which provides short length curves. If Rmin > Rmax then the subproblem (SP) is assumed to have no solution. One might also define separately restrictions for left and right turns if left (port) and right (starboard) turn radii are different (which indeed is true for
190
most ships, see [2]). In that case, the lower bound Rmin has to be replaced by l r two different values Rmin and Rmin , according to the turn direction.
It remains to create a speed change model. During turn, the drag of the rudder and the sideways drift of the ship result in a progressive loss of speed, even though the engine revolutions are maintained constant (see e.g. [2]). Ex195
perience shows that most ships with full rudder deflection lose about one-third of their original speed after a 90-degree turn. We assume that the ship travels with cruising speed v0 along the straight line segments. The ship decelerates during the first clothoid where the curvature increases, achieves its minimum speed v1 at the end of the first clothoid
200
(at symmetry point S), and accelerates during the second clothoid where the curvature decreases, until it again achieves cruising speed v0 . We suppose a velocity curve with average speed (v0 + v1 )/2 and define v1 as a function of the turning radius r and turning angle α by v1 = v(α, r) = v0 (1 −
1 α Rmin ) 3 π/2 r
where α is restricted by 90◦ . The largest decrease in speed occurs for maximal 205
angle and minimal radius where π 2 v1 = v( , Rmin ) = v0 . 2 3 Since the rate of deceleration varies between different types of ship, one might have to adapt the model accordingly. 10
2.3. Threat model The operational area is a marine area which is passable by ships. The area 210
possibly includes islands, reefs, shallows and nonnatural impassable sites such as shipwrecks and other obstacles. The operational area depends also on the draught of the ship, the tides and other environmental conditions. Modeling is done by defining appropriate (non-convex) polygons. A large outer polygon defines the shoreline or the boundary of the mission area, cut-out inner polygons
215
represent the impassable regions. The threat areas consist of solitary mines and minefields. With the former we mean single mines for which accurate information is available. The latter are aggregations of mines with imprecise or vague information. For a classification of mine types, mine counter measures and methods of mine delivery, see e.g.
220
[16] and [22]. The detection and mapping of sea mines is discussed e.g. in [3] and [24]. Solitary mines are characterized by their location and damage radius. The damage radius is mainly determined by the explosive charge of the mine. However, it also depends on technical ship data including signature and vulnerability.
225
The signature of a ship consists of those emitted signals (e.g. acoustic, magnetic, visual, infrared) which can be used by sensors to detect the ship’s position. The vulnerability of a ship is based on its construction, hull material, armour, etc. These data determine the effective radius of a mine which, under no circumstances, may be entered. Otherwise the ship will sustain damage. The radius
230
can be increased by additional safety margins to take into account uncertainty of the position of the mine and the navigational error of the ship. In contrast to the mapped mines described above, there is only inaccurate information for minefields. The individual positions (and possibly the types) of the mines are not precisely known. Minefields are modeled by polygons which
235
include the supposed positions of the mines with an appropriate safety margin. They can be considered as forbidden areas where the ship is not allowed to pass through. Risk models for mines have been discussed in several papers, e.g. [5, 19, 21]. 11
They may rely on information from different sensors, available actuation and 240
damage curve data and other sources. Mine detection, classification uncertainty, and damage are then represented by suitable probability distributions. Since a full treatment of all relevant data would lead to highly complex models, some simplifications seem inevitable. The point of closest approach to the closest mine is used in [5] to determine
245
the risk of a path. With D denoting the minimal distance, the risk is defined proportional to any user selectable inverse power of D, such as 1/D or 1/D2 . A similar point risk model is used in [19]. The probability of being damaged at each location is an inverse function of the distance to the location of the mine. The risk associated with each path is then estimated by summing the
250
point risks along the path. A more sophisticated model is introduced in [21] where the location of a mine is represented by a normal distribution. This is combined with a triggering model for the mine which is given by a Gaussian-like distribution. 2.4. Optimization problem
255
We call a navigation route passable if it proceeds completely within the operational area. Such a route is clear of obstacles and can be traveled without the danger of collision. A route is called feasible if it considers the technical and operational limitations of the ship, i.e. the ship is able to follow the route. Safety will be expressed by the level of risk that the mission planner or head
260
of operations is willing to bear. It is reflected by a safety margin which is added to the effective radius of each mine. For that purpose, an appropriate risk function f (D) has to be defined which is inversely proportional with distance D to the mine. If pmax denotes the maximum acceptable risk then the solution of the equation f (D) = pmax provides the radius around the mine which can
265
be considered as a forbidden area for the ship. We choose the risk function proportional to 1/D, noting that it can easily be replaced by a more sophisticated model (see [5, 21]). A similar idea is realized for minefields. The risk for the ship reduces if the distance to the minefield increases. The associated 12
polygon is enlarged by a buffer zone (see [4]) until the maximum acceptable risk 270
is achieved. The ship routing problem discussed in this paper is to determine a passable and feasible route of minimum length or minimum traveling time from an initial configuration (starting point with direction) to a final configuration (destination point with direction) which avoids regions with risk higher than a given
275
threshold. The above formulation aims to minimize one criterion (length or time) under a hard bound on the other criterion (the risk). An alternative approach would be to determine a route which minimizes the risk, or a combination of risk and length (time). Both problem definitions have been extensively discussed e.g. in
280
[18, 19]. Let us mention that our algorithm, although designed for the constraintbased approach, can easily be extended to find a passable and feasible route with smallest accumulated risk or combination of risk and length (time). For the calculation of the accumulated risk along a path we refer to [19].
285
3. Route planning The basic idea of our approach is to generate a network over the water surface in such a way that each path in this network represents a possible navigation route. The shortest navigation route in the network is then determined by applying standard methods of combinatorial optimization. The quality of the
290
solution increases with the density of the network. Network-based solutions are also proposed in other publications, e.g. in [5, 15, 18, 19, 21]. The common procedure is to evenly distribute vertices over an area and link close pairs of vertices by straight lines. This provides a discretization of the area in the form of a regular grid network. A drawback is that the
295
ship’s characteristic (turning radius, maneuverability, etc.) is not sufficiently considered. The resulting navigation route is a polygonal line which might be difficult, or even impossible to navigate. Moreover, if the number of vertices in
13
the network is small then the route might be far away from an optimal solution. With increasing number of vertices the route might become an annoying zigzag 300
course. Our network considers the technical and operational limitations of the ship. The network is not regular since its generation is affected by some random process. The vertices are not connected by straight lines but by curved lines consisting of two straight line segments and two symmetric clothoids in between
305
those segments. This provides a smooth track which allows convenient control of the ship. In contrast to the situation described above the track remains smooth, independent of the density of the network. The following algorithm scheme describes the generation of the network and the route planning.
310
Algorithm Find-Ship-Track Start configuration P , target configuration Q
Input:
Operational area, performance data of ship, threat areas Output:
315
Ship track γ(t) from P to Q
1
Create vertices vP and vQ for start configuration P and target configuration Q
2
Initialize vertex set V := {vP , vQ } and edge set E := ∅
3
Generate a huge number of configurations in the operational area by randomly choosing position and direction
4
For each configuration A Create a new vertex vA , V := V ∪ {vA }
5 320
6
end
7
For each pair of configurations (A, B)
8
Generate a feasible track γAB (t) from A to B consisting of clothoids and straight line segments (if existing) by solving the corresponding (SP)
9
Check whether γAB (t) is passable, i.e. the track is entirely within the operational area
325
10
Check whether γAB (t) is safe, i.e. the track avoids the threat areas
11
If γAB (t) is passable and safe then
12
Create a new edge (vA , vB ) with cost c(vA , vB ) :=
13
E := E ∪ {(vA , vB )}
14
|γAB (t)| dt,
Figure 4: Schematic representation of network
330
14
end
15 end 16 Find a shortest path in the network G = (V, E, c) from vP to vQ 17 Concatenate the tracks belonging to the path to a P - Q track γ(t)
The vertices of the network represent randomly generated configurations. Each 335
(directed) edge represents a track segment connecting two configurations. Each path is a concatenation of track segments. The construction of the track segments, described in the previous section as the solutions of subproblems (SP), guarantees that the transition between any two of them is smooth. More precisely, this means that each path corresponds to a twice continuously differen-
340
tiable curve. As a result of the checks performed in lines (9) and (10), each path in the network represents a passable, feasible and safe navigation route. The cost of an edge is the length of the associated track segment. It is √ obtained by adding twice the arc length s = a π T of a clothoid to the lengths of the straight line segments. Alternatively, if a route with smallest traveling
345
time has to be found, the cost of an edge is the time required to pass the track. This is twice the time 2s/(v0 + v1 ) to pass one clothoid plus the time for the straight line segments. The idea of the network generation is sketched in Figure 4. It shows a number of configurations (arrows) with connecting track segments. The scenario 15
350
contains solitary mines (circles), minefields (polygons) and other obstacles. For operational use, the algorithm Find-Ship-Track should deliver a solution within a couple of seconds. Unfortunately, a straightforward implementation is too slow. Hence further efforts are needed to speed up the algorithm.
4. Efficient implementation 355
4.1. Intersection procedures The main workload of the algorithm occurs in lines (8) – (10) where a huge number of tracks is generated and analyzed. The generation follows the description given in Sections 2.1 and 2.2. According to the construction, the tracks are feasible. Checking whether tracks are passable and safe is time-critical.
360
The former means to find out whether a track lies completely within the outer (polygonal) boundary of the operational area and does not collide with an inner (polygonal) obstacle, the latter whether a track passes through a (circular) solitary mine or through a (polygonal) minefield. Hence, we have to develop efficient procedures to check the intersection of a
365
track with circles and polygons. The first problem is addressed by solving the following more general question: Given an arbitrary curve γ(t) with varying speed (describing the motion of a vehicle) and a set of circles with varying radii, does the curve intersect with at least one circle ?
370
We propose the following procedure. Procedure Intersect-Curve-with-Circles Input:
Curve γ(t) with 0 ≤ t ≤ T , minimal time step Δt Maximal velocity or upper bound vmax Set of circles Ck with center Mk and radius rk
375
Output: 1 2
Message, if γ(t) intersects with a circle
Let A := γ(0) and B := γ(T ) be the start and end point of the curve, respectively Calculate the length of the curve l := 0T |γ (t)| dt
16
380
3
Initialize the set of candidate circles C := ∅
4
For each circle Ck If distance(A, Mk ) + distance(Mk , B) ≤ l + 2rk then
5
C := C ∪ {Ck }
6 7
385
390
end
8
end
9
If C = ∅ then
10
t := 0, i := 1
11
While t ≤ T
12
Pi := γ(t)
13
Find the circle Ck ∈ C which is closest to Pi
14
and let di := distance(Pi , Mk ) − rk
15
If di < 0 then output ”INTERSECTION” and STOP end
16
ti := di /vmax
17
If ti < Δt then ti := Δt end
18
t := t + ti , i := i + 1
19 395
end
20 end
A preselection test is performed in lines (4) – (8) of the procedure. We construct an ellipse with focal points located at the start and end points of the curve and with major diameter equal to the length l of the curve. If a circle lies outside the ellipse then it certainly cannot intersect with the curve (see Figure 5 (a)). 400
The other circles are subject to an intersection test with the curve. This test comprises lines (9) – (20) of the procedure. The idea is illustrated in Figure 5 (b). In the first step, the smallest distance d1 from the start point P1 = γ(0) of the curve to all candidate circles is determined. Then t1 = d1 /vmax is the smallest possible time to reach that circle, presumed a straight line connection
405
from P1 to the circle and maximal velocity. We can proceed along the curve to the point P2 = γ(t1 ) without intersecting a circle. This is repeated until the end of the curve is reached or a point lies inside a circle. In case of intersection, termination of the procedure is ensured by introducing a minimal time step Δt. Note that a (potentially occurring) discretization error is neglectable in our
17
(a)
(b)
Figure 5: Checking the intersection of a curve with circles
410
application. Intersection tests with polygons can be constructed in a similar way. Due to lack of space we omit the details. 4.2. Data preprocessing A further measure to improve the efficiency of the algorithm is the separation
415
into preprocessing and online planning. The preprocessing covers all preparatory work, including the generation of a basic network which takes into account all a priori available information. In a first step we create a network with a huge number of configurations (vertices) and feasible tracks (edges). All tracks are eliminated which intersect
420
with obstacles within the operational area or with its outer boundary. The result is a network where all tracks are passable and feasible. Since threat areas are not considered in this phase, the tracks are not necessarily safe. In a second step, the network is thinned out in order to remove non-relevant vertices and edges. Only those vertices must be preserved which lie on a path
425
from the source vertex vP to the target vertex vQ . All other vertices are worthless and can be deleted from the network, together with the incident edges. Hence, a vertex vA remains in the network if (a) the vertex vA can be reached from source vertex vP , and (b) the target vertex vQ can be reached from vertex vA .
430
We determine the set of vertices that are reachable from the source vertex vP
18
and, in the reversed network, the set of vertices that are reachable from the target vertex vQ . The reversed network is obtained by reversing the directions of all edges (i.e. if a network G contains an edge (u, v) then the reverse of G contains the edge (v, u) and vice versa). The vertices which belong to both 435
sets remain in the network, all other vertices are deleted. Finding the set of reachable vertices in a network is a basic problem in graph theory. It can be solved very efficiently by simple graph-searching algorithms such as depth-first search or breadth-first search (see e.g. [6]). The basic networks for different operational areas can be generated well in
440
advance and stored in a database. 4.3. Online planning The online planning considers up-to-date information about threats. In contrast to the preprocessing phase, where computer running time is not crucial, the online planning phase should be very fast.
445
A straightforward idea would be to adjust the basic network to the current scenario, and then find a shortest path in the updated network. This means to eliminate all edges (tracks) which are not safe, i.e. which intersect with threat areas. However, investigating all edges is very time-consuming. Instead, only those edges are inspected which are visited during the execution
450
of the shortest-path algorithm. A variety of shortest-path algorithms have been proposed in the literature (see e.g. [6]). Particularly well suited for our purpose is the A∗ -algorithm [11] which can be considered as an extension of the famous algorithm of Dijkstra [7]. By performing a target-directed search, a significantly smaller number of vertices and edges are visited than by other search algorithms,
455
making it much more efficient. We adapt the A∗ -algorithm in two ways: by integrating the intersection test of tracks with threats and by using a specific heuristic distance function of the vertices. The input for the algorithm is a network with passable and feasible edges (tracks). The result is a shortest path consisting of safe edges only.
19
460
Modified A∗ -Algorithm Network G = (V, E, c) with c : E → R+
Input:
Source vertex vP and target vertex vQ Estimated distance h(v) from v to vQ , for each v ∈ V Shortest path of safe edges from vP to vQ
Output: 465
1
For each v ∈ V d(v) := ∞ , π(v) := N IL
2
470
475
3
end
4
S := {vP } , d(vP ) := 0
5
While S = ∅
6
Find u ∈ S with f (u) = d(u) + h(u) minimal
7
If u = vQ then STOP end
8
S := S \ {u}
9
For each v ∈ V with (u, v) ∈ E
10
Check whether edge (u, v) is safe
11
If (u, v) is not safe then Delete the edge from the network, E := E \ {(u, v)}
12 13
480
else If d(v) > d(u) + c(u, v) then
14
d(v) := d(u) + c(u, v)
15
π(v) := u
16
If v ∈ S then S := S ∪ {v} end
17
end
18 19
end end
20 end
485
The modification covers lines (10) - (12). Checking whether an edge is safe takes place by the intersection procedures for circles and polygons presented above. d(v) is the distance of a vertex v to the source vertex vP . π(v) is the predecessor of v in the shortest path from vP to v. The shortest path from vP to vQ is obtained by iteratively building the sequence vQ , π(vQ ), π(π(vQ )), ...
490
until the source vertex vP is reached, and reversing the order. The length of the path is d(vQ ).
20
Figure 6: Safe route in a scenario with minefields and solitary mines
The A∗ -algorithm traverses the network by visiting first those vertices which are expected to lead to the target very quickly. In order to find the most promising vertex, each vertex v is assigned a value f (v). It estimates the length 495
of the path from the source to the target via vertex v. The algorithm always selects the vertex with smallest value of f . This realizes a target-directed search. f (v) consists of the exact distance d(v) already traveled from the source vertex to v and an estimate h(v) of the distance from v to the target vertex. The heuristic function h(v) must not overestimate the distance. The Euclidean
500
distance is often used as a suitable heuristic since the actual distance in the network cannot be shorter than the straight line connection. Note that, the lower h(v) is, the more vertices are visited which makes the algorithm slower. If h(v) is 0, then only d(v) plays a role, and the A∗ -algorithm turns into the algorithm of Dijkstra.
505
We use an estimate which is more accurate than the Euclidean distance, namely the exact distance in the basic network (where threats have not yet been considered). In the extreme, if there is no threat on the shortest path, the algorithm follows the optimal path and never expands anything else. This means that only edges on the shortest path, or edges incident to vertices of the
21
Figure 7: Safe routes obtained by different route planning algorithms
510
shortest path are investigated, which keeps the effort extremely small. The distances are calculated in the preprocessing phase by applying Dijkstra’s algorithm to the reversed network, with vQ being the source vertex. As mentioned before, Dijkstra’s algorithm is obtained from the (original) A∗ algorithm by setting h(v) = 0 for all vertices v ∈ V . If the termination condition
515
in line (7) is omitted, then we obtain shortest paths from the source vertex to all other vertices of the network.
5. Computational Results The algorithm has been implemented with Matlab R2013b on a standard PC with 3.5 GHz Intel Xeon Processor E3-1270v2 and 8 GB of RAM. Figure 520
6 illustrates a safe route in a typical scenario (a sea passage to a harbor) with two minefields and several solitary mines with different safety radii. Figure 7 shows two routes in a small scenario with ten solitary mines and two obstacles. The ship is a merchant vessel with turning circle diameter 1090 m. The short route is obtained by our algorithm, the other one by an algorithm
525
derived from Ari et al. [4] which relies on a regular grid (lattice) discretization. The ship turn constraints are realized in the following way: the track must first follow a straight line of about 450 m length (called the first leg) along the grid lines or along the diagonals, then performs a 45◦ turn, and then again follows a straight line of the same length (called the second leg). The polygonal line is 22
Figure 8: Safe routes with high risk threshold (red) and low risk threshold (blue)
530
smoothed after the optimization by replacing the inner halves of the two legs by the arc segment of the minimum turning circle. The resolution of the grid can be refined by subdividing the leg lengths into a larger number of edges. However, independent of the selected resolution, the method cannot find a route which is close to the optimum.
535
The test case in Figure 8 shows the results of our route planning algorithm for two different risk thresholds, with the objective to minimize traveling time. The scenario covers an operational area of 4 km × 3 km with five obstacles and a large number of solitary mines. The small circles correspond to a high acceptable risk, the large circles to a lower risk. The associated basic network has
540
4 000 vertices (configurations) and approximately 300 000 edges (tracks). The computation time of the preprocessing phase is about 2.2 minutes, the online planning requires less than 1 second. Experience shows that the calculated routes tend to get more winding if the number of mines and obstacles increases, or if the acceptable risk of the
545
threats is reduced (which inflates the threat areas). From a practical point of view, permanent changes of the ship’s course are not desired. One reason is that routes with fewer turns are easier to execute accurately, and navigational accuracy is critical when crossing an area affected by mines. Another reason is
23
Figure 9: Safe routes of short length (blue) and with few maneuvers (red)
simply not to overstress ship and crew. Hence, in those cases it seems advisable 550
to modify the optimization objective. Instead of minimizing the route length or traveling time, we minimize the number of maneuvers. Since each edge (track) implies one single maneuver, the number of edges in a path from the start to the target is identical to the number of maneuvers in the route. We change the costs of the edges in the network to the constant value 1
555
and determine a path with the smallest number of edges. This is again done by applying the modified A∗ -algorithm from above. The heuristic function h(v) is the minimal number of edges in a path from v to the target in the basic network. h(v) is calculated in the preprocessing phase using a simple breadth-first search. Figure 9 shows the obtained routes with minimum length and with minimal
560
number of maneuvers in a scenario with solitary mines and minefields. The basic network consists of 4 000 vertices and over 500 000 edges. The preprocessing time is about 3.1 minutes, the online planning time little more than 1 second. To evaluate the performance of the algorithm, we tested it on randomly generated scenarios. The size of the operational area is 6 km× 4.5 km. A varying
565
number of mines is randomly distributed within this area. The safety radius of the mines is 200 m. The ship is a naval vessel with cruising speed 23 knots and turning circle diameter 580 m.
24
Number of vertices
2 000
4 000
6 000
8 000
10 000
Number of edges
471 868
1 904 895
4 267 466
7 602 529
11 878 617
Route length [m]
8 370
8 143
7 765
7 687
7 642
Preprocessing [min]
1.3
5.2
12.5
22.1
35.7
Online planning [sec]
0.8
2.5
5.0
8.3
16.2
Table 1: Computational results for increasing network size
The impact of the network size on the quality of the solutions and the required calculation times are listed in Table 1. The results are the average values 570
for 10 random scenarios, each with 50 randomly distributed mines. The number of edges refers to the basic networks which are generated in the preprocessing phase (i.e. the edges correspond to tracks which are passable and feasible, but not necessarily safe). The objective is to minimize the route length. As expected, the length decreases if the density of the network increases. However,
575
at some point the effort grows overproportionally for a very modest improvement of the solution. As a simple rule-of-thumb, 300 vertices per square kilometer seem to be a reasonable compromise between accuracy and processing time. Table 2 contains results for different optimization objectives. The associated basic networks have 10 000 vertices, the number of mines is increasing. We
580
minimized the route length (see the left section of the table) and the number of maneuvers (right section). It becomes evident that, the shorter the route, the
Number
Route
Number of
Route
Number of
of mines
length [m]
maneuvers
length [m]
maneuvers
50
7 642
6.8
9 454
3.2
60
7 746
7.7
9 479
4.1
70
7 810
8.0
9 707
4.6
Table 2: Minimizing route length vs. minimizing number of maneuvers
25
Figure 10: Safe routes with minimized length (blue), minimized number of maneuvers (red) and minimized linear combination (black)
more maneuvers are included. On the other hand, a route with few maneuvers is usually rather long. A compromise between both objectives can be achieved by optimizing a 585
linear combination of length and number of maneuvers. The usual approach is to calculate the mean length (mean cost) of all edges
1 c(vA , vB ) c¯ = |E| (vA ,vB )∈E
and define the length ratio (cost ratio) of an edge as c∗ (vA , vB ) =
c(vA , vB ) . c¯
The weighted edge cost used for optimization is α · c∗ (vA , vB ) + (1 − α) · 1 The parameter α ∈ [0, 1] allows a balance between the objectives and is selected 590
by the mission planner. α = 1 means to minimize the route length, α = 0 is to minimize the number of maneuvers. Figure 10 shows a test case with a route of minimum length, a route with a minimal number of maneuvers, and a route obtained by setting α = 0.5. 26
6. Conclusion 595
In this paper we presented an algorithm to find safe navigation routes in an environment with threats and obstacles. The former are given in terms of solitary mines with sufficiently reliable information and minefields with very vague information. The focus of our approach is on finding a trajectory that considers the technical and operational restrictions of the ship’s movement. We
600
are able to determine a route of minimum length, minimum traveling time, with smallest number of maneuvers, or some compromise. The routes allow a large variety of maneuvers of different strengths and without turn degree limitations. A realistic ship track model ensures that the ship will be able to precisely follow the route - an indispensable requirement in an environment threatened by sea
605
mines. In contrast to other network-based methods, our algorithm is not restricted to regular grid graphs but depends on randomized networks where each path corresponds to a continuous-curvature track. In this way, high-quality solutions can be achieved if the networks are sufficiently large. Moreover, no path smooth-
610
ing is needed after optimization. The planning process can be performed very efficiently by separation into a preprocessing and an online planning phase. The algorithm is flexible enough as to allow the integration of more sophisticated risk models known from the literature. Future work could refine the shape of the navigation routes. The calculated
615
routes are concatenations of tracks, each consisting of straight line segments and two symmetric clothoid curves. One could allow more complex tracks, e.g. by combining a clothoid, a circular arc and a clothoid, or joining several nonsymmetric clothoids. References
620
References [1] M. Abramowitz, I.A. Stegun, eds. Handbook of Mathematical Functions with Formulas, Graphs, and Mathematical Tables. Dover, New York, 1972. 27
[2] American Bureau of Shipping, Guide for Vessel Maneuverability, Houston, TX, 2006. 625
[3] M. An, J.T. Cobb, B. Shenefelt, R. Tolimieri, Advances in group filter applications to sea mine detection, Proceedings of the IEEE OCEANS 2006 Boston Conference and Exhibition, 1-5, 2006. [4] I. Ari, V. Aksakalli, V. Aydoˇ gdu, S. Kum, Optimal ship navigation with safety distance and realistic turn constraints, European Journal of Opera-
630
tional Research 229 (2013) 707-717. [5] J.F. Bekker, J.P. Schmid, Planning the safe transit of a ship through a mapped minefield, ORION 22:1 (2006) 1-18. [6] T.H. Cormen, C.E. Leiserson, R.L. Rivest, C. Stein, Introduction to Algorithms. 3rd edn. MIT Press, Cambridge, 2009.
635
[7] E.W. Dijkstra, A note on two problems in connexion with graphs, Numerische Mathematik 1 (1959) 269-271. [8] K. Fagerholt, S. Heimdal, A. Loktu, Shortest path in the presence of obstacles: an application to ocean shipping, Journal of the Operational Research Society 51 (2000) 683-688.
640
[9] D. Ferguson, M. Likhachev, A. Stentz, A guide to heuristic-based path planning, Proceedings of ICAPS Workshop on Planning under Uncertainty for Autonomous Systems, AAAI, 9-18, 2005. [10] C. Goerzen, Z. Kong, B. Mettler, A survey of motion planning algorithms from the perspective of autonomous UAV guidance, Journal of Intelligent
645
& Robotic Systems 51 (2010) 65-100. [11] P.E. Hart, N.J. Nilsson, B. Raphael, A formal basis for the heuristic determination of minimum cost paths, IEEE Transactions on System Science and Cybernetics 4:2 (1968) 100-107.
28
[12] Y. Kanayama, N. Miyake, Trajectory generation of mobile robots, Robotics 650
Research 3 (1985) 334-340. [13] J.C. Latombe, Robot Motion Planning, Kluwer Academic Publishers, Boston, 1991. [14] S.M. LaValle, Planning Algorithms, Cambridge University Press, 2006. [15] P.C. Li, Planning the optimal transit for a ship through a mapped minefield,
655
Thesis, Naval Postgraduate School, Monterey, California, 2009. [16] H.C. Lu, Using expert systems in mine warfare, Thesis, Naval Postgraduate School, Monterey, California, 1991. [17] D. Mohovi´c, R. Mohovi´c, I. Rudan, Ship track and speed model in case of steering gear breakdown with rudder remaining fix at non zero angle,
660
Brodogradnja 63:2 (2012) 117-124. [18] C.D. Piatko, C.P. Diehl, P. McNamee, C.L. Resch, I.J. Wang, Stochastic Search and Graph Techniques for MCM Path Planning, Proceedings of the SPIE Volume 4742, Detection and Remediation Technologies for Mines and Minelike Targets, 583-593, 2002.
665
[19] C.D. Piatko, C. Priebe, L. Cowen, I.J. Wang, P. McNamee, Path Planning for Mine Countermeasures Command and Control, Proceedings of the SPIE Volume 4394, Detection and Remediation Technologies for Mines and Minelike Targets, 836-843, 2001. [20] W.H. Press, S.A. Teukolsky, W.T. Vetterling, B.P. Flannery, Numerical
670
Recipes. The Art of Scientific Computing. Cambridge University Press, 3nd edition, 2007. [21] C. Resch, C. Piatko, F.J. Pineda, J. Pistole, I.J. Wang, Path planning for mine countermeasures, Proceedings of the SPIE Volume 5089, Detection and Remediation Technologies for Mines and Minelike Targets, 1279-1286,
675
2003. 29
[22] J.J. Rios, Naval mines in the 21st century: can NATO navies meet the challenge? Thesis, Naval Postgraduate School, Monterey, California, 2005. [23] M. Zabarankin, S. Uryasev, P. Pardalos, Optimal risk path algorithms, In: R. Murphey, P. Pardalos, Cooperative control and optimization, vol. 66, p. 680
271-300, Kluwer Academic Publishers, Dordrecht, 2002. [24] B. Zerr, J. Fawcett, D. Hopkin, Adaptive algorithm for sea mine classification, 3rd international conference and exhibition on underwater acoustic measurements: technologies and results, 319-326, 2009.
30
Highlights x x x x x
We investigate the problem of finding a safe transit of a ship through an area threatened by sea mines. We present a network-based algorithm which takes into account the technical and operational restrictions of the ship’s movement. The network is especially tailored to the maneuverability of the ship. Each path in this network represents a continuous-curvature track based on combinations of clothoids and straight line segments. The algorithm allows minimizing the route length, the traveling time, the number of maneuvers, or other objectives.
S0377-2217(14)00618-3 http://dx.doi.org/10.1016/j.ejor.2014.08.002 EOR 12444
To appear in:
European Journal of Operational Research
Received Date: Accepted Date:
17 March 2014 2 August 2014
Please cite this article as: Babel, L., Zimmermann, T., Planning safe navigation routes through mined waters, European Journal of Operational Research (2014), doi: http://dx.doi.org/10.1016/j.ejor.2014.08.002
This is a PDF file of an unedited manuscript that has been accepted for publication. As a service to our customers we are providing this early version of the manuscript. The manuscript will undergo copyediting, typesetting, and review of the resulting proof before it is published in its final form. Please note that during the production process errors may be discovered which could affect the content, and all legal disclaimers that apply to the journal pertain.
Planning safe navigation routes through mined waters Luitpold Babel∗ , Thomas Zimmermann Institut f¨ ur Mathematik und Informatik, Fakult¨ at Betriebswirtschaft, Universit¨ at der Bundeswehr M¨ unchen, D-85579 Neubiberg, Germany
Abstract In this paper we investigate the problem of finding a safe transit of a ship through areas threatened by sea mines. The aim is to provide decision-making support by a tool that can be integrated into a naval command and control system. We present a route finding algorithm which avoids regions of risk higher than a given threshold. The algorithm takes into account the technical and operational restrictions of the ship’s movement. It allows to minimize the route length, the traveling time, the number of maneuvers, or other objectives. The basic idea is to embed a network in the operational area and compute a least-cost path. Instead of using a regular grid graph which strongly restricts the types of maneuvers and necessitates a path smoothing after optimization, we design a network which is especially tailored to the maneuverability of the vessel. Each path in this network represents a continuous-curvature track based on combinations of clothoids and straight line segments. The approach allows a large variety of maneuvers, hence high-quality solutions are achievable provided a sufficiently dense network. Keywords: Graph theory, Path planning, Navigation routes, Sea mines, Clothoid curves
∗ Corresponding
author. Tel.: +49 89 6004 3267. Email address: [email protected]
Preprint submitted to European Journal of Operational Research
August 6, 2014
1. Introduction The rapid increase in the proliferation of sea mines presents a unique challenge to maritime security. Sea mines are self-contained explosive devices placed in water to destroy surface ships or submarines. Since mines are inexpensive 5
they form an effective means to deny the use of large areas. For over a century, they have been used as a strategic weapon in military conflicts to blockade ports, waterways or other maritime zones. The naval forces must be prepared to counter the danger posed by mines. This, in particular, includes to find safe passages in order to sustain commercial and military shipping operations.
10
Planning a navigation route through an area threatened by sea mines is closely related to the classical problem of finding a collision-free path of a vehicle in an environment with obstacles. The latter, also known as path planning or trajectory planning problem, has long been the subject of intensive research. We refer to the reviews of Latombe [13] and LaValle [14].
15
Numerous ideas have been investigated to resolve the path planning problem. The most promising include the potential field method, cell decomposition, the roadmap method, the mass-spring-damper method, mixed integer linear programs, and several network-based methods with discretizations of the space by visibility graphs, Voronoi diagrams or regular grids. Newer approaches include
20
the rapidly-exploring random tree method, model predictive control methods, and mathematical programming methods. Comprehensive reviews are given by Goerzen et al. [10] and Ferguson et al. [9]. While previous research mainly focused on vehicles moving on ground (mobile robots, autonomous vehicles) or in air (unmanned aerial vehicles, drones,
25
autonomous helicopters), less work has been published for sea vehicles. Although the problems are closely related the proposed algorithms are not directly applicable since, from a technological point of view, ships behave quite differently than ground vehicles or aircraft. Moreover, shipping may imply specific operational requirements.
30
An important fact is that continuous methods can only find locally optimal
2
solutions which may be far from globally optimal. On the other hand, networkbased methods give an approximate problem, but means that one can find a globally optimal solution (see e.g. Zabarankin et al. [23]). Furthermore, networkbased methods appear more promising for realizing specific technological and 35
operational constraints. Fagerholt et al. [8] present an algorithm to find a shortest route for a ship located at sea to a destination port, in the presence of polygonal obstacles. It relies on the partial computation of a visibility graph and solves the problem using two variants of the algorithm of Dijkstra.
40
In a series of papers, Piatko et al. [18, 19, 21] present path planning techniques to find low-risk paths of a ship through a minefield. They discretize the operational area by evenly distributing vertices, with eight edges emerging from each vertex (along the axes and the diagonals), hence generating a regular grid graph. The authors develop probabilistic risk models and identify different path
45
quality attributes. The problem is then solved by defining appropriate costs on the edges of the graph and calculating a least-cost path. A similar method has been suggested by Bekker and Schmid [5]. They also use a regular grid graph, and remove all edges which lie within some safety radius around the mines. The remaining edges get as a cost either the length,
50
or a risk value depending on the distance to the mines, or a mixture of both. Again a least-cost path has to be determined. Their approach has been extended by Li [15] who added longer edges to the graph in order to improve the solution. A drawback of the above methods is that important factors including the characteristics of the ship (maneuverability, turning circle, etc.) or restrictions
55
which are implied by operational reasons (e.g. no sharp maneuvers due to passenger or cargo safety) are not considered. Each method provides a polygonal line, possibly with sharp turns or extreme maneuvers, which the ship may not be able to follow. Turn constraints have been successfully implemented in a recent research of
60
Ari et al. [4] where the goal is to find a shortest path through an area with obstacles. Their approach is based on a regular grid as above, however, they 3
use a vertex replication technique where each vertex is split into copies labeled by the direction the ship is coming from. This incorporates navigation history into the present location and allows to impose realistic ship turn constraints. 65
Moreover, the operational area can be discretized at any desired resolution. Indeed, there is still room for improvement since ship maneuvers are limited to 45◦ turns with full rudder deflection (which might imply that the solution is far from optimal) and a path smoothing is required after optimization (which might cause collision with an obstacle if not provided with a buffer zone).
70
The aim of our work is to develop an algorithm which finds a safe route through a mined area. Particular attention will be given to the technical capabilities and operational restrictions of the ship’s movement. Moreover, the maritime environment including shallows, islands, obstacles, etc. will be considered. The algorithm is intended to be integrated as a decision-making tool into
75
a naval command and control system. We propose a network optimization approach, focusing on the generation of a special randomized network which is adapted to the maneuverability of the ship. The paths in the network correspond to smooth (continuous-curvature) tracks consisting of clothoids and straight line segments. The approach allows to set
80
restrictions on the track including (non-symmetric) turn radius constraints, thus avoiding too sharp maneuvers, but also other restrictions such as the exclusion of long-lasting weak maneuvers (which might be undesirable due to operational reasons). The change of velocity during turning maneuvers is taken into account. The random structure of the network allows a large variety of maneuvers of
85
different strengths and without turn degree limitations. By applying standard shortest-path methods, we can determine routes which avoid regions of risk higher than a given threshold. With appropriate costs defined on the edges, different objectives can be considered such as minimizing route length, traveling time, the number of maneuvers, or some compromise of
90
these. Like all network-based methods, our approach is of a heuristic nature. On the other hand, the quality of the solutions can be continually improved by increasing the size of the network, as illustrated in Section 5. 4
In the next section, a mathematical model of the ship’s movement is presented, along with a brief description of the threat model and the optimization 95
problem. Section 3 contains the basic idea of the network generation and the route planning algorithm. In Section 4 we discuss how to speed up the algorithm by using efficient intersection procedures, and by separating the planning into preprocessing and online planning. Computational results are presented in Section 5. We conclude with final remarks in Section 6.
100
2. Mathematical modeling 2.1. Ship track model Our first task is to establish a realistic model of the movement of a ship. From an operational point of view, a planned ship track should fulfill certain requirements. First, the track should be easy to navigate, i.e. it should allow easy
105
handling and control of the ship by the crew. To the extent possible, long-lasting maneuvers and extreme course alterations should be avoided. Finally, the track’s curvature should be continuous in order to guarantee a smooth navigation. This is because following a planned track should not unduly burden the ship and the crew. These aspects become crucial in transporting people or dangerous goods.
110
Comfort and safety can only be ensured for continuous-curvature trajectories. A clothoid, also known as Euler or Cornu spiral, is a curve whose curvature changes linearly with its length [1]. It is well known that a vehicle following the curve at constant speed will have a constant change of centripetal acceleration. Hence, the clothoid is a perfect representation of a vehicle traveling at constant
115
speed with a constant rate of change of the steering angle [12]. The linear steering motion makes such a track comfortable for vehicle and crew. It has been shown in [17] that a clothoid very well approximates the curve of a turning ship up to a turning angle of about 100◦ to 120◦ , although the speed is not constant. In the first phase of the maneuver, when the rudder
120
is deflected, the track is close to a straight line since the ship starts turning very slowly. Then the curvature increases and the track gradually turns into a
5
Figure 1: Phases of ship turning, turning circle and clothoid
circle. During this part of the turn, the speed of the ship decreases. If the angle exceeds 120◦ the ship moves with constant speed in a circle with constant radius. The characteristic phases of a turning ship and the approximating clothoid are 125
sketched in Figure 1. Our goal is to realize the ship track as a combination of straight line segments and clothoids. The transitions between the track segments have to be smooth. A clothoid can be expressed as a parametrized curve ⎞ ⎛ x(t) ⎠ , 0≤t≤T γ(t) = ⎝ y(t) with domain [ 0, T ] and t √ πτ 2 x(t) = a π dτ , cos 2 0
130
√ y(t) = a π
t
sin 0
πτ 2 dτ . 2
It is well known and straightforward to verify that the arc length and the curvature of the curve are √ s = a πt ,
√ κ =
π t a
respectively. Hence κ =
1 s a2
i.e. the curvature is proportional to the arc length. The parameter a represents a scale factor which defines the rate of change of the curvature and hence the 6
135
size of the clothoid. The Fresnel integrals t t πτ 2 πτ 2 dτ , S(t) = dτ C(t) = cos sin 2 2 0 0 do not have an analytical solution. According to the power series expansions of cosine and sine cos(t) = 1 −
t2 t4 t6 + − + ... 2! 4! 6!
sin(t) = t −
t5 t7 t3 + − + ... 3! 5! 7!
they can be solved numerically using 1 π 4 t9 1 π 6 t13 1 π 2 t5 + − + ... 2! 2 5 4! 2 9 6! 2 13 1 π 3 t7 1 π 5 t11 1 π 7 t15 π t3 − + − + ... S(t) = 2 3 3! 2 7 5! 2 11 7! 2 15
C(t) = t − 140
For more details, particularly with regard to convergence, numerical stability and efficient computation, we refer to the standard text book [20]. We define a configuration as a combination of a position and a direction in the plane. For the route planning algorithm that will be presented in Section 3, 145
we need to solve a huge number of subproblems (SP) of the following type: Given two configurations, find a ship track from the first (initial) (SP)
to the second (final) configuration. The velocity at both endpoints should equal the cruising speed v0 .
To construct such a track from a configuration A to a configuration B, we first connect A and B by two line segments and a circle segment. We confine ourselves to circle segments with angle smaller than 180◦ . Otherwise 150
(SP) is assumed to have no solution. Hence the ray originating from A intersects with the ray ending in B. A circle with some predefined radius r is placed tangential to the line segments between A and B (see Figure 2). Note that the transition from a straight line path to a tangential circular path implies a sudden centripetal acceleration (which a ship is not able to follow). This is
155
avoided by replacing the circle segment by two symmetric clothoid curves. They
7
Figure 2: Connecting two configurations with two lines and a circle segment
are adjusted in such a way that the acceleration is continuous along the whole track and equal to zero at both ends. Consider a circle segment with radius r and angle 2α as depicted in Figure 3. The first clothoid starts in the origin of the coordinate system tangential to 160
the x-axis and ends in point S, perpendicular to the angle bisector. The second clothoid is mirrored at the bisector line. We denote the first clothoid by γ(t). The parameter a and the range T of γ(t) is then calculated as follows. The tangent vector of γ(t) is
⎞
⎞ ⎛ πt2 cos √ 2 ⎠ ⎠ = a π⎝ γ (t) = ⎝ . 2 y (t) sin πt2 ⎛
x (t)
Since the slope of γ(t) in S equals α we get ⎛ ⎞ ⎞ ⎛ 2 cos πT2 cos α ⎝ ⎠ . ⎠ = ⎝ 2 sin α sin πT2 165
This implies
T =
2α . π
The point S lies on the bisector line y = − cot α · x + r. Hence cot α · x(T ) + y(T ) = r and, by substituting x(T ) and y(T ),
T T √ πτ 2 πτ 2 a π cot α dτ + dτ = r. cos sin 2 2 0 0 8
Figure 3: Replacing a circle segment by two symmetric clothoids
The expression in brackets, abbreviated with I, is evaluated numerically. Then r a = √ πI which determines the parameter of the clothoid. 170
2.2. Ship characteristics The curvature of the curve is maximal at connection point S of the clothoids. At that point
√ κmax =
π 1√ T = 2π α I . a r
Since the ship must be able to follow the track, κmax must not exceed the curvature of the minimum turning circle (which is related to the maximal rudder 175
deflection). If Rturn denotes the minimum turning radius then 1√ 1 2π α I ≤ r Rturn must hold. We define Rmin =
√
2π α I Rturn .
For the construction from above, we restrict the radius r of the circle segment by Rmin ≤ r ≤ Rmax . The upper bound Rmax is the radius of the largest circle that can be placed 180
tangential to the line segments. Note that this bound might be used to forbid 9
long-lasting weak maneuvers, if those are not desired due to operational reasons. The radius r is selected randomly within the bounds in order to enable a large variety of maneuvers. Alternatively, to control the type of the occuring maneuvers, one might also define the radius deterministically, e.g. by always selecting 185
the smallest possible value which provides short but strong maneuvers, or the largest possible value which provides short length curves. If Rmin > Rmax then the subproblem (SP) is assumed to have no solution. One might also define separately restrictions for left and right turns if left (port) and right (starboard) turn radii are different (which indeed is true for
190
most ships, see [2]). In that case, the lower bound Rmin has to be replaced by l r two different values Rmin and Rmin , according to the turn direction.
It remains to create a speed change model. During turn, the drag of the rudder and the sideways drift of the ship result in a progressive loss of speed, even though the engine revolutions are maintained constant (see e.g. [2]). Ex195
perience shows that most ships with full rudder deflection lose about one-third of their original speed after a 90-degree turn. We assume that the ship travels with cruising speed v0 along the straight line segments. The ship decelerates during the first clothoid where the curvature increases, achieves its minimum speed v1 at the end of the first clothoid
200
(at symmetry point S), and accelerates during the second clothoid where the curvature decreases, until it again achieves cruising speed v0 . We suppose a velocity curve with average speed (v0 + v1 )/2 and define v1 as a function of the turning radius r and turning angle α by v1 = v(α, r) = v0 (1 −
1 α Rmin ) 3 π/2 r
where α is restricted by 90◦ . The largest decrease in speed occurs for maximal 205
angle and minimal radius where π 2 v1 = v( , Rmin ) = v0 . 2 3 Since the rate of deceleration varies between different types of ship, one might have to adapt the model accordingly. 10
2.3. Threat model The operational area is a marine area which is passable by ships. The area 210
possibly includes islands, reefs, shallows and nonnatural impassable sites such as shipwrecks and other obstacles. The operational area depends also on the draught of the ship, the tides and other environmental conditions. Modeling is done by defining appropriate (non-convex) polygons. A large outer polygon defines the shoreline or the boundary of the mission area, cut-out inner polygons
215
represent the impassable regions. The threat areas consist of solitary mines and minefields. With the former we mean single mines for which accurate information is available. The latter are aggregations of mines with imprecise or vague information. For a classification of mine types, mine counter measures and methods of mine delivery, see e.g.
220
[16] and [22]. The detection and mapping of sea mines is discussed e.g. in [3] and [24]. Solitary mines are characterized by their location and damage radius. The damage radius is mainly determined by the explosive charge of the mine. However, it also depends on technical ship data including signature and vulnerability.
225
The signature of a ship consists of those emitted signals (e.g. acoustic, magnetic, visual, infrared) which can be used by sensors to detect the ship’s position. The vulnerability of a ship is based on its construction, hull material, armour, etc. These data determine the effective radius of a mine which, under no circumstances, may be entered. Otherwise the ship will sustain damage. The radius
230
can be increased by additional safety margins to take into account uncertainty of the position of the mine and the navigational error of the ship. In contrast to the mapped mines described above, there is only inaccurate information for minefields. The individual positions (and possibly the types) of the mines are not precisely known. Minefields are modeled by polygons which
235
include the supposed positions of the mines with an appropriate safety margin. They can be considered as forbidden areas where the ship is not allowed to pass through. Risk models for mines have been discussed in several papers, e.g. [5, 19, 21]. 11
They may rely on information from different sensors, available actuation and 240
damage curve data and other sources. Mine detection, classification uncertainty, and damage are then represented by suitable probability distributions. Since a full treatment of all relevant data would lead to highly complex models, some simplifications seem inevitable. The point of closest approach to the closest mine is used in [5] to determine
245
the risk of a path. With D denoting the minimal distance, the risk is defined proportional to any user selectable inverse power of D, such as 1/D or 1/D2 . A similar point risk model is used in [19]. The probability of being damaged at each location is an inverse function of the distance to the location of the mine. The risk associated with each path is then estimated by summing the
250
point risks along the path. A more sophisticated model is introduced in [21] where the location of a mine is represented by a normal distribution. This is combined with a triggering model for the mine which is given by a Gaussian-like distribution. 2.4. Optimization problem
255
We call a navigation route passable if it proceeds completely within the operational area. Such a route is clear of obstacles and can be traveled without the danger of collision. A route is called feasible if it considers the technical and operational limitations of the ship, i.e. the ship is able to follow the route. Safety will be expressed by the level of risk that the mission planner or head
260
of operations is willing to bear. It is reflected by a safety margin which is added to the effective radius of each mine. For that purpose, an appropriate risk function f (D) has to be defined which is inversely proportional with distance D to the mine. If pmax denotes the maximum acceptable risk then the solution of the equation f (D) = pmax provides the radius around the mine which can
265
be considered as a forbidden area for the ship. We choose the risk function proportional to 1/D, noting that it can easily be replaced by a more sophisticated model (see [5, 21]). A similar idea is realized for minefields. The risk for the ship reduces if the distance to the minefield increases. The associated 12
polygon is enlarged by a buffer zone (see [4]) until the maximum acceptable risk 270
is achieved. The ship routing problem discussed in this paper is to determine a passable and feasible route of minimum length or minimum traveling time from an initial configuration (starting point with direction) to a final configuration (destination point with direction) which avoids regions with risk higher than a given
275
threshold. The above formulation aims to minimize one criterion (length or time) under a hard bound on the other criterion (the risk). An alternative approach would be to determine a route which minimizes the risk, or a combination of risk and length (time). Both problem definitions have been extensively discussed e.g. in
280
[18, 19]. Let us mention that our algorithm, although designed for the constraintbased approach, can easily be extended to find a passable and feasible route with smallest accumulated risk or combination of risk and length (time). For the calculation of the accumulated risk along a path we refer to [19].
285
3. Route planning The basic idea of our approach is to generate a network over the water surface in such a way that each path in this network represents a possible navigation route. The shortest navigation route in the network is then determined by applying standard methods of combinatorial optimization. The quality of the
290
solution increases with the density of the network. Network-based solutions are also proposed in other publications, e.g. in [5, 15, 18, 19, 21]. The common procedure is to evenly distribute vertices over an area and link close pairs of vertices by straight lines. This provides a discretization of the area in the form of a regular grid network. A drawback is that the
295
ship’s characteristic (turning radius, maneuverability, etc.) is not sufficiently considered. The resulting navigation route is a polygonal line which might be difficult, or even impossible to navigate. Moreover, if the number of vertices in
13
the network is small then the route might be far away from an optimal solution. With increasing number of vertices the route might become an annoying zigzag 300
course. Our network considers the technical and operational limitations of the ship. The network is not regular since its generation is affected by some random process. The vertices are not connected by straight lines but by curved lines consisting of two straight line segments and two symmetric clothoids in between
305
those segments. This provides a smooth track which allows convenient control of the ship. In contrast to the situation described above the track remains smooth, independent of the density of the network. The following algorithm scheme describes the generation of the network and the route planning.
310
Algorithm Find-Ship-Track Start configuration P , target configuration Q
Input:
Operational area, performance data of ship, threat areas Output:
315
Ship track γ(t) from P to Q
1
Create vertices vP and vQ for start configuration P and target configuration Q
2
Initialize vertex set V := {vP , vQ } and edge set E := ∅
3
Generate a huge number of configurations in the operational area by randomly choosing position and direction
4
For each configuration A Create a new vertex vA , V := V ∪ {vA }
5 320
6
end
7
For each pair of configurations (A, B)
8
Generate a feasible track γAB (t) from A to B consisting of clothoids and straight line segments (if existing) by solving the corresponding (SP)
9
Check whether γAB (t) is passable, i.e. the track is entirely within the operational area
325
10
Check whether γAB (t) is safe, i.e. the track avoids the threat areas
11
If γAB (t) is passable and safe then
12
Create a new edge (vA , vB ) with cost c(vA , vB ) :=
13
E := E ∪ {(vA , vB )}
14
|γAB (t)| dt,
Figure 4: Schematic representation of network
330
14
end
15 end 16 Find a shortest path in the network G = (V, E, c) from vP to vQ 17 Concatenate the tracks belonging to the path to a P - Q track γ(t)
The vertices of the network represent randomly generated configurations. Each 335
(directed) edge represents a track segment connecting two configurations. Each path is a concatenation of track segments. The construction of the track segments, described in the previous section as the solutions of subproblems (SP), guarantees that the transition between any two of them is smooth. More precisely, this means that each path corresponds to a twice continuously differen-
340
tiable curve. As a result of the checks performed in lines (9) and (10), each path in the network represents a passable, feasible and safe navigation route. The cost of an edge is the length of the associated track segment. It is √ obtained by adding twice the arc length s = a π T of a clothoid to the lengths of the straight line segments. Alternatively, if a route with smallest traveling
345
time has to be found, the cost of an edge is the time required to pass the track. This is twice the time 2s/(v0 + v1 ) to pass one clothoid plus the time for the straight line segments. The idea of the network generation is sketched in Figure 4. It shows a number of configurations (arrows) with connecting track segments. The scenario 15
350
contains solitary mines (circles), minefields (polygons) and other obstacles. For operational use, the algorithm Find-Ship-Track should deliver a solution within a couple of seconds. Unfortunately, a straightforward implementation is too slow. Hence further efforts are needed to speed up the algorithm.
4. Efficient implementation 355
4.1. Intersection procedures The main workload of the algorithm occurs in lines (8) – (10) where a huge number of tracks is generated and analyzed. The generation follows the description given in Sections 2.1 and 2.2. According to the construction, the tracks are feasible. Checking whether tracks are passable and safe is time-critical.
360
The former means to find out whether a track lies completely within the outer (polygonal) boundary of the operational area and does not collide with an inner (polygonal) obstacle, the latter whether a track passes through a (circular) solitary mine or through a (polygonal) minefield. Hence, we have to develop efficient procedures to check the intersection of a
365
track with circles and polygons. The first problem is addressed by solving the following more general question: Given an arbitrary curve γ(t) with varying speed (describing the motion of a vehicle) and a set of circles with varying radii, does the curve intersect with at least one circle ?
370
We propose the following procedure. Procedure Intersect-Curve-with-Circles Input:
Curve γ(t) with 0 ≤ t ≤ T , minimal time step Δt Maximal velocity or upper bound vmax Set of circles Ck with center Mk and radius rk
375
Output: 1 2
Message, if γ(t) intersects with a circle
Let A := γ(0) and B := γ(T ) be the start and end point of the curve, respectively Calculate the length of the curve l := 0T |γ (t)| dt
16
380
3
Initialize the set of candidate circles C := ∅
4
For each circle Ck If distance(A, Mk ) + distance(Mk , B) ≤ l + 2rk then
5
C := C ∪ {Ck }
6 7
385
390
end
8
end
9
If C = ∅ then
10
t := 0, i := 1
11
While t ≤ T
12
Pi := γ(t)
13
Find the circle Ck ∈ C which is closest to Pi
14
and let di := distance(Pi , Mk ) − rk
15
If di < 0 then output ”INTERSECTION” and STOP end
16
ti := di /vmax
17
If ti < Δt then ti := Δt end
18
t := t + ti , i := i + 1
19 395
end
20 end
A preselection test is performed in lines (4) – (8) of the procedure. We construct an ellipse with focal points located at the start and end points of the curve and with major diameter equal to the length l of the curve. If a circle lies outside the ellipse then it certainly cannot intersect with the curve (see Figure 5 (a)). 400
The other circles are subject to an intersection test with the curve. This test comprises lines (9) – (20) of the procedure. The idea is illustrated in Figure 5 (b). In the first step, the smallest distance d1 from the start point P1 = γ(0) of the curve to all candidate circles is determined. Then t1 = d1 /vmax is the smallest possible time to reach that circle, presumed a straight line connection
405
from P1 to the circle and maximal velocity. We can proceed along the curve to the point P2 = γ(t1 ) without intersecting a circle. This is repeated until the end of the curve is reached or a point lies inside a circle. In case of intersection, termination of the procedure is ensured by introducing a minimal time step Δt. Note that a (potentially occurring) discretization error is neglectable in our
17
(a)
(b)
Figure 5: Checking the intersection of a curve with circles
410
application. Intersection tests with polygons can be constructed in a similar way. Due to lack of space we omit the details. 4.2. Data preprocessing A further measure to improve the efficiency of the algorithm is the separation
415
into preprocessing and online planning. The preprocessing covers all preparatory work, including the generation of a basic network which takes into account all a priori available information. In a first step we create a network with a huge number of configurations (vertices) and feasible tracks (edges). All tracks are eliminated which intersect
420
with obstacles within the operational area or with its outer boundary. The result is a network where all tracks are passable and feasible. Since threat areas are not considered in this phase, the tracks are not necessarily safe. In a second step, the network is thinned out in order to remove non-relevant vertices and edges. Only those vertices must be preserved which lie on a path
425
from the source vertex vP to the target vertex vQ . All other vertices are worthless and can be deleted from the network, together with the incident edges. Hence, a vertex vA remains in the network if (a) the vertex vA can be reached from source vertex vP , and (b) the target vertex vQ can be reached from vertex vA .
430
We determine the set of vertices that are reachable from the source vertex vP
18
and, in the reversed network, the set of vertices that are reachable from the target vertex vQ . The reversed network is obtained by reversing the directions of all edges (i.e. if a network G contains an edge (u, v) then the reverse of G contains the edge (v, u) and vice versa). The vertices which belong to both 435
sets remain in the network, all other vertices are deleted. Finding the set of reachable vertices in a network is a basic problem in graph theory. It can be solved very efficiently by simple graph-searching algorithms such as depth-first search or breadth-first search (see e.g. [6]). The basic networks for different operational areas can be generated well in
440
advance and stored in a database. 4.3. Online planning The online planning considers up-to-date information about threats. In contrast to the preprocessing phase, where computer running time is not crucial, the online planning phase should be very fast.
445
A straightforward idea would be to adjust the basic network to the current scenario, and then find a shortest path in the updated network. This means to eliminate all edges (tracks) which are not safe, i.e. which intersect with threat areas. However, investigating all edges is very time-consuming. Instead, only those edges are inspected which are visited during the execution
450
of the shortest-path algorithm. A variety of shortest-path algorithms have been proposed in the literature (see e.g. [6]). Particularly well suited for our purpose is the A∗ -algorithm [11] which can be considered as an extension of the famous algorithm of Dijkstra [7]. By performing a target-directed search, a significantly smaller number of vertices and edges are visited than by other search algorithms,
455
making it much more efficient. We adapt the A∗ -algorithm in two ways: by integrating the intersection test of tracks with threats and by using a specific heuristic distance function of the vertices. The input for the algorithm is a network with passable and feasible edges (tracks). The result is a shortest path consisting of safe edges only.
19
460
Modified A∗ -Algorithm Network G = (V, E, c) with c : E → R+
Input:
Source vertex vP and target vertex vQ Estimated distance h(v) from v to vQ , for each v ∈ V Shortest path of safe edges from vP to vQ
Output: 465
1
For each v ∈ V d(v) := ∞ , π(v) := N IL
2
470
475
3
end
4
S := {vP } , d(vP ) := 0
5
While S = ∅
6
Find u ∈ S with f (u) = d(u) + h(u) minimal
7
If u = vQ then STOP end
8
S := S \ {u}
9
For each v ∈ V with (u, v) ∈ E
10
Check whether edge (u, v) is safe
11
If (u, v) is not safe then Delete the edge from the network, E := E \ {(u, v)}
12 13
480
else If d(v) > d(u) + c(u, v) then
14
d(v) := d(u) + c(u, v)
15
π(v) := u
16
If v ∈ S then S := S ∪ {v} end
17
end
18 19
end end
20 end
485
The modification covers lines (10) - (12). Checking whether an edge is safe takes place by the intersection procedures for circles and polygons presented above. d(v) is the distance of a vertex v to the source vertex vP . π(v) is the predecessor of v in the shortest path from vP to v. The shortest path from vP to vQ is obtained by iteratively building the sequence vQ , π(vQ ), π(π(vQ )), ...
490
until the source vertex vP is reached, and reversing the order. The length of the path is d(vQ ).
20
Figure 6: Safe route in a scenario with minefields and solitary mines
The A∗ -algorithm traverses the network by visiting first those vertices which are expected to lead to the target very quickly. In order to find the most promising vertex, each vertex v is assigned a value f (v). It estimates the length 495
of the path from the source to the target via vertex v. The algorithm always selects the vertex with smallest value of f . This realizes a target-directed search. f (v) consists of the exact distance d(v) already traveled from the source vertex to v and an estimate h(v) of the distance from v to the target vertex. The heuristic function h(v) must not overestimate the distance. The Euclidean
500
distance is often used as a suitable heuristic since the actual distance in the network cannot be shorter than the straight line connection. Note that, the lower h(v) is, the more vertices are visited which makes the algorithm slower. If h(v) is 0, then only d(v) plays a role, and the A∗ -algorithm turns into the algorithm of Dijkstra.
505
We use an estimate which is more accurate than the Euclidean distance, namely the exact distance in the basic network (where threats have not yet been considered). In the extreme, if there is no threat on the shortest path, the algorithm follows the optimal path and never expands anything else. This means that only edges on the shortest path, or edges incident to vertices of the
21
Figure 7: Safe routes obtained by different route planning algorithms
510
shortest path are investigated, which keeps the effort extremely small. The distances are calculated in the preprocessing phase by applying Dijkstra’s algorithm to the reversed network, with vQ being the source vertex. As mentioned before, Dijkstra’s algorithm is obtained from the (original) A∗ algorithm by setting h(v) = 0 for all vertices v ∈ V . If the termination condition
515
in line (7) is omitted, then we obtain shortest paths from the source vertex to all other vertices of the network.
5. Computational Results The algorithm has been implemented with Matlab R2013b on a standard PC with 3.5 GHz Intel Xeon Processor E3-1270v2 and 8 GB of RAM. Figure 520
6 illustrates a safe route in a typical scenario (a sea passage to a harbor) with two minefields and several solitary mines with different safety radii. Figure 7 shows two routes in a small scenario with ten solitary mines and two obstacles. The ship is a merchant vessel with turning circle diameter 1090 m. The short route is obtained by our algorithm, the other one by an algorithm
525
derived from Ari et al. [4] which relies on a regular grid (lattice) discretization. The ship turn constraints are realized in the following way: the track must first follow a straight line of about 450 m length (called the first leg) along the grid lines or along the diagonals, then performs a 45◦ turn, and then again follows a straight line of the same length (called the second leg). The polygonal line is 22
Figure 8: Safe routes with high risk threshold (red) and low risk threshold (blue)
530
smoothed after the optimization by replacing the inner halves of the two legs by the arc segment of the minimum turning circle. The resolution of the grid can be refined by subdividing the leg lengths into a larger number of edges. However, independent of the selected resolution, the method cannot find a route which is close to the optimum.
535
The test case in Figure 8 shows the results of our route planning algorithm for two different risk thresholds, with the objective to minimize traveling time. The scenario covers an operational area of 4 km × 3 km with five obstacles and a large number of solitary mines. The small circles correspond to a high acceptable risk, the large circles to a lower risk. The associated basic network has
540
4 000 vertices (configurations) and approximately 300 000 edges (tracks). The computation time of the preprocessing phase is about 2.2 minutes, the online planning requires less than 1 second. Experience shows that the calculated routes tend to get more winding if the number of mines and obstacles increases, or if the acceptable risk of the
545
threats is reduced (which inflates the threat areas). From a practical point of view, permanent changes of the ship’s course are not desired. One reason is that routes with fewer turns are easier to execute accurately, and navigational accuracy is critical when crossing an area affected by mines. Another reason is
23
Figure 9: Safe routes of short length (blue) and with few maneuvers (red)
simply not to overstress ship and crew. Hence, in those cases it seems advisable 550
to modify the optimization objective. Instead of minimizing the route length or traveling time, we minimize the number of maneuvers. Since each edge (track) implies one single maneuver, the number of edges in a path from the start to the target is identical to the number of maneuvers in the route. We change the costs of the edges in the network to the constant value 1
555
and determine a path with the smallest number of edges. This is again done by applying the modified A∗ -algorithm from above. The heuristic function h(v) is the minimal number of edges in a path from v to the target in the basic network. h(v) is calculated in the preprocessing phase using a simple breadth-first search. Figure 9 shows the obtained routes with minimum length and with minimal
560
number of maneuvers in a scenario with solitary mines and minefields. The basic network consists of 4 000 vertices and over 500 000 edges. The preprocessing time is about 3.1 minutes, the online planning time little more than 1 second. To evaluate the performance of the algorithm, we tested it on randomly generated scenarios. The size of the operational area is 6 km× 4.5 km. A varying
565
number of mines is randomly distributed within this area. The safety radius of the mines is 200 m. The ship is a naval vessel with cruising speed 23 knots and turning circle diameter 580 m.
24
Number of vertices
2 000
4 000
6 000
8 000
10 000
Number of edges
471 868
1 904 895
4 267 466
7 602 529
11 878 617
Route length [m]
8 370
8 143
7 765
7 687
7 642
Preprocessing [min]
1.3
5.2
12.5
22.1
35.7
Online planning [sec]
0.8
2.5
5.0
8.3
16.2
Table 1: Computational results for increasing network size
The impact of the network size on the quality of the solutions and the required calculation times are listed in Table 1. The results are the average values 570
for 10 random scenarios, each with 50 randomly distributed mines. The number of edges refers to the basic networks which are generated in the preprocessing phase (i.e. the edges correspond to tracks which are passable and feasible, but not necessarily safe). The objective is to minimize the route length. As expected, the length decreases if the density of the network increases. However,
575
at some point the effort grows overproportionally for a very modest improvement of the solution. As a simple rule-of-thumb, 300 vertices per square kilometer seem to be a reasonable compromise between accuracy and processing time. Table 2 contains results for different optimization objectives. The associated basic networks have 10 000 vertices, the number of mines is increasing. We
580
minimized the route length (see the left section of the table) and the number of maneuvers (right section). It becomes evident that, the shorter the route, the
Number
Route
Number of
Route
Number of
of mines
length [m]
maneuvers
length [m]
maneuvers
50
7 642
6.8
9 454
3.2
60
7 746
7.7
9 479
4.1
70
7 810
8.0
9 707
4.6
Table 2: Minimizing route length vs. minimizing number of maneuvers
25
Figure 10: Safe routes with minimized length (blue), minimized number of maneuvers (red) and minimized linear combination (black)
more maneuvers are included. On the other hand, a route with few maneuvers is usually rather long. A compromise between both objectives can be achieved by optimizing a 585
linear combination of length and number of maneuvers. The usual approach is to calculate the mean length (mean cost) of all edges
1 c(vA , vB ) c¯ = |E| (vA ,vB )∈E
and define the length ratio (cost ratio) of an edge as c∗ (vA , vB ) =
c(vA , vB ) . c¯
The weighted edge cost used for optimization is α · c∗ (vA , vB ) + (1 − α) · 1 The parameter α ∈ [0, 1] allows a balance between the objectives and is selected 590
by the mission planner. α = 1 means to minimize the route length, α = 0 is to minimize the number of maneuvers. Figure 10 shows a test case with a route of minimum length, a route with a minimal number of maneuvers, and a route obtained by setting α = 0.5. 26
6. Conclusion 595
In this paper we presented an algorithm to find safe navigation routes in an environment with threats and obstacles. The former are given in terms of solitary mines with sufficiently reliable information and minefields with very vague information. The focus of our approach is on finding a trajectory that considers the technical and operational restrictions of the ship’s movement. We
600
are able to determine a route of minimum length, minimum traveling time, with smallest number of maneuvers, or some compromise. The routes allow a large variety of maneuvers of different strengths and without turn degree limitations. A realistic ship track model ensures that the ship will be able to precisely follow the route - an indispensable requirement in an environment threatened by sea
605
mines. In contrast to other network-based methods, our algorithm is not restricted to regular grid graphs but depends on randomized networks where each path corresponds to a continuous-curvature track. In this way, high-quality solutions can be achieved if the networks are sufficiently large. Moreover, no path smooth-
610
ing is needed after optimization. The planning process can be performed very efficiently by separation into a preprocessing and an online planning phase. The algorithm is flexible enough as to allow the integration of more sophisticated risk models known from the literature. Future work could refine the shape of the navigation routes. The calculated
615
routes are concatenations of tracks, each consisting of straight line segments and two symmetric clothoid curves. One could allow more complex tracks, e.g. by combining a clothoid, a circular arc and a clothoid, or joining several nonsymmetric clothoids. References
620
References [1] M. Abramowitz, I.A. Stegun, eds. Handbook of Mathematical Functions with Formulas, Graphs, and Mathematical Tables. Dover, New York, 1972. 27
[2] American Bureau of Shipping, Guide for Vessel Maneuverability, Houston, TX, 2006. 625
[3] M. An, J.T. Cobb, B. Shenefelt, R. Tolimieri, Advances in group filter applications to sea mine detection, Proceedings of the IEEE OCEANS 2006 Boston Conference and Exhibition, 1-5, 2006. [4] I. Ari, V. Aksakalli, V. Aydoˇ gdu, S. Kum, Optimal ship navigation with safety distance and realistic turn constraints, European Journal of Opera-
630
tional Research 229 (2013) 707-717. [5] J.F. Bekker, J.P. Schmid, Planning the safe transit of a ship through a mapped minefield, ORION 22:1 (2006) 1-18. [6] T.H. Cormen, C.E. Leiserson, R.L. Rivest, C. Stein, Introduction to Algorithms. 3rd edn. MIT Press, Cambridge, 2009.
635
[7] E.W. Dijkstra, A note on two problems in connexion with graphs, Numerische Mathematik 1 (1959) 269-271. [8] K. Fagerholt, S. Heimdal, A. Loktu, Shortest path in the presence of obstacles: an application to ocean shipping, Journal of the Operational Research Society 51 (2000) 683-688.
640
[9] D. Ferguson, M. Likhachev, A. Stentz, A guide to heuristic-based path planning, Proceedings of ICAPS Workshop on Planning under Uncertainty for Autonomous Systems, AAAI, 9-18, 2005. [10] C. Goerzen, Z. Kong, B. Mettler, A survey of motion planning algorithms from the perspective of autonomous UAV guidance, Journal of Intelligent
645
& Robotic Systems 51 (2010) 65-100. [11] P.E. Hart, N.J. Nilsson, B. Raphael, A formal basis for the heuristic determination of minimum cost paths, IEEE Transactions on System Science and Cybernetics 4:2 (1968) 100-107.
28
[12] Y. Kanayama, N. Miyake, Trajectory generation of mobile robots, Robotics 650
Research 3 (1985) 334-340. [13] J.C. Latombe, Robot Motion Planning, Kluwer Academic Publishers, Boston, 1991. [14] S.M. LaValle, Planning Algorithms, Cambridge University Press, 2006. [15] P.C. Li, Planning the optimal transit for a ship through a mapped minefield,
655
Thesis, Naval Postgraduate School, Monterey, California, 2009. [16] H.C. Lu, Using expert systems in mine warfare, Thesis, Naval Postgraduate School, Monterey, California, 1991. [17] D. Mohovi´c, R. Mohovi´c, I. Rudan, Ship track and speed model in case of steering gear breakdown with rudder remaining fix at non zero angle,
660
Brodogradnja 63:2 (2012) 117-124. [18] C.D. Piatko, C.P. Diehl, P. McNamee, C.L. Resch, I.J. Wang, Stochastic Search and Graph Techniques for MCM Path Planning, Proceedings of the SPIE Volume 4742, Detection and Remediation Technologies for Mines and Minelike Targets, 583-593, 2002.
665
[19] C.D. Piatko, C. Priebe, L. Cowen, I.J. Wang, P. McNamee, Path Planning for Mine Countermeasures Command and Control, Proceedings of the SPIE Volume 4394, Detection and Remediation Technologies for Mines and Minelike Targets, 836-843, 2001. [20] W.H. Press, S.A. Teukolsky, W.T. Vetterling, B.P. Flannery, Numerical
670
Recipes. The Art of Scientific Computing. Cambridge University Press, 3nd edition, 2007. [21] C. Resch, C. Piatko, F.J. Pineda, J. Pistole, I.J. Wang, Path planning for mine countermeasures, Proceedings of the SPIE Volume 5089, Detection and Remediation Technologies for Mines and Minelike Targets, 1279-1286,
675
2003. 29
[22] J.J. Rios, Naval mines in the 21st century: can NATO navies meet the challenge? Thesis, Naval Postgraduate School, Monterey, California, 2005. [23] M. Zabarankin, S. Uryasev, P. Pardalos, Optimal risk path algorithms, In: R. Murphey, P. Pardalos, Cooperative control and optimization, vol. 66, p. 680
271-300, Kluwer Academic Publishers, Dordrecht, 2002. [24] B. Zerr, J. Fawcett, D. Hopkin, Adaptive algorithm for sea mine classification, 3rd international conference and exhibition on underwater acoustic measurements: technologies and results, 319-326, 2009.
30
Highlights x x x x x
We investigate the problem of finding a safe transit of a ship through an area threatened by sea mines. We present a network-based algorithm which takes into account the technical and operational restrictions of the ship’s movement. The network is especially tailored to the maneuverability of the ship. Each path in this network represents a continuous-curvature track based on combinations of clothoids and straight line segments. The algorithm allows minimizing the route length, the traveling time, the number of maneuvers, or other objectives.