Ocean Engineering 153 (2018) 215–224

Contents lists available at ScienceDirect

Ocean Engineering
journal homepage: www.elsevier.com/locate/oceaneng

Ship voyage optimization for safe and energy-efficient navigation:

A dynamic programming approach
R. Zaccone a, *, E. Ottaviani b, M. Figari a, M. Altosole a
a
DITEN - Department of Electrical, Electronic, Telecommunications Engineering and Naval Architecture, Polytechnic School, University of Genoa, Genoa, Italy
b
OnAIR s.r.l., Genoa, Italy

A R T I C L E I N F O A B S T R A C T

Keywords: The paper presents a 3D dynamic programming based ship voyage optimization method, aiming to select the
Ship voyage optimization optimal path and speed profile for a ship voyage on the basis of weather forecast maps. The optimization is
Dynamic programming performed in accordance to a minimum fuel consumption strategy taking also into account ship motions and
Weather routing comfort. The optimization is carried out in a discretized space-time domain: the ship voyage is parametrized as a
Ship propulsion multi-stage decision process in order to formulate a dynamic programming optimization problem. Waves and
Ship motions
wind conditions are estimated for each route segment by weather forecasting maps then seakeeping related in-
dexes and fuel oil consumption are computed taking into account wave-induced ship motions and added resis-
tance. The best routing solution is thus selected by a dynamic programming algorithm developed and
implemented by the authors. Results and discussion of the proposed method are presented for a merchant ship
application in a test case voyage through the Northern Atlantic Ocean and compared to the constant speed great
circle solution.

1. Introduction rather than on external conditions (Chen, 2013), in order to better fit

different ship types, shapes and dimensions: different ships have different
In the recent years the continuously increasing availability of reliable responses in the same weather and speed conditions. Ship and human life
weather forecast data has significantly improved the safety of the ship safety, fuel consumption, energy efficiency, crew and passengers com-
voyages, helping the operators to select proper routes to avoid rough fort, voyage time management, control of delays are possible tasks which
weather and have an estimate of the time of arrival (ETA) and the voyage can be pursued by voyage optimization.
cost. Moreover, increased attention is put nowadays to seakeeping abil- Weather routing services and codes available on the market are
ities of ships in order to improve passage safety even with rough sea. usually not supported by public domain scientific papers due to confi-
However, most of the time medium intensity weather conditions are dentiality reasons. In author's knowledge very often they are based on the
encountered; these conditions do not affect ship safety but influence fuel principle of ‘storm avoidance’. The typical approach is to implement a set
consumption and comfort on board. In this framework optimization al- of generic speed reduction curves in function of sea state parameters (for
gorithms can provide a significant support to the decision making process example the Beaufort Wind Force Scale) and a number of global ship
in order to select the best choice in sight of one or more objectives. parameters (for example, displacement and length). The fuel consump-
The selection of the optimal route combines a number of objective tion is estimated on the basis of empirical relationships, while ship hull
functions as well as various constraints. In principle a ship voyage is geometry, seakeeping abilities and propulsion system features are
characterized by a starting point, an arrival point, a constrained arrival neglected. This approach can lead to unnecessary diversions to avoid
time window, eventually a number of fixed way points. Geographical rough weather, badly affecting the fuel consumption. The ship speed
(static) constraints need to be considered as well: for example, the shore decrease in rough weather can occur as a consequence of two reasons:
line, traffic separation schemes, restricted areas, bathymetry. Realistic involuntary speed reduction due to the additional resistance induced by
modeling of the ship behaviour in relation to weather conditions is wind and waves, and/or voluntary speed reduction to avoid navigation
crucial to correctly estimate and compare the ship performance in hazards and excessive ship motions which would result in propeller
different conditions. Decision making needs to be based on ship response, racing, slamming or green water. Additional constraints related to safety

* Corresponding author.
E-mail address: raphael.zaccone@edu.unige.it (R. Zaccone).

https://doi.org/10.1016/j.oceaneng.2018.01.100
Received 16 March 2017; Received in revised form 19 October 2017; Accepted 25 January 2018
Available online 4 February 2018
0029-8018/© 2018 Elsevier Ltd. All rights reserved.
R. Zaccone et al. Ocean Engineering 153 (2018) 215–224

and/or comfort might thus be considered to limit ship motions or wave computational efficiency. The implemented algorithm performs an effi-
induced forces. A proper ship response modeling is crucial in order to cient systematic exploration of the domain of the solutions, defined by a
correctly estimate these phenomena on a case by case basis. Furthermore three dimensional space time grid. The result of the optimization is a
a model of the actual propulsion system is essential in order to predict the sequence of waypoints and intermediate times determining the ship
propulsion performances in rough sea and to take into account the pro- trajectory and speed profile. With respect to present state of the art
pulsion system limits when steaming into stormy conditions. Moreover, weather routing and voyage optimization methods, the presented
ship speed changes need to be managed by the optimization algorithm in approach presents a twofold benefit: it is based on a detailed description
order to use the speed adjustment as a rough weather avoidance of the propulsion system allowing a realistic evaluation of the fuel con-
parameter in addition to course deviation but within the propulsion sumption and exhaust emissions versus ship speed and it implements an
system thrust/power capability. innovative 3D dynamic programming optimization routine which allows
The weather routing problem has been addressed by many authors in to manage both ship route and speed profile.
the past and different approaches have been proposed. The first pioneer To this end, the main tasks to be solved, in order to manage the route
works were centred on finding the minimum time of arrival on a voyage optimization problem, are:
(James, 1957; Zoppoli, 1972; Papadakis and Perakis, 1990) however
most of the authors neglected the ship response behaviour in rough sea.
weather forecast data availability and management;
More recently, voyage optimization has been approached in a wider
ship propulsion performance and hydrodynamic response modeling;
sense, taking into account ship motions and/or fuel consumption. The
optimization problem formulation and solution.
most used techniques include multi objective genetic algorithms (Marie
et al., 2009; Maki et al., 2011; Vettor and Guedes Soares, 2016; Zaccone A short overview is given to the first point, which is considered as
et al., 2016), deterministic enumerative algorithms (Lin et al., 2013; Fang input data. The paper is centred on the second and third points which are
and Lin, 2015; Shao et al., 2012), brute force optimization (Lu et al., going to be described in detail. Finally, a case study is analysed and re-
2015), local search algorithms (Safaei et al., 2015). Simulation of ship sults are shown in order to highlight the features of the presented
behaviour in rough sea have been proposed in the past by several au- approach in a realistic case.
thors: Journ
ee (1976); Journ
ee and Meijers (1980) proposed some
guidelines to develop a ship model to take into account either propulsion 2. Weather conditions
and ship motions. More recently weather routing oriented propulsion
simulation, without optimization, have been published by Coraddu et al. Two weather actions are considered: waves and wind. Wind forecast
(2013). Crew and passengers comfort evaluation techniques for human data is provided in terms of wind speed components ðu; vÞ as a function of
bodies exposed to mechanical vibrations and accelerations are suggested geographical coordinates ðλ; φÞ and time t. Forecasted wave significant
in the ISO 2631-1: 1997, based on the methods proposed by O'Hanlon parameters are provided as well in function of the same variables, in
and McCauley (1973), and Lawther and Griffin (1988). The authors particular significant wave height H13 , wave period T1 , and direction θ.
gained extensive experience in steady state and transient modeling of The wave parameters are used to fit a parametric spectral formulation.
ship propulsion systems (Altosole et al., 2008, 2014, 2016, 2017; Martelli ITTC 84 JONSWAP spectral formulation with a cosine square spreading
et al., 2014a,b) mainly for control design purposes. function is used:
In the here presented work ship voyage optimization problem is
tackled by means of a problem-specific algorithm based on 3D dynamic H12

944 2
programming (3DDP) coupled with a dynamic ship propulsion model and Sðω; θÞ ¼ 155 3
exp  4 4 γ Y cos2 θ (1)
ω
5T4
1 ω T1 π
with weather forecast data. Fuel oil consumption is considered as the
objective function to minimize while ship motions and expected time of being:
arrival (ETA) are used as constraints. In particular, ship motions based


constraints are imposed on probability of slamming and deck wetness 0:191ωT1  1
Y ¼ exp  pffiffiffi (2)
(Journ
ee and Meijers, 1980), motions sickness index (MSI) (O'Hanlon 2σ
and McCauley, 1973), and lateral forces (Perez, 2006). The constraints
allow to implicitly take into account any voluntary speed reduction in where ω is the wave circular frequency, and γ ¼ 3:3 for Northern Atlantic
presence of rough weather, while involuntary speed reductions are applications. The choice of the spectral formulation can significantly
simulated by considering the effect of the ship resistance increase due to affect the results, as it directly influences the ship motions and resistance.
waves (Journ
ee, 1976) on the propulsion system. Wave and wind fore- In particular, a correct selection of parameter γ should be made in
cast data are considered in function of time and geographical co- function of the geographical area in which the optimization is performed.
ordinates: wave significant height, mean period, direction and wind Moreover, real and parametric spectra may differ significantly in pres-
speed components are obtained via space-time interpolation in the ence of crossing seas (Mentaschi et al., 2015; Spentza et al., 2017): these
forecast maps. The obtained data are used to estimate ship motions and aspects will be deeply investigated in the future research, while the
added resistance in order to assess the required engine power and fuel presented work is mainly focused on investigating the potential of the
consumption via steady state ship propulsion simulation. The fuel con- optimization procedures in a voyage optimization framework.
sumption is used as the cost function to evaluate each route segment and
compute the optimal solution. DP approach fits well to ‘best path search’ 3. Ship model
problems (Bellman, 1958), so is one of the classic choices to tackle ship
voyage optimization tasks. In particular, the problem is solved by The ship numerical model allows to estimate the fuel consumption in
exhaustively exploring a discretization of the search domain while very given geographical positions, weather conditions and ship speeds. The
little restrictions are put on either the objective function and optimiza- model is structured in four main blocks:
tion constraints because no derivatives are calculated. The deterministic
nature of the algorithm may be seen as an advantage with respect to
hydrodynamic resistance calculation
heuristic global search algorithms. Nevertheless, in a weather routing
ship motions and comfort assessment
framework, the benefit is partially limited by the fact that input data is
propeller performance prediction
affected by significant uncertainties. DP method has some drawbacks in
engine performance prediction
terms of compuatation speed. However, the presented implementation
includes proper problem-specific pruning strategies to boost The model arrangement is shown in Fig. 1. At any time and

216
R. Zaccone et al. Ocean Engineering 153 (2018) 215–224

Fig. 1. A scheme of the ship model arrangement.

 2
geographical coordinates the weather conditions (wave and wind data) Sηi ¼ RAOηi  Sζ (6)
are estimated by interpolating into weather forecast data then the hy-
drodynamic resistance as well as ship motions are computed. Finally the where ω is the wave circular frequency, μ is the encounter angle, Sζ ðω; μÞ
required propeller thrust and torque are evaluated using propeller open is the directional wave energy spectral density, and RAOηi is the complex
water diagrams in order to finally assess the engine power and fuel response amplitude operator of motion ηi , with i ¼ 1 to 6 represent surge,
consumption. sway, heave, roll, pitch and yaw motions respectively. The probability
that motion ηi exceeds the value ηi is given by the following relationship:
3.1. Resistance "  2#
 η
P ηi > ηi ¼ exp  i (7)
The hydrodynamic resistance of the ship, Rtot , is decomposed in 2m0ηi
accordance with the following equation:
where m0ηi is the spectral moment of order zero of motion ηi . The joint
Rtot ¼ Rt þ Raw þ Rwind (3) probability that motions ηi and ηj exceed the values ηi and ηj respec-
tively, is obtained under the hypothesis of independence between ηi and
where Rt is the still water hydrodynamic resistance, Raw is the added
ηj , and is given by the following formula:
resistance due to the rough sea effect, and Rwind is the wind resistance.
The added wave resistance Raw is estimated according to the following 8 2 2
39
>
<   2 >
=
expression (Lewis, 1989): 6 η η
j 7
P ηi > ηi ; ηj > ηj ¼ exp  4 i þ 5 (8)
>
: 2m0ηi 2m0ηj >;
∞ 2π
Raw ¼ ∫ 0 ∫ 0 Φaw Sζ d ωdμ (4)
The spectral moment of order n per squared wave height of motion ηi
where ω is the wave circular frequency, Sζ is the directional wave energy
is given by:
spectral density, Φaw ðω; θ; VÞ is the wave added resistance pseudo-
response amplitude operator (which expresses the longitudinal drift ∞ 2π
mn;1 ¼ ∫ 0 ∫ 0 ωn Sηi d ωdμ (9)
forces per squared wave amplitude and wave circular frequency) and μ is
the encounter angle. Lateral force estimator (LFE) is used as a global index to take into
The wind resistance Rwind is expressed through a function of the account either lateral accelerations and displacements (Perez, 2006). The
relative wind speed Vw;r : LFE at the point P of cartesian coordinates ðxP ; yP ; zP Þ is given by the
following equation:
1 2
Rwind ¼ ρair AF Vw;r cX ðγ w Þ (5)  
2 LFEP ¼ €η2  €η4 zP  zg þ €η6 xP  xg  gη4 (10)

where ρair is the air density, AF is the ship's above water front projection where €ηi represent the acceleration of motion ηi . The above mentioned
area and cX is the wind resistance coefficient depending on the encounter equation is an expression of side component of the apparent acceleration
angle γ w . The effects of rough weather are thus modelled in terms of acting on a body on a ship, expressed in its relative system of coordinates,
resistance increase rather than speed reductions: the resistance increase simplified for small amplitude motions. In particular, The first term
is actually the cause of the involuntary speed reduction or the power represents the acceleration induced by the sway motion, the second and
output increase to achieve the required speed (within engine limits). third express the contributes due to roll and yaw angular accelerations,
multiplied by the respective arms, while the last term is the lateral
3.2. Ship motions and comfort component of the gravity acceleration originated by the roll angle. A
more detailed description of the mathematical proof of equation (10) can
The energy spectrum of absolute motion ηi is evaluated in accordance be found in Perez (2006).
with the following equation (Lewis, 1989): Ship passengers comfort evaluation has been discussed by various

217
R. Zaccone et al. Ocean Engineering 153 (2018) 215–224

authors, which propose proper parameters depending mainly on vertical the gearbox efficiency respectively considered as constant values. The
accelerations mean value and frequency. The most used parameters are engine fuel mass flow rate m_ f is evaluated by using the engine perfor-
the Motion Sickness Incidence (MSI), proposed by O'Hanlon and mance map modelled using the engine MCR power and revolution. By the
McCauley (1973), and the Vomiting Incidence (VI), by Lawther and performance map each working conditions are checked for not exceeding
Griffin (1986, 1988). In this study, MSI has been used to estimate crew the engine limits (Fig. 2). The feasibility or infeasibility of any combi-
comfort. The MSI at point P for a reference exposition time of 2 h is given nation of speed and weather conditions as well as the associated cost in
by the following equation (Cepowski, 2012; Piscopo and Scamardella, terms of fuel consumption are thus a result of the propulsion simulation.
2015), in dependence of the spectral moments of the absolute vertical The procedure is summarized in Algorithm 1.
motion η3P :
Algorithm 1. Fuel mass flow rate computation procedure.


log10 agv  μMSI
MSI ¼ 100 0:5 þ erf (11)
0:4

where g is the gravitational acceleration, erfðxÞ is the error function, and

av and μMSI are defined as follows:
pffiffiffiffiffiffiffiffiffiffiffi
av ¼ 0:798 m4;η3P
 rffiffiffiffiffiffiffiffiffiffiffi 2 (12)
m4;η3P
μMSI ¼ 0:819 þ 2:32 log10
m2;η3P

3.3. Propulsion

The propeller performance simulation is carried out by using the open

water diagrams reporting thrust and torque coefficients, KT and KQ
respectively, versus advance coefficient J respectively defined as follows:

Va
J¼ (13)
nD

T
KT ¼ (14)
ρn2 D4

Q
KQ ¼ (15)
ρn2 D5

where ρ is the density of water, T and Q are the propeller open water
thrust and torque, Va is the propeller advance speed, n is the propeller
revolution speed and D is the propeller diameter. The propeller revolu-
tions are evaluated imposing the equilibrium between required vs
delivered thrust by using the auxiliary variable KJ 2T is used. The engine
power is computed in accordance to the following equation: The mass of fuel burned mf of each route segment is evaluated by
numerical integration of the fuel mass flow rate with respect to time as
2πρH2 O n3 D5 KQ shown in equation (17):
PB ¼ (16)
ηr ηsηg
t
mf ðXi ;ti Þ→ðXj ;tj Þ ¼ ∫ tji m_ f ðtÞdt (17)
where ηr , ηs, ηg are the relative rotative efficiency, the shaft efficiency and
where is supposed to ship travel from point Xi at time ti to point Xj at time
tj .
Expression 17 allows to compute the cost function at each DP step.
The cost of the whole voyage is obtained by summation of all the
contributes.

4. Route parametrization

A parametric definition of the ship voyage is needed in order to

perform the optimization. The proposed parametrization is based on
great circle geometry in order to easily identify the shortest solution: in
an optimal solution any deviation from the great circle is thus caused by
rough weather.
The minimum angular distance ΘAB , measured in radiants, between
two points A and B satisfies the following relationship:

cosΘAB ¼ sinφA sinφB þ

(18)
þcosφA cosφB cosðλB  λA Þ

Fig. 2. Feasibility check of propulsion system working points. where φA , φB , λA , λB , denote latitude and longitude of the points.

218
R. Zaccone et al. Ocean Engineering 153 (2018) 215–224

The course α is defined as the angle between the ship heading and the 5. Dynamic programming optimization
true north direction, it is measured clockwise and it changes continuously
along a segment of great circle. The following relationships provide the Dynamic programming is an optimization strategy that aims to solve a
initial and final course angles, respectively αi and αf , of the great circle problem by using a multi-stage approach in which the final result is the
segment linking points A and B: consequence of a number of separate decisions. Optimal solution is found
on the basis of Bellmans principle of optimality:
sinφB  sinφA cosΘAB
cosαi ¼ (19) “an optimal policy has the property that whatever the initial state and
cosφA sinΘAB
initial decision are, the remaining decisions must constitute an
sinφB cosΘAB  sinφA optimal policy with regard to the state resulting from the first deci-
cosαf ¼ (20) sion” (Bellman, 1954).
cosφB sinΘAB
The geographical coordinates of a point X laying on the great circle The structure of the method well adapts to any multi stage decision
segment linking A and B and the local course angle are given by the process and among all ‘best path search’ type problems (Bellman, 1958),
following equations: (Bellman, 1962), as they can be reduced as a sequence of decisions in
dependance of which a return is obtained.
sinφX ¼ sinφA þ cosðmAB xÞþ
(21)
þcosφA sinðΘAB xÞcosαi 5.1. Formalization

cosðΘAB xÞ The proposed formulation of the ship voyage optimization problem is

cosðλX  λA Þ ¼ 
cosφA cosφX (22) solved in a three dimensional domain: start and arrival points and time
þtanφA tanφX instants are assigned a-priori, then a two dimensional domain is identi-
fied by discretizing parameters x and δ in equations (21)–(25) (Fig. 4). A
sinφX cosðΘAB xÞ  sinφA discretization of the time axis (Fig. 5) is as well performed. The route is
cosαX ¼ (23) modelled as finite set of Nþ1 way points through which the ship passes at
cosφX sinðΘAB xÞ
assigned times. The ship is supposed to travel at constant speed between
where x ¼ ΘΘAX
AB
is an advance parameter. each pair of waypoints. The route r is thus given by the following sets:
The coordinates of points Y are given as a function of the deviation
X ¼ fX0 ; X1 ; …; XN g
parameter δ by the following expressions: r: (26)
t ¼ ft0 ; t1 ; …; tN g
φY ¼ φX  δsinαX (24)
where Xi ¼ ðλi ; φi Þ, X0 ¼ A, XN ¼ B, t0 ¼ 0, tN ¼ ETA.
Equation (17) allows to compute the fuel mf ðXi ;ti Þ→ðXj ;tj Þ consumed
λY ¼ λX þ σ δcosαX (25)
from point Xi at time ti to point Xj at time tj steaming at constant speed.
where σ ¼ 1 going eastward and σ ¼ 1 going westward. Fig. 3 sum- The objective of the problem is to minimize the total fuel consump-
marizes the above presented sets of equations. tion:
In the proposed dynamic programming framework the ship route is X
N 1
modelled as a finite set of way points fX0 ; X1 ; …; Xn g such that Xi ¼ mf ¼ mf ðXi ;ti Þ→ðXiþ1 ;tiþ1 Þ (27)
ðλi ; φi Þ and a corresponding set of arrival times at the way points, i¼0

ft0 ; t1 ; …; tn g, where t0 ¼ 0 and tn ¼ ETA. Each segment Xi →Xiþ1 , is sailed The problem can be written in the form of the Bellman's equation as
at constant speed Vi . follows. Let fi ðX; tÞ be the fuel consumed at stage i from the start,
The described scheme allows to identify unequivocally the great following an optimal policy, being the ship in position X at time t: the
circle i.e. minimum distance solution which is the most convenient in target of the search is thus fN ðXN ; ETAÞ, i.e. the mass of fuel consumed at
calm sea and obstacle-free conditions. the end of the voyage identified by the arrival point XN and estimated
time of arrival ETA, following an optimal policy. In accordance to Bell-
man's Principle of Optimality, the following functional relationship can
be written:

fi ðX; tÞ ¼ min mf ðξ;τÞ→ðX;tÞ þ fi1 ðξ; τÞ (28)
ξ;τ

Fig. 3. Deviation δ of point Y from the corresponding point X laying on great

circle segment linking A and B.
Fig. 4. A schematic representation of the geographical search grid.

219
R. Zaccone et al. Ocean Engineering 153 (2018) 215–224

proposed by Holtrop (1984). Wind resistance coefficient cX has been

computed according to the method proposed by Blendermann (1994).
Finally, open source strip theory seakeeping code PDSTRIP (Bertram et.
al, 2006) has been used to evaluate the pseudo response amplitude
operator Φaw and ship motions absolute RAOs. The optimization is per-
formed with reference to a voyage representing a North Atlantic crossing.

6.1. Grid setup

Two optimization modes are compared: a speed profile optimization,

in which no deviations from the great circle route are allowed (referred
as speed optimization, SO), and a speed and course optimization
(referred as voyage optimization, VO) in which both speed and course
adjustments are performed in order to achieve the optimum. Moreover, a
constant speed great circle solution is presented for comparison.
A 10 stages by 11 steps grid with 1∘ maximum deviation is
Fig. 5. A representation of the times of arrival at each grid node. considered: the resulting waypoints are about 300 miles apart. Times of
arrival are discretized every 60 min, while integration time step is set to
30 min. The computation is carried out considering a time of arrival
where ξ and τ represent the possible positions and times of arrival at stage window between 200 and 400 h from departure.
i  1.
6.2. Constraints
5.2. Computational aspects
Ship safety and comfort constraints are imposed as well, in particular,
DP is an exhaustive approach, i.e. it guarantees to identify the optimal slamming and deck wetness probabilities are checked (Journ
ee, 1976;
solution associated to the discretization of the problem. Nevertheless, Journ
ee and Meijers, 1980) in accordance with the following
this advantage over heuristic global search algorithms is partially limited inequalities:
by the fact that input data is affected by significant uncertainties. DP
presents some drawbacks in terms of compuatation speed: Bellman-Ford pffiffiffiffiffiffi
P η3r;0:9L > T; η_ 3r;0:9L > 0:093 gL < 0:03 (29)
implementation on a graph (Bellman, 1958) runs in OðjVjjEjÞ, where jVj
and jEj are the number of graph and edges respectively. In the present 
work, the breadth-first-search structure of the Bellman-Ford algorithm P η3r;B > fb < 0:07 (30)
has been maintained, however a significant speed-up has been made by

implementing problem-specific modifications to the algorithm. P €η3;B > 0:4g < 0:07 (31)
In the proposed solution scheme the cost of each route segment
linking each couple of nodes (all geographical positions and all time in- Being η3r;0:9L and η3r;B the relative vertical motion at 90% of the ship
stants) is evaluated in the forward phase. The search domain is pro- length and at bow respectively, L, T and fb the ship length, draft and free
gressively explored following a breadth-first approach. Nodes are sorted board at the bow respectively.
by priority to minimize the number of segments estimations: a compu- In addition, the following inequality constraints are imposed:
tationally cheap lower bound of the cost of each segment is first
LFERMS < LFEmax (32)
computed in order to early prune the less promising solutions. The
compliance of each node with the time of arrival window is checked as
MSI < MSImax (33)
early as possible. Nodes visibility based on geographical constraints is
checked by using an efficient nearest neighbour search algorithm based Where MSImax has been set to 10%, LFERMS is the root mean square of the
on kd-tree structures (Friedman et al., 1977) before evaluating the cost of lateral force estimator, and LFEmax has been selected in accordance with
each segment. The cost associated to each segment is computed by nu- (Perez, 2006), who proposes 0:1g as an upper bound for general
merical integration of the ship model output. purposes.
At the end of the calculation a number of nodes, characterized by a
time of arrival (i.e the ETA), and a cost associated (fuel consumption) are
6.3. Weather forecasts
identified. The final nodes identify global optima in terms of fuel con-
sumption at fixed time or the fastest solutions at fixed fuel consumption.
The forecasts used for the test case are based on the National Weather
The final set of nodes can thus be seen as the Pareto Frontier of a two
Service (NOAA) Wave Watch III model which provides forecast maps
objective optimization problem which aims to minimize both time of
every 6 h. The files are available in GRIB format, they are 162 h deep,
arrival and fuel consumption. Each node contains the information of the
they have a 1:25∘  1∘ . grid resolution and a 6 h time resolution (NOAA
previous one, so each solution is extracted by back-tracking.
MMAB Operational Wave Model). The weather forecasts have been
The optimization routine has been developed and implemented by
downloaded from: http://www.globalmarinenet.com.
the authors in Cþþ language, significantly relying on Cþþ Standard
Template Library. Armadillo linear algebra library (Sanderson, 2010) has
Table 1
been used for basic matrix manipulation.
Main characteristics of the case study ship.

Ship type Bulk Carrier

6. Case study
Length b. p. 182 m
Beam 31 m
The above described optimization code has been tested on a bulk Draft 9m
carrier to evaluate the effectiveness of the adopted optimization strategy. Displacement 41577 t
The main data of the considered vessel are summarized in Table 1. The Propeller 1  FPP, 4 blades, Diameter ¼ 6 m
Main engine 1  10 MW, 113 RPM 2 stroke diesel engine
still water resistance has been calculated in accordance with the method

220
R. Zaccone et al. Ocean Engineering 153 (2018) 215–224

In particular, the case study optimization is carried out in a typical

Northern Atlantic winter rough weather condition, using the forecast
maps of day 2016/01/21, 12:00 a.m. An example screenshot of the storm
occurring is shown in Fig. 6.

6.4. Results

Fig. 7 shows the fuel consumption vs. ETA trade off curve. Each point
of the curve is the minimum fuel consumption solution at the corre-
sponding ETA, i.e. the curve can be seen as the Pareto Frontier of the
minimum ETA and minimum fuel consumption optimization problem.
VO and SO solutions are compared and the results achievable by sailing
at constant speed are reported for comparison. Note that the course and
speed optimization curve (VO) has a significant lower consumption with
respect to the other two, however also speed optimization (SO) provides
good results especially for high ETAs.
Note in addition that both the optimization strategies allow to reduce
the negative effect of involuntary speed reduction on the ETA, allowing
thus to achieve lower times of arrival by adjusting speed and course to
the weather conditions, i.e. catching up the lost time when the weather
conditions are better. Moreover, each of the suggested solutions is the Fig. 6. An example map of the weather forecast used as a case study.
most fuel saving strategy to perform the voyage in the assigned time.
Finally, note that Fig. 7 allows to easily estimate the cost of each hour of
anticipation or delay in terms of tons of fuel consumed.
The trajectories associated to different ETAs (220, 250 and 300 h) are
presented in Figs. 8–10. The corresponding speed profiles are presented
in Figs. 11–13. It is worth noting that the optimization algorithm VO acts
drastically on ship speed to search the optimum, finding articulated
speed profiles.
The route and/or speed adjustments allow a global reduction of the
fuel mass flow rate required by the propulsion system (Figs. 14–16)
resulting in a reduction of the total energy requirement and CO2 emis-
sions. The fuel consumption profiles are presented in Figs. 17–19.
The solutions presented until now are referred to a westward voyage:
this direction is associated to head or bow sea mainly. If the opposite
voyage in the same day of departure is considered, the solutions change
radically: the fuel consumption to ETA curves of the two voyages are
compared in Fig. 20. The eastward voyage presents significantly lower
fuel consumption as expected: the following sea condition is predomi-
nant, thus the added resistance is lower. Moreover, lower ETAs are
achievable as the weather conditions are less severe: involuntary speed
reduction is thus less significant. In addition, the eastward optimized
solutions do not present significant fuel saving with respect to the con- Fig. 7. Fuel consumption to ETA trade off for the considered case
stant speed solution, if compared to the westward ones due to the fact study voyage.
that the course and speed changes have small effect on the hydrodynamic
resistance in following seas. Fig. 21 presents the trajectories in the two
cases, associated to an ETA of 250 h: note that the suggested routes are
significantly different, due to the different the conditions encountered.
Note finally that sailing in following sea may lead to significant lateral
motions. The constraint on lateral forces allows to limit these effects,
however a check of parametric roll (Maki et al., 2011) might be included
in a future version.

6.5. Computational time

The above presented computations have been carried out on a stan-

dard laptop which main data are summarized in Table 2. The computa-
tion time is significantly influenced by the chosen discretization
parameters and by the set up of the constraints. The average computation Fig. 8. Optimal routes for ETA ¼ 220 h.
times achieved are presented in Table 3. The computation times required
in case a unique ETA with no delay is selected are also reported: note that 7. Conclusions
the computation in this case is significantly faster. On the contrary the
computation time increases as the ETA window gets wider. In general, A ship voyage optimization method has been developed and pre-
the more the problem is constrained, the less nodes are processed, sented, aiming to find the minimum fuel consumption voyage, subject to
resulting in a faster computation.

221
R. Zaccone et al. Ocean Engineering 153 (2018) 215–224

Fig. 9. Optimal routes for ETA ¼ 250 h.

Fig. 12. Optimal speed profiles for ETA ¼ 250 h.

Fig. 10. Optimal routes for ETA ¼ 300 h.

Fig. 13. Optimal speed profiles for ETA ¼ 300 h.

Fig. 11. Optimal speed profiles for ETA ¼ 220 h.

safety and comfort constraints, by using 3D Dynamic Programming

optimization. A steady state ship model has been described. The model
computes the seakeeping and propulsion performances of the ship in
assigned wind and wave conditions. A parametric route model, based on
great circle navigation, has been introduced, in order to give a dynamic
programming formulation of the problem. The solution is found by
solving the Bellman's equation which cost function, i.e. the total mass of
burned fuel, is evaluated by means of a ship steady sate simulation
model. In order to take into account ship safety and crew comfort in the
optimization, slamming, deck wetness and high bow acceleration prob-
abilities, as well as motion sickness incidence and lateral forces in the Fig. 14. Fuel mass flow rate profiles for ETA ¼ 220 h.
crew accommodations are checked not to exceed threshold values. An

222
R. Zaccone et al. Ocean Engineering 153 (2018) 215–224

Fig. 15. Fuel mass flow rate profiles for ETA ¼ 250 h. Fig. 18. Fuel consumption profiles for ETA ¼ 250 h.

Fig. 16. Fuel mass flow rate profiles for ETA ¼ 300 h.
Fig. 19. Fuel consumption profiles for ETA ¼ 300 h.

Fig. 17. Fuel consumption profiles for ETA ¼ 220 h. Fig. 20. Fuel consumption to ETA trade off of eastward vs. westward voyages.

223
R. Zaccone et al. Ocean Engineering 153 (2018) 215–224

Bellman, R., 1962. Dynamic programming treatment of the travelling salesman problem.
J. ACM 9 (1), 61–63.
Bertram, V., Veelo, B, Soding, H., 2006. Program PDSTRIP: Public Domain Strip Method.
Hamburg, unpublished.
Blendermann, W., 1994. Parameter identification of wind loads on ships. J. Wind Eng.
Ind. Aerod. 51 (3), 339–351.
Cepowski, T., 2012. The Prediction of the Motion Sickness Incidence Index at the Initial
Design Stage. Zeszyty Naukowe/Akademia Morska w Szczecinie, pp. 45–48.
Chen, H., 2013. Voyage optimization supersedes weather routing. Jeppesen Commerc.
Mar. J.
Coraddu, A., Figari, M., Savio, S., Villa, D., Orlandi, A., 2013. Integration of seakeeping
and powering computational techniques with meteo-marine forecasting data for in-
service ship energy assessment. In: Developments in Maritime Transportation and
Exploitation of Sea Resources: IMAM 2013, p. 93.
Fang, M.C., Lin, Y.H., 2015. The optimization of ship weather-routing algorithm based on
the composite influence of multi-dynamic elements (ii): optimized routings. Appl.
Ocean Res. 50, 130–140.
Fig. 21. Trajectories of eastward vs. westward voyages for ETA ¼ 250 h. Friedman, J.H., Bentley, J.L., Finkel, R.A., 1977. An algorithm for finding best matches in
logarithmic expected time. ACM Trans. Math Software (TOMS) 3 (3), 209–226.
Holtrop, J., 1984. Statistical re-analysis of resistance and propulsion data. Int. Shipbuild.
Prog. 31 (363), 272–276.
Table 2 ISO 2631-1: 1997, 1997. Mechanical Vibration and Shock: Evaluation of Human Exposure
Main data of the PC used for computation. to Whole-body Vibration. Part 1, General Requirements: International Standard ISO
Processor Intel Core i7-6500U 2631-1: 1997 (E). ISO.
James, R.W., 1957. Application of Wave Forecasts to Marine Navigation.
Memory 16 GB
Journ
ee, J.M.J., 1976. Prediction of Speed and Behaviour of a Ship in a Sea-way. Delft
OS Ubuntu 16.04 64-bit
University of Technology.
C/Cþþ Compiler GCC Journ
ee, J.M.J., Meijers, J.H.C., 1980. Ship Routeing for Optimum Performance. Delft
University of Technology.
Lawther, A., Griffin, M.J., 1986. The motion of a ship at sea and the consequent motion
sickness amongst passengers. Ergonomics 29 (4), 535–552.
Table 3 Lawther, A., Griffin, M.J., 1988. A survey of the occurrence of motion sickness amongst
Time required by computations. passengers at sea. Aviat Space Environ. Med. 59 (5), 399–406.
Lewis, E.V., 1989. Principles of Naval Architecture. second revision. In: Motions in Waves
Computation Time Num. sol. and Controllability, iii. The Society of Naval Architects and Marine Engineers, Jersey
Speed opt. (no delay) <100 1 City, USA.
Lin, Y.H., Fang, M.C., Yeung, R.W., 2013. The optimization of ship weather-routing
Voy. opt. (no delay) 10 2000 1
algorithm based on the composite influence of multi-dynamic elements. Appl. Ocean
Speed opt. (Pareto) 600 134
Res. 43, 184–194.
Voy. opt. (Pareto) 100 5000 134
Lu, R., Turan, O., Boulougouris, E., Banks, C., Incecik, A., 2015. A semi-empirical ship
operational performance prediction model for voyage optimization towards energy
efficient shipping. Ocean Eng. 110, 18–28.
algorithm has been developed and implemented by the authors in Cþþ Maki, A., Akimoto, Y., Nagata, Y., Kobayashi, S., Kobayashi, E., Shiotani, S., Ohsawa, T.,
language to solve the proposed task. The presented approach has been Umeda, N., 2011. A new weather-routing system that accounts for ship stability based
on a real-coded genetic algorithm. J. Mar. Sci. Technol. 16 (3), 311–322.
applied to a case study voyage from France to New York along the Marie, S., Courteille, E., et al., 2009. Multi-objective optimization of motor vessel route.
Atlantic Ocean: the results associated to a constant speed voyage, a speed Mar. Navig. Saf. Sea Transport. 411.
profile only optimization, and a full speed and course optimization have Martelli, M., Figari, M., Altosole, M., Vignolo, S., 2014a. Controllable pitch propeller
actuating mechanism, modelling and simulation. Proc. Inst. Mech. Eng. M J. Eng.
been compared: the presented results include the voyage time to fuel Marit. Environ. 228 (1), 29–43.
consumption trade-off curves, the optimal trajectories and speed profiles Martelli, M., Viviani, M., Altosole, M., Figari, M., Vignolo, S., 2014b. Numerical modelling
obtained for different voyage times. Moreover, fuel mass flow rate and of propulsion, control and ship motions in 6 degrees of freedom. Proc. Inst. Mech.
Eng. M J. Eng. Marit. Environ. 228 (4), 373–397.
fuel consumption profiles associated to the considered solutions have
Mentaschi, L., Besio, G., Cassola, F., Mazzino, A., 2015. Performance evaluation of
been discussed. The sensitivity of the optimization procedure to different wavewatch iii in the mediterranean sea. Ocean Model. 90, 82–94.
weather conditions has been tested by optimizing the reverse voyage in O'Hanlon, J.F., McCauley, M.E., 1973. Motion Sickness Incidence as a Function of the
Frequency and Acceleration of Vertical Sinusoidal Motion, Technical Report, DTIC
the same days. In light of the obtained results a notable fuel saving po-
Document.
tential is observed, in particular in the voyage optimization (VO) case, in Papadakis, N.A., Perakis, A.N., 1990. Deterministic minimal time vessel routing. Oper.
which the algorithm takes full advantage from the ship propulsion model Res. 38 (3), 426–438.
features. Perez, T., 2006. Ship Motion Control: Course Keeping and Roll Stabilisation Using Rudder
and Fins. Springer Science & Business Media.
Piscopo, V., Scamardella, A., 2015. The overall motion sickness incidence applied to
References catamarans. Int. J. Naval Architect. Ocean Eng. 7 (4), 655–669.
Safaei, A., Ghassemi, H., Ghiasi, M., 2015. Voyage Optimization for a Very Large Crude
Altosole, M., Benvenuto, G., Campora, U., Laviola, M., Zaccone, R., 2017. Simulation and Carrier Oil Tanker: a Regional Voyage Case Study. Zeszyty Naukowe/Akademia
performance comparison between diesel and natural gas engines for marine Morska w Szczecinie.
applications. Proc. Inst. Mech. Eng. M J. Eng. Marit. Environ. 231 (2), 690–704. Sanderson, C., 2010. Armadillo: an Open Source Cþþ Linear Algebra Library for Fast
Altosole, M., Benvenuto, G., Figari, M., Campora, U., Bagnasco, A., D'Arco, S., Prototyping and Computationally Intensive Experiments.
Giuliano, M., Giuffra, V., Spadoni, A., Zanichelli, A., Michetti, S., Ratto, M., 2008. Shao, W., Zhou, P., Thong, S.K., 2012. Development of a novel forward dynamic
Real time simulation of the propulsion plant dynamic behaviour of the aircraft carrier programming method for weather routing. J. Mar. Sci. Technol. 17 (2), 239–251.
“cavour”. In: Proceedings of the Institute of Marine Engineering, Science and Spentza, E., Besio, G., Mazzino, A., Gaggero, T., Villa, D., 2017. A ship weather-routing
Technology - INEC 2008: Embracing the Future. tool for route evaluation and selection: influence of the wave spectrum. In:
Altosole, M., Figari, M., Ferrero, C., Giuffra, V., Piva, L., 2014. Propulsion retrofitting of Proceedings of the 17th International Congress of the International Maritime
the tall ship amerigo vespucci: automation design by simulation. In: 2014 Association of the Mediterranean, IMAM 2017.
International Symposium on Power Electronics, Electrical Drives, Automation and Vettor, R., Guedes Soares, C., 2016. Development of a ship weather routing system. Ocean
Motion, SPEEDAM 2014, pp. 313–318. Eng. 123, 1–14.
Altosole, M., Piastra, F., Canepa, E., 2016. Performance analysis of a motor-sailing Zaccone, R., Figari, M., Altosole, M., Ottaviani, E., 2016. Fuel saving-oriented 3D dynamic
propulsion system for control design purposes. Ships Offshore Struct. 11 (7), programming for weather routing applications. In: Proceedings of the 3rd
688–699. International Conference on Maritime Technology and Engineering, MARTECH 2016,
Bellman, R., 1954. The Theory of Dynamic Programming, Technical Report, DTIC p. 183189.
Document. Zoppoli, R., 1972. Minimum-time routing as an n-stage decision process. J. Appl.
Bellman, R., 1958. On a routing problem. Q. Appl. Math. 87–90. Meteorol. 11 (3), 429–435.

224