VEO Weather Routing - an innovative algorithm for calculating routes based on weather, map and polar plots of yachts

Introduction

Mapping the route of a sailing yacht based on weather data has long been an important issue in the practice of sailing. So far, this task has been carried out mainly on the basis of planning "in the head" of the yacht captain, i.e. weather data and the sailor's knowledge / intuition regarding the use of wind and sea currents by a specific ship model. The mathematical, algorithmic version of determining sea routes was known more in the racing world, where every mistake results in a worse result in the race. So far, isochrone algorithms have dominated in regattas and tourist applications. These algorithms are executable on low power computers but have a number of functional limitations. Namely, they do not bypass the complex coastline - and if they have workarounds for this problem, they lose the quality of the designated routes or stop working altogether.

In the near future, the use of algorithms for determining sea routes in the standard navigation activities of the captain of the yacht will be popularized. As was the case with Google Maps and land navigation, the mathematical basis of this process will be Dijkstra's algorithms. The Dijkstra algorithm, however, must be adapted to the issues of maritime navigation. Such an algorithm is VEO Weather Routing presented in this article. The logic of Dijkstra's procedure in the VEO system allows much better avoidance of land and zones excluded from water traffic, while gaining very high accuracy in terms of the size of the geographic square to which the ship is heading. This is at the expense of the computing power of computers, in this case high-performance graphics processing unit computers.

Theoretical basis

The main tool used in determining the routes is the so-called Dijkstra's algorithm [1]. It is used to determine the shortest path in a graph including all nodes - starting from the first node (source) or to determine the shortest path between two nodes.

Dijkstra's algorithm divides the set of graph vertices into three sets: V, Q, S. The first of them V is the set of all nodes in the graph. Q is the set of all vertices sharing an edge with the source node. The node from the set Q with the shortest path from the center goes to the set S, and then it transforms into the middle node, and the original source node is no longer included in the further operation of the algorithm. The algorithm ends when the set Q is empty or when it reaches a previously defined node.

Dijkstra's algorithm for determining the route between two nodes in the graph can be described on the following example - Figure 1.

Figure 1 Dijkstry’s algorithm procedure

The vertex of S is the source node and the vertex F is the ending node. Each edge in the graph is weighted (non-negative). The operation of the Dijkstra algorithm is to find the vertex connected to the source node by the lowest weight common edge, which in the figure is the connection of the source node S to node A. Then the node connected to the source node takes over its role (A). This iteration continues until the endpoint (S-> A-> C-> E-> F) is reached. The figure shows (with a bold line) the path with the lowest sum of weights (4 + 3 + 3 + 2 = 12). Their sum was 12 and represents the optimal route (minimum cost to get there) between nodes S and F.

The Dijkstra algorithm is commonly used in single vehicle routing programs. On the basis of this algorithm works, among others route mapping system in Google Maps.

Marine adaptation of Dijkstra's algorithm

In order to adapt the Dijkstra algorithm for marine navigation, it is necessary to consider the following aspects:

1. Provide a method of weighing sections of the route that are adequate to the specificity of maritime navigation, thus:

a. Take into account that some edges will be faster than others due to the optimal matching of the sailing characteristics of the yacht and the wind that blows in a given geographic square at a given time

b. Take into account that reaching the boundary of the geographic square requires providing new weather data for the next period of time

2. Provide a map with areas where it is possible to swim safely without reaching the land / shallow water

3. Provide a map with areas where it is not safe to swim:

a. sailing is risky (e.g. shallow water, forecasted storm in a given time window)

b. sailing is uncomfortable (high wave, strong swing)

c. sailing is forbidden (fishing grounds, military zones, etc.)

4. Take into account the effect of sea currents coinciding with certain edges of the voyage and adding to the speed of the yacht

5. Take into account the drift of the yacht, so the need to set a different sailing angle to correct the wind and wave drift of the yacht

6. Provide such detailed map squaring that it is possible to precisely avoid small water obstacles - eg objects 100 meters wide.

7. Provide a routing procedure based on an interactive approach - first rough topology of edge points, then correction for detailed map and weather conditions.

The inclusion of some of these guidelines is described in the Implementation section. All of them will be fully implemented in the next development of the algorithm.

Implementation

The solution was implemented as a set of tables and stored procedures of the PostgreSQL database with the PostGIS spatial extension and the PG_Routing routing machine.

For efficiency reasons, the operation of the algorithm has been divided into two phases:

• In the preliminary phase, all global and time-stable data are prepared, performing more time-consuming operations, such as calculating the land content in the mesh, preparing the polar plots of yachts, determining the edges for initial routing or determining the vertices of the topological graph used in later steps to determine the final routes

• In the proper phase, only the operations necessary to calculate the route are performed and only on the part of the space in which the route will be mapped. This allows for a significant acceleration of the algorithm and does not generate unnecessary operations on parts of the space that will not be visited.

During the server operation, cyclic operations such as downloading weather data and updating them in the database are also performed.

The calculation process consists of two phases: the environment preparation phase and the route calculation phase.

The environment preparation phase

Land tables preparation

Initial table with lands is created from Natural Earth [3] data, which is imported into the database using QGIS [7].

Preparation of a table with excluded areas

Excluded areas are created on the basis of data from the OpenStreetMap project [8] by the community gathered around the OpenSeaMap project [9]. The data extract is downloaded and filtered using the osmosis application, and the resulting file is imported into the database with the osm2pgsql processor. The imported data is rewritten to the tables compliant with the solution scheme.

Preparation of a table with the speed of the yachts in relation to the wind

Polar plot is a chart linking the sailing characteristics of the yacht with the strength and direction of the wind that affects it, allowing you to read what speed the yacht will reach in a given direction and strength of the wind.

An example of a polar plot is presented below:

Figure 2 Polar plot example [4]

The table with the speed of yachts in relation to the wind is created on the basis of the data available in the repository [5]. As part of the project, a dedicated tool was created that allows you to download this data and import it to the database in a structure consistent with the solution diagram.

Preparation of a table with weather data

For the needs of the solution, a service has been created that periodically downloads current weather data in GRIB format from NOAA GFS [6] sources, extracts the values necessary for the correct operation of the algorithm, and then imports them into the database.

Creation of a table with 0.5 degree squares

Based on the conducted experiments, we determined that the compromised maximum size of a tile is 0.5 degrees (in the future, we will reduce this tile at least 10 times as we increase the processing power). This size is used to create the initial subdivision grid.

Figure 3 0.5 degrees tiles example

Inserting the land fragments intersecting with the square into the table

For each square of the mesh generated in the previous step, a piece of land overlapping with it is separated. The geometry of the land in the later steps is used to determine the land content in the mesh, and when dividing the mesh further, it is used as the land source for the smaller mesh sizes. This approach significantly reduced the time needed to generate the initial data.

Determination of the land content in a square

Based on the comparison of the mesh area of the mesh and the area of the land it contains and the excluded areas, the level of the mesh filling with land is determined in three values:

• f - full - these stitches will be ignored during further processing

• p - partial - these meshes will still be divided in order to isolate the landless parts from them

• n - none - these stitches are considered regular and are not further divided

Figure 4 Identification of tiles with land, with land and water, with only water

Recursive dynamic division of squares with partial land content

The next step is to extract from the meshes, partially land-containing, non-land fragments. Thanks to this, we separate fragments of space through which the yacht can safely move. All the meshes of the semi-land mesh and the exclusion areas are subdivided, the full land meshes are discarded and the landless meshes are left valid. The process takes place in a recursive loop until an experimentally determined compromise size of the minimum mesh size is reached. This step involves a compromise between the initial resource generation time, the route generation time, and the accuracy of the route and the possibility of determining it in areas densely covered with islands or in narrow straits. In subsequent phases of application development, this tradeoff may be shifted towards accuracy as higher performance is achieved in the routing stage.

Figure 5 Recursive separation of squares with land and water into smaller elements

Supplementing the relation to the table with weather data

In the next step, based on the spatial relation, the relation of the foreign key to the table containing weather data is created. The mesh identifier of the weather table is rewritten to all derived elements (edges of the topological graph), which prevents multiple spatial joins and significantly affects the efficiency of the final routing procedure.

Figure 6 Adding a link with weather data to the tiles

Creation of a table with rough edges

In the next step, the topological graph used in the first step of the routing procedure to determine the shortest route, disregarding the weather conditions, is determined between the centroids of the mesh that do not contain land.

Creation of the main mesh

In the next step, on the basis of the coarse mesh, the main mesh is created by dividing the multilines that make up the mesh mesh into individual lines, which in the next steps will be divided and used to create a topological graph used to calculate the target route taking into account weather conditions.

Creation of topology vertices

As the last step of the data preparation phase, points are created on the created mesh, which are used in the next steps of mapping the route as vertices of the topological graph after which the routing will take place. Points are created by dividing each edge of the mesh into 5 parts (experimentally selected parameter).

Figure 7 Topology of the edge points of the tiles

Route calculation phase

Finding the shortest path along a rough topology

The first step in mapping out a route is to find the mesh meshes that intersect with the start and end points, and calculate a rough route. The conducted experiments allowed to determine that the increase in the length of the route due to weather conditions is profitable only in the range of 30%, therefore the final route must directly correspond to the shortest route. Only the edge length was used as the cost of traveling the edge, and the route is determined using a function implementing the Dijkstra algorithm in the PG_Routing extension.

Figure 8 A simplified route along the rough topology of edge points

Creation a table with edges

At the present stage, we already know which way roughly leads the shortest route between the given points. In the next step, a dense topological graph is determined on the area of the mesh meshes covered by the shortest route, after which we will determine the target route, taking into account weather conditions. The vertices of the topology on the meshes intersecting the coarse route, created during the data preparation phase, are joined by the edges.

Figure 9 Detailed topology of marginal points related to the rough route

Attaching the start and end points of the route to the topology

As the previous calculations were made on the basis of the centroids of the mesh meshes traveled, it is necessary to connect the correct start and end points of the route to the topological graph to perform the target route search. For this purpose, the edges of the mesh meshes intersecting with the start and end points are removed, and new edges are created, which start at the route start point and end at the previously prepared topology vertices.

Extending the graph with adjacent tiles

In order to be able to determine the fastest route corresponding to the coarse route, the topology with tiles adjacent must be extended to those on which this route runs. A routine that extends an existing edge set is prepared for recursive operation and can extend the topology multiple times. This point is another compromise point - each expansion of the graph extends the time of calculating the final route, but at the same time increases the probability of finding a much more effective route at a greater distance from the coarse route. In the next steps of the algorithm development, this step can also be made dependent on the calculation result - if the destination route is not found, you can return to it, extend the topology and try to calculate the route again.

Figure 10 Detailed topology of marginal points related to the rough route

Adding the current and wind information to the edge

The information about currents and winds is added to the designated edges based on the downloaded weather forecasts. The weather information is shifted in time based on the estimated time of reaching each edge. In the next step, the azimuth of each edge and its angles in relation to the wind are calculated both for traveling in and against its direction - these parameters will be necessary to determine the time to pass the edge.

Calculating the time it takes for each edge to travel in both directions

The last operation necessary to prepare the data for the route calculation is the calculation of the time to travel each edge in both directions. Edge crossing times are interpolated on the basis of a table of average yacht speeds for a given yacht type, wind force and direction, which was imported into the system during the data preparation phase.

In the presented picture, the color of the edge reflects the time in which the edge can be traversed.

Figure 11 Generation of different routes and times marked with color

At this stage, it is best to notice the differences with the isochronous algorithm popular on the market. The routes determined by the Dijkstra algorithm are located in the tiles associated with the rough route, on the edge points of high density, while in isochronous algorithms these routes are located in funnels generated at a certain angle of inclusion, without being related to geographic tiles and their edge points. An example of an isochronous funnel in the MAX SEA TIMEZERO program:

Figure 12 Isochronous algorithm for a similar route in the MAX SEA TIMEZERO program  

Search for the fastest route along the edges

The data preparation process culminates in mapping the target route on the prepared edges. To calculate the route - as in the case of the coarse route, a procedure implementing the Djikstra algorithm from the PG_Routing extension is used. In this case, the cost of traversing the edge is the time it takes to travel it, calculated in the previous step.

The picture below shows the calculated roughing route against the edge with the calculated time of their travel.

Figure 13 Determining the optimal route (light green)

The comparison of the coarse route and the destination route against the weather data is as follows:

Figure 14 The coarse route and the optimized route (light green) calculated automatically

Examples of calculations

Below are 3 calculations for exemplary routes, for the Salona 44 yacht, for an exemplary weather for a given date.

The "Simplified cruise time" parameter represents approximately the cruise time for a manually marked route, based on manually designated waypoints, converted to the typical speed of a tourist yacht of this class. The parameter "Time of the cruise on the optimized route" shows the effect of the VEO Weather Routing algorithm for map quadratisation of the order of 0.5 degrees (including recursive correction for a given map area). Ultimately, the squaring of maps will reach 0.005 degrees with at least 20 times faster calculations thanks to the use of GPU processors with very high processing parameters. There is further potential to reduce the squares of the map being analyzed and further speed up the algorithms (increasing the computing power of the host computer).

The parameter "Algorithm calculation time in the medium PC class environment" should be treated as a comparative information, giving only an idea of ​​the typical computing load from the perspective of average personal computers.

Sea route from-to

Cruise time on simplified route (coarse topology)

Cruise time on an optimized route

(detailed topology)

Time of calculating the algorithm in the medium PC class environment

Gdynia (marina) - Łeba (port) – ok. 80 mil

16,5h

14,3h

12,8s

Ystad (Sweden) - Rosklide (Denmark) – ok. 150 mil

31,2h

26,1h

15,4s

Krk (town, KRK Island Croatia) - Gujak, Kornat (Croatia) – ok. 80 mil

16,8h

13,4h

12,8s

Table 1: Comparison of results gathered with the algorithm usage

Podsumowanie

The article demonstrates the usefulness of the Dijkstra algorithm in the marine adaptation built by the authors called VEO Weather Routing. This mechanism is used to solve the problem of determining the sea route in terms of weather data, navigable characteristics of the yacht - with automatic avoidance of land (and zones excluded from navigation) with high efficiency and precision. The full calculation procedure in both phases of the algorithm's operation was presented. As a result, a significant time gain was proved on the route marked with the use of the constructed mechanism in relation to the route marked by hand.

Bibliography 

1. Dijkstra E. W., A note on two problems in connexion with graphs, No. 1 (1), pp. 269-271, Numerische mathematik, 1959.

2. Tan Y., Wu D., Research on Optimization of Distribution Routes for Fresh Agricultural Products Based on Dijkstra Algorithm, No. 336, pp. 2500-2503, Applied Mechanics and Materials, 2013.

3. Website: https://www.naturalearthdata.com/

4. Wikipedia website: https://en.wikipedia.org/wiki/Polar_diagram_(sailing)

5. Website http://jieter.github.io/orc-data/site/

6. Website https://www.ncei.noaa.gov/products/weather-climate-models/global-forecast

7. Website and QGIS system: https://www.qgis.org/pl/site/

8. Website of the OpenStreetMaps project: https://www.openstreetmap.org/

9. The website of the OpenSeaMaps project: https://www.openseamap.org/index.php?id=openseamap&no_cache=1

Summary

Title

VEO Weather Routing - an innovative algorithm for calculating routes based on weather, map and navigable data of yachts

Abstract (j. ang)

The article demonstrates the usefulness of the Dijkstra algorithm in the marine adaptation built by the authors called VEO Weather Routing. This mechanism is used to solve the problem of determining the sea route in terms of weather data, navigable characteristics of the yacht - with automatic avoidance of land with high efficiency and precision. The full calculation procedure in both phases of the algorithm's operation was presented. As a result, a significant time gain was proved on the route marked with the use of the constructed mechanism in relation to the route marked by hand.

Keywords: Weather routing, sea navigation, yacht navigation