Error Bounds for Discrete-Continuous Free Flight Trajectory Optimization

Abstract

Two-stage methods addressing continuous shortest path problems start local minimization from discrete shortest paths in a spatial graph. The convergence of such hybrid methods to global minimizers hinges on the discretization error induced by restricting the discrete global optimization to the graph, with corresponding implications on choosing an appropriate graph density. A prime example is flight planning, i.e., the computation of optimal routes in view of flight time and fuel consumption under given weather conditions. Highly efficient discrete shortest path algorithms exist and can be used directly for computing starting points for locally convergent optimal control methods. We derive a priori and localized error bounds for the flight time of discrete paths relative to the optimal continuous trajectory, in terms of the graph density and the given wind field. These bounds allow designing graphs with an optimal local connectivity structure. The properties of the bounds are illustrated on a set of benchmark problems. It turns out that localization improves the error bound by four orders of magnitude, but still leaves ample opportunities for tighter error bounds by a posteriori estimators.

1 Introduction

There are applications of numerical optimization that call for the computation of global instead of just local optima. One example is free flight planning, an instance of airborne navigation, where travel time is optimized subject to a given wind field (the travel time T between origin and destination is almost proportional to fuel consumption, CO$_2$ emission, and cost [22]). Going left or right around obstacles or adverse wind situations gives rise to locally optimal trajectories with considerably different costs [19], and airlines are naturally interested in the best of those.

Various approaches to global optimization have been proposed: stochastic ones like multistart or simulated annealing [6], biologically inspired metaheuristics like genetic algorithms or particle swarms [14], and rigorous ones based on objective bounds and branching [2, 9]. The former approaches converge to a global minimizer only almost surely at increasing computational costs, but provide no guarantees for finiteness. The latter ones usually require in-depth structural knowledge of the objective, or the use of interval arithmetics, and quickly suffer from the curse of dimensionality for practically relevant problems.

In this paper we consider a two-stage multistart approach along the following lines. It (i) defines a sufficiently large set of possible starting points, (ii) selects few promising candidates, and (iii) performs local optimization starting from those candidates. For this to be computationally feasible, the representation and selection of starting points needs to be highly efficient even for large and high-dimensional design spaces. This is, of course, problem-dependent. Some problems allow a discretization in terms of discrete network optimization problems such as minimum cost flow and, in particular, shortest path problems, which can be solved efficiently to global optimality in theory and practice [13]. If such discrete problems are close to their continuous counterpart, their solutions might provide promising starting points for local optimization to converge to a nearby global optimizer.

Obviously, flight planning and discrete shortest path search are related in this way and can hence serve as examples to substantiate the general idea. The starting points covering the design space of trajectories between origin and destination can be implicitly described as paths in a graph covering the spatial domain. In this discrete approximation of the problem, the selection of promising candidate points can be efficiently performed using Dijkstra’s algorithm or its $A^*$ variants. This leads to a hybrid discrete-continuous algorithm combining discrete global optimization methods with continuous local optimal control methods.

A first successful step in this direction has been taken with the development of the hybrid algorithm DisCOptER [7] that has been proposed by the authors of this paper for free flight planning. Even though the potential applications are manifold, no other similar method has been proposed so far. Note that the discrete stage alone is traditionally used as a standalone optimizer for practical flight planning on given airway networks [4, 5], but gets quickly inefficient when the airway networks need to be refined significantly to exploit the potential benefits of free flight [22]. Similarly, for robot path planning, rapidly exploring random graphs and trees (RRT) are used for sampling the trajectory design space at many discrete points [23].

Guaranteed convergence of a hybrid two-stage algorithm to a global minimizer hinges on the one hand on a sufficiently dense sampling of possible starting points in the design space, and on the other hand on the ability of the local optimizer to converge reliably to a nearby local optimum when started from one of these candidate points. The present paper investigates the first aspect, i.e., we derive bounds on the required resolution of the discretization. To this purpose, we introduce a continuous problem formulation that allows a direct comparison of continuous and discrete 2D flight paths (Sect. 2), and derive bounds for the flight duration deviation $T(\xi )-T(\xi _C)$ between different paths $\xi$ and $\xi _C$ in terms of spatial distance $\Vert \xi -\xi _C\Vert$, angular distance $\Vert (\xi -\xi _C)_\tau \Vert$, and bounds on the stationary wind and its derivatives. Based on the (h, l) graph density property from [7], we obtain corresponding flight duration bounds for discrete optimal trajectories (Sect. 3), which also yield a theoretically optimal ratio $h=\mathcal {O}(l^2)$ of vertex distance and characteristic edge length. We derive two types of error bounds: an a priori bound $T(\xi _G)\le T(\xi _C) + \kappa l^2$ depending only on problem quantities but not on a particular solution, and a local bound based on bounds for the wind in a neighborhood of the optimal trajectory. Taking more detailed information into account, the latter one improves on the former one by several orders of magnitude. The theoretical predictions are confirmed by numerical examples for a set of benchmark problems with varying wind complexity (Sect. 4), which reveal that there is still ample room for improvement by using a posteriori error estimators. The flight planning application leaves its imprint on the nature and derivation of these bounds, but the general idea should work for similar applications that have a discrete-continuous nature.

2 Shortest Flight Planning: Continuous & Discrete

For simplicity of presentation, we consider flight planning in the Euclidean plane. We aim at minimizing the travel time T between an origin $x_O$ and a destination $x_D$, with a fixed departure at $t=0$ and a constant airspeed $\overline{v}> 0$, thus neglecting start and landing phase. Moreover we assume a spatially heterogeneous, twice continuously differentiable wind field w to be given, with a bounded magnitude $\Vert w\Vert _{L^\infty (\mathbb {R}^2)} < \overline{v}$. Focusing on free flight areas, we also neglect any traffic flight restrictions.

2.1 Continuous: Optimal Control

In free flight areas, the flight trajectory is not restricted to a predefined set of airways. Instead, we consider any Lipschitz-continuous path $x:[0,T]\rightarrow \mathbb {R}^2$ in the Sobolev space $H^1([0,T])$, connecting origin $x_O$ and destination $x_D$, as a valid trajectory if it satisfies the following ODE almost everywhere,

$$\begin{aligned} x_t(t) = v(t) + w(x(t),t), \end{aligned}$$

(1)

where $x_t$ denotes the derivative of x with respect to t and is obtained by adding the vectors of airspeed and wind. The airspeed $v\in L^2([0,T]): [0,T]\mapsto \mathbb {R}^2$ lives in the Lebesgue space of square integrable functions. Among those trajectories, we need to find one with minimal flight duration T, since that is essentially proportional to fuel consumption [22]. This classic of optimal control is known as Zermelo’s navigation problem [24].

The optimal control problem of finding a trajectory (x, v) that minimizes travel time $T\in \mathbb {R}$ while obeying the dynamics described in (1) and travelling at constant airspeed $\overline{v}$ now reads

$$\begin{aligned} \min _{T,x,v} T \quad \text {s.t.}\quad c(T,x,v) = \begin{bmatrix} x(0) - x_O \\ x(T) - x_D \\ x_t(t) - (v(t) + w(x(t),t)) \\ v(t)^\textrm{T} v(t) - \overline{v}^2 \end{bmatrix} = 0 \end{aligned}$$

(2)

with $c:\mathbb {R}\times H^1([0,T])^2 \times L^2([0,T])^2 \rightarrow \mathbb {R}^2 \times \mathbb {R}^2 \times L^2([0,T])^2 \times L^2([0,T])$. Note that due to $v^\textrm{T}v = \overline{v}^2$, the airspeed v (and therefore also the ground speed $v+w$) is bounded almost everywhere, such that any feasible trajectory x is Lipschitz-continuous, i.e., $x \in C^{0,1}([0,T])$ [21, Thm. 1.36]. Moreover, it is immediately clear that there is an ellipse $\varOmega \subset \mathbb {R}^2$ with focal points $x_O$ and $x_D$, in which any trajectory with minimal flight duration is contained.

Problem (2) can be numerically solved efficiently with either direct methods using a discretization of the variables to formulate a finite-dimensional nonlinear programming problem [10], or with indirect methods relying on Pontryagin’s maximum principle, leading to a boundary value problem for ordinary differential equations [12, 16,17,18, 20, 24]. These approaches have also been considered explicitly for free flight planning [3, 11].

While the optimal control formulation (2) is convenient for numerically solving the optimization problem, we will consider a different formulation defining trajectories on the unit interval that is better suited for direct comparison with graph-based approaches here. Assume the flight trajectory $x:[0,T]\rightarrow \varOmega$ is given by a strictly monotonously increasing parametrization $t(\tau )$ on [0, 1] as $x(t(\tau )) = \xi (\tau )$, and $\xi \in H^1([0,1]): [0,1]\rightarrow \varOmega$ being a Lipschitz-continuous path with $\xi (0)=x_O$, $\xi (1)=x_D$. Due to Rademacher’s theorem, its derivative $\xi _\tau$ exists almost everywhere, and we assume it not to vanish. Then, $t(\tau )$ is defined by the state equation $x_t \overset{(1)}{=}\ v+w \ne 0$ and the airspeed constraint $\Vert v\Vert = \overline{v}$, since

$$\begin{aligned} \overline{v}= \Vert x_t - w \Vert \quad \text {and}\quad x_t t_\tau = \xi _\tau \ne 0 \end{aligned}$$

imply

$$\begin{aligned}&(t_\tau ^{-1}\xi _\tau -w)^\textrm{T}(t_\tau ^{-1}\xi _\tau -w) = \overline{v}^2 \\&\quad \Leftrightarrow t_\tau ^{-2}\xi _\tau ^\textrm{T} \xi _\tau -2 t_\tau ^{-1}\xi _\tau ^\textrm{T} w + w^\textrm{T}w - \overline{v}^2 = 0 \\&\quad \Leftrightarrow (\overline{v}^2 - w^\textrm{T}w) t_\tau ^2 + 2\xi _\tau ^\textrm{T} w t_\tau - \xi _\tau ^\textrm{T} \xi _\tau = 0 \end{aligned}$$

due to $t_\tau > 0$. Solving the quadratic equation yields

$$\begin{aligned} \quad&t_\tau = \frac{-\xi _\tau ^\textrm{T}w + \sqrt{(\xi _\tau ^\textrm{T}w)^2+(\overline{v}^2 - w^\textrm{T}w)(\xi _\tau ^\textrm{T} \xi _\tau )}}{\overline{v}^2 - w^\textrm{T}w} =: f(t,\xi ,\xi _\tau ). \end{aligned}$$

(3)

The flight duration T is then given by integrating the ODE (3) from 0 to 1 as $T=t(1)$. For the ease of presentation let us assume that the wind w is stationary, i.e., independent of t, and thus $f(t,\xi ,\xi _\tau ) = f(\xi ,\xi _\tau )$. Doing so, we avoid the more complicated work with an ODE. Instead, we obtain

$$\begin{aligned} T(\xi ) = \int _0^1 f\big (\xi (\tau ),\xi _\tau (\tau )\big )\, \hbox {d}\tau . \end{aligned}$$

(4)

We, however, strongly expect our results to directly carry over to the more complex case. Since the flight duration T as defined in (4) is based on a reparametrization $x(t) = \xi (\tau (t))$ of the path such that $\Vert x_t(t)-w(x(t))\Vert = \overline{v}$, the actual parametrization of $\xi$ is irrelevant for the value of T. Calling two paths $\xi ,\hat{\xi }$ equivalent if there exists a Lipschitz-continuous bijection $r:[0,1]\rightarrow [0,1]$ such that $\hat{\xi }(r(\tau )) = \xi (\tau )$, we can restrict the optimization to equivalence classes $[\xi ]$. Thus, the admissible set is

$$\begin{aligned} X = \{[\xi ] \mid \xi \in C^{0,1}([0,1],\varOmega ), \; \xi (0)=x_O, \; \xi (1)=x_D\}. \end{aligned}$$

(5)

Example 2.1

This is illustrated in Fig. 1. Consider the case that the wind gets stronger the farther the airplane proceeds to the right with $w(\xi (\tau ), \tau ) = \begin{bmatrix}{\bar{w}} &{}0\\ 0&{}0\end{bmatrix}\xi (\tau )$. Obviously, the optimal route is the straight line. This route can be represented differently, depending on the choice of the time parametrization $\tau (t)$. E.g., with $\tau (t) = t/T$ constant air speed is maintained. Consequently, the distance travelled in a certain time step increases over the course of the flight, due to the increasing wind speed (Fig. 1a). Alternatively, one may choose $\tau (t)$ such that the ground speed is constant (Fig. 1b).

Fig. 1

Illustration of the effect of the time parametrization $\tau (t)$. Blue dots represent position at equidistant time steps in the normalized time $\tau \in [0,1]$. Wind is blowing to the right with increasing speed

Since every equivalence class contains a representative with constant ground speed $\Vert \xi _\tau (\tau )\Vert$, we will subsequently often assume $\Vert \xi _\tau (\tau )\Vert =\textrm{const}$ without loss of generality, such that $\Vert \xi _\tau \Vert$ is just the length of the flight trajectory. For convenience, let us define the set of representatives with constant ground speed as

$$\begin{aligned} {\hat{X}} = \{ \xi \mid [\xi ]\in X, \Vert \xi _\tau \Vert =\textrm{const} \text { f.a.a. } \tau \in [0,1] \}. \end{aligned}$$

(6)

The reduced minimization problem, equivalent to (2), now reads

$$\begin{aligned} \min _{[\xi ] \in X} T(\xi ), \quad \text {or, equivalently,} \quad \min _{\xi \in {\hat{X}}} T(\xi ). \end{aligned}$$

(7)

Remark 2.1

Let us interpret this representation of flight duration. In the absence of wind, i.e., $\Vert w\Vert =0$, (3) yields $t_\tau = \xi _\tau / \overline{v}$. Integrating over [0, 1] yields just the total path length divided by the velocity (airspeed and ground speed coincide). For low wind, i.e., $\Vert w\Vert \ll \overline{v}$, we obtain $t_\tau \approx (\xi _\tau -\xi _\tau ^\textrm{T} w) / \overline{v}$ from (3), and hence a reduction of flight duration due to the tail wind component $\xi _\tau ^\textrm{T} w$ (or an increase in the case $\xi _\tau ^\textrm{T}w<0$ of head wind). For $\Vert w\Vert \rightarrow \overline{v}$, we obtain $t_\tau \rightarrow \xi _\tau / (2 \Vert \xi _\tau \Vert ^{-1} |\xi _\tau ^\textrm{T} w|)$ from (3) in case of a tailwind component $\xi _\tau ^\textrm{T} w>0$ and $t_\tau \rightarrow \infty$ otherwise. In any case, flight duration scales linearly with the length of the path.

In contrast to the optimal control formulation (2), the reduced formulation (7) allows a direct comparison of continuous and discrete flight trajectories, and is therefore the ideal tool for deriving error bounds in Sect. 3. We point out, however, that it is less suited for actually computing an optimal solution.

2.2 Discrete: Airway Networks

If flight trajectories are restricted to certain airways connecting predefined waypoints, flight planning is a special kind of shortest path problem on a graph. Let $V\subset \mathbb {R}^2$ be a finite set of waypoints including $x_O$ and $x_D$, and $E\subset V\times V$ a set of airways such that $G=(V,E)$ is a connected directed graph. A discrete flight path is a finite sequence $(x_i)_{0\le i \le n}$ of waypoints with $(x_{i-1},x_i)\in E$ for $i=1,\dots ,n$, connecting $x_0=x_O$ with $x_n=x_D$.

We define a mapping $\varXi :(x_i)_{0\le i\le n} \mapsto [\xi ]\in X$ of discrete flight paths to continuous paths by piecewise linear interpolation

$$\begin{aligned} \xi (\tau ) = x_{\lfloor n\tau \rfloor } + (n\tau -\lfloor n\tau \rfloor )(x_{\lceil n\tau \rceil }-x_{\lfloor n\tau \rfloor }) \end{aligned}$$

(8)

resulting in polygonal chains, which are Lipschitz-continuous with piecewise constant derivative. We denote its image $\mathop \textrm{im} \varXi \subset X$, i.e., the set of flight trajectories in the Euclidean plane that can be realized by adhering to the airway network, by $X_G$. The discrete flight planning problem then reads

$$\begin{aligned} \min _{[\xi ]\in X_G} T(\xi ), \end{aligned}$$

(9)

and differs from its continuous counterpart (7) only by the admissible set, effectively acting as a particular discretization.

Shortest path problems on static graphs with non-negative weights are usually solved with the $A^*$ variant of Dijkstra’s algorithm [15].

3 Approximation Error Bounds

Having established a setting in which discrete and continuous flight trajectories can be directly compared, we are interested in bounding the suboptimality, i.e., the increase of flight duration T relative to the continuous optimum, due to restricting the flight path to predefined airways. In particular, we aim at relating this approximation error to the airway network density.

3.1 A Posteriori Error

For estimating the flight time deviation, we start with a Taylor-based bound in terms of the actual path deviation $\delta \xi$ from a minimizer. This bound will serve as the basis for computable bounds in Sect. 3.3 and, in addition, provide a quantitative idea of the efficiency of a posteriori error estimators using computable estimates of $\Vert \delta \xi \Vert$.

At this point we want to point out that $\xi$, $\xi _\tau$, $\delta \xi$, and $\delta \xi _\tau$ are in general functions of $\tau$. In favor of a more compact notation we will usually omit the argument $\tau$ in the remainder of the paper.

Lemma 3.1

For any $p\in \varOmega$ let $c_0(p)=\Vert w(p)\Vert$, $c_1(p)=\Vert w_x(p)\Vert$, and $c_2(p)=\Vert w_{xx}(p)\Vert$, and assume $c_0 \le \overline{v}/ \sqrt{5}$. Moreover, let $\xi \in {\hat{X}}$, $L:=\Vert \xi _\tau \Vert >0$ and ${v}^2(p):= \overline{v}^2 - c_0^2(p)$. Then the second total directional derivative of f as defined in (3) is bounded by

$$\begin{aligned} |f''(\xi ,\xi _\tau ) [\delta \xi ,\delta \xi _\tau ]^2| \le \;&\alpha _0(\xi ) \Vert \delta \xi \Vert ^2 + \alpha _1(\xi ) \Vert \delta \xi \Vert \, \Vert \delta \xi _\tau \Vert + \alpha _2(\xi ) \Vert \delta \xi _\tau \Vert ^2 \end{aligned}$$

(10)

for almost all $\tau \in [0,1]$, with $\alpha _i: \varOmega \rightarrow \mathbb {R}^+$, $i=0,\dots ,2$, given as

$$\begin{aligned} \alpha _0(p)&= \frac{L}{\underline{v}^{3}(p)} \left( 12 c_1^2(p) + 4 \underline{v}(p) c_2(p)\right) \,, \\ \alpha _1(p)&= \frac{8 c_1(p)}{\underline{v}^2(p)}\,,\\ \alpha _2(p)&= \frac{2}{ L \underline{v}(p)}. \end{aligned}$$

The proof of this lemma is, though not difficult, rather technical and lengthy calculus and is provided in the “Appendix”.

Remark 3.1

The assumption of $c_0\le \overline{v}/\sqrt{5}$ covers the usually experienced wind velocities, but not the possible extremes.

Theorem 3.1

Let $\xi _C\in {\hat{X}}$ be a minimizer of (7) and $\delta \xi := \xi -\xi _C$. Then there is a constant $r>0$ depending on $\xi _C$ and w, such that the a posteriori bound

$$\begin{aligned} T(\xi ) \le T(\xi _C) + \int _0^1 \big (&\alpha _0(\xi _C) \Vert \delta \xi \Vert ^2+ \alpha _1(\xi _C) \Vert \delta \xi \Vert \Vert \delta \xi _\tau \Vert + \alpha _2(\xi _C) \Vert \delta \xi _\tau \Vert ^2 \big ) \, \hbox {d}\tau \end{aligned}$$

(11)

holds for all paths $\xi \in {\hat{X}}$ with $\Vert \xi -\xi _C\Vert _{C^{0,1}([0,1])} \le r$ and $\alpha _i$ as defined in Lemma 3.1.

Proof

We note that $T:C^{0,1}([0,1])^2 \rightarrow \mathbb {R}$ as defined in (4) is twice continuously Fréchet-differentiable at $\xi _C\in {\hat{X}}$ due to $\Vert (\xi _C)_\tau \Vert \underset{(6)}{=} L>0$ for almost all $\tau$. By Lemma 3.1, there are functions $\alpha _0,\alpha _1,\alpha _2$ depending on the local wind w and its derivatives as well as the overall trajectory length L, such that

$$\begin{aligned} f''(\xi _C,(\xi _C)_\tau )[(\delta \xi ,\delta \xi _\tau ),{} & {} (\delta \xi ,\delta \xi _\tau )] \\{} & {} \le \alpha _0(\xi _C) \Vert \delta \xi \Vert ^2 + \alpha _1(\xi _C) \Vert \delta \xi \Vert \Vert \delta \xi _\tau \Vert + \alpha _2(\xi _C) \Vert \delta \xi _\tau \Vert ^2 \end{aligned}$$

holds for almost all $\delta \xi ,\delta \xi _\tau \in \mathbb {R}^2$ and $\tau \in [0,1]$. Integrating over $\tau$ yields the bound

$$\begin{aligned} T''(\xi _C)[\delta \xi ,\delta \xi ] \le \int _0^1 \big (&\alpha _0(\xi _C) \Vert \delta \xi \Vert ^2 + \alpha _1 (\xi _C)\Vert \delta \xi \Vert \Vert \delta \xi _\tau \Vert + \alpha _2(\xi _C) \Vert \delta \xi _\tau \Vert ^2\big ) \, \hbox {d}\tau \end{aligned}$$

for second directional derivatives of the flight duration T in direction $\delta \xi \in C^{0,1}([0,1])^2$ with $\delta \xi (0)=\delta \xi (1)=0$. Due to continuity of $T''$, there exists a neighborhood $B_r(\xi _C)$ of radius $r>0$, such that $T''(\tilde{\xi })[\delta \xi ,\delta \xi ] \le 2\int _0^1 \alpha _0 \Vert \delta \xi \Vert ^2 + \alpha _1 \Vert \delta \xi \Vert \Vert \delta \xi _\tau \Vert + \alpha _2 \Vert \delta \xi _\tau \Vert ^2 \hbox {d}\tau$ for all $\tilde{\xi }\in B_r(\xi _C)$. Consequently, by Taylor’s theorem and using $\delta \xi =\xi -\xi _C$, we can bound

$$\begin{aligned} T(\xi )&= T(\xi _C) + \underbrace{T'(\xi _C)[\delta \xi ]}_{=0} + \int _0^1 (1-\nu ) T''(\xi _C + \nu \delta \xi )[\delta \xi ,\delta \xi ] \, \hbox {d}\nu \\&\le T(\xi _C) + \int _0^1 \big ( \alpha _0(\xi _C) \Vert \delta \xi \Vert ^2 + \alpha _1(\xi _C) \Vert \delta \xi \Vert \Vert \delta \xi _\tau \Vert + \alpha _2(\xi _C) \Vert \delta \xi _\tau \Vert ^2\big ) \,\hbox {d}\tau \end{aligned}$$

due to $\xi _C$ being a minimizer. $\square$

3.2 Trajectory Approximation in Locally Dense Graphs

The approximation error of the optimal discrete flight path $\xi _G$ according to (9) relative to the continuous optimum $\xi _C$ of (7) due to the smaller admissible set $X_G \subset X$ depends on the density of the airway network. The discussion will be limited to a certain class of locally dense digraphs as defined in [7].

Definition 3.1

A digraph $G=(V,E)$ is said to be (h, l)-dense in a convex set $\varOmega \subset \mathbb {R}^2$ for $h,l\ge 0$, if it satisfies the following conditions:

• 1. containment: $V\subset \varOmega$

• 2. vertex density: $\forall p\in \varOmega : \exists v\in V: \Vert p-v\Vert \le h$

• 3. local connectivity: $\forall v,w\in V, \Vert v-w\Vert \le l+2h: (v,w)\in E$

An example for such an airway digraph is shown in Fig. 2. Note that, even for $l\rightarrow 0$, the minimum local connectivity length of $2\,h$ guarantees that a vertex is connected to its neighbors. It is easy to show that any (h, l)-dense digraph is connected, such that a path from origin to destination exists.

Fig. 2

A locally densely connected digraph with Cartesian structure. The center node (dark blue) is connected to all nodes in a circular neighborhood of radius $2h+l$ (light blue) with edges in both directions

Let $\xi _C\in X$ be a global minimizer of the continuous problem formulation (7), and $\xi _G\in X_G$ be a shortest discrete path in the (h, l)-dense airway digraph G satisfying (9). For establishing a bound for the excess flight duration in terms of the airway density, we first construct a particular discrete path $\xi _R(\xi _C)\in X_G$ using a rounding procedure, and derive a bound for $T(\xi _R)-T(\xi _C)\le e(h,l)$, from which the actual error bound $T(\xi _G)-T(\xi _C)\le e(h,l)$ immediately follows from optimality of $\xi _G$.

For defining $\xi _R(\xi )\in X_G$, with an (h, l)-dense digraph G with $l>0$, for a given continuous path $\xi \in {\hat{X}}$ with $\xi _\tau =\textrm{const}$, we first choose an equidistant grid $\tau _i = i/n$, $n=\lceil \xi _\tau / l \rceil$, for $i=0,\dots ,n$. By construction, the distance of the corresponding trajectory points is bounded by $\Vert \xi (\tau _i)-\xi (\tau _{i+1})\Vert \le l$. For each i, there is some $v_i \in V$ with $\Vert v_i - \xi (\tau _i)\Vert \le h$, such that $\Vert v_i - v_{i+1}\Vert \le l+2\,h$. Consequently, $(v_i,v_{i+1})\in E$, and $(v_i)_{0\le i \le n}$ is a valid discrete path, for which we define $[\xi _R] = \varXi (v_i)_i$.

It is intuitively clear—and rigorously confirmed below—that the excess flight duration $T(\xi _R)-T(\xi _C)$ is affected by both, the spatial distance between $\xi _R$ and $\xi _C$, e.g., taking a longer detour or flying through an area with adverse wind conditions, and the angular deviation, e.g., a zigzag path tends to take longer than a straight trajectory. In order to capture these effects, we will first bound the spatial distance $\Vert \xi _R-\xi _C\Vert _{L^\infty [0,1]}$ and the angular deviation $\Vert (\xi _R-\xi _C)_\tau \Vert _{L^\infty [0,1]}$, and equip the space of Lipschitz-continuous functions with the norm $\Vert f\Vert _{C^{0,1}([0,1])}:= \Vert f\Vert _{L^\infty ([0,1])} + \Vert f_\tau \Vert _{L^\infty ([0,1])}$.

Theorem 3.2

Assume $\xi \in {\hat{X}} \cap C^{1,1}([0,1],\varOmega )$ has bounded curvature, i.e., there is some $\bar{\sigma }$ with $\Vert \xi _\tau (a)-\xi _\tau (b)\Vert \le \overline{\sigma }|a-b|$ for $a,b\in [0,1]$, and denote the length of the trajectory by $L=\Vert \xi _\tau \Vert$. Then, the following bounds hold for the discrete approximation $\xi _R(\xi )$ in an (h, l)-dense digraph:

$$\begin{aligned}&\mathrm{(distance~error)}&\Vert \xi _R(\xi )-\xi \Vert _{L^\infty ([0,1])}&\le \frac{\overline{\sigma }l^2}{8L^2} + h\,, \end{aligned}$$

(12)

$$\begin{aligned}&\quad \mathrm{(angular~error)}&\Vert (\xi _R(\xi )-\xi )_\tau \Vert _{L^\infty ([0,1])}&\le \frac{\sqrt{2}\overline{\sigma }l}{L} + 2h\left( \frac{L}{l} +1\right) . \end{aligned}$$

(13)

If $l\le L$, we obtain the total error bound

$$\begin{aligned} \Vert \xi _R(\xi )-\xi \Vert _{C^{0,1}([0,1])}&\le \left( \frac{1}{8}+\sqrt{2}\right) \bar{\sigma }\frac{l}{L} + 2h\frac{L}{l} + 3h \nonumber \\&\le 2 \bar{\sigma }\frac{l}{L} + 2h\frac{L}{l} + 3h. \end{aligned}$$

(14)

Proof

Let $\hat{\xi }(\tau ) = \xi (\tau _{\lfloor n\tau \rfloor }) + (n\tau -\lfloor n\tau \rfloor )(\xi (\tau _{\lceil n\tau \rceil })-\xi (\tau _{\lfloor n\tau \rfloor }))$ be the linear interpolant of the continuous trajectory $\xi$ on $n=\lceil L / l\rceil$ equisized intervals. Standard interpolation error estimates yield

$$\begin{aligned} \Vert (\hat{\xi }- \xi )(\tau )\Vert \le \overline{\sigma }/ (8n^2) \le \frac{\overline{\sigma }l^2}{8L^2} \end{aligned}$$

for all $\tau$ [1, Ch. 3.1, p. 93 ff.]. Moreover, with $\alpha = (n\tau -\lfloor n\tau \rfloor ) \in [0,1]$,

$$\begin{aligned} \hat{\xi }(\tau )-\xi _R(\tau )&= \xi (\tau _{\lfloor n\tau \rfloor }) - x_{\lfloor n\tau \rfloor } + \alpha \Big ( \xi (\tau _{\lceil n\tau \rceil })-x_{\lceil n\tau \rceil }-\xi (\tau _{\lfloor n\tau \rfloor })+x_{\lfloor n\tau \rfloor } \Big ) \nonumber \\&= (1-\alpha )\big ( \xi (\tau _{\lfloor n\tau \rfloor })- x_{\lfloor n\tau \rfloor } \big ) + \alpha \big ( \xi (\tau _{\lceil n\tau \rceil })-x_{\lceil n\tau \rceil } \big ) \end{aligned}$$

(15)

implies $\Vert (\hat{\xi }- \xi _R)(\tau )\Vert \le h$, which yields the distance error bound (12) by triangle inequality.

Let $\phi = (\hat{\xi }- \xi )_k$, $k\in \{1,2\}$, be one of the two components of the difference between continuous trajectory and linear interpolant. By the mean value theorem, there is a point $\hat{\tau }\in \mathopen ]\tau _i,\tau _{i+1}\mathclose [$ with $\phi _\tau (\hat{\tau }) = 0$. Thus,

$$\begin{aligned} |\phi _\tau (\tau )| = |\phi _\tau (\tau ) -\phi _\tau (\hat{\tau })| \le \frac{\overline{\sigma }}{n} \quad \forall \tau \in [\tau _i,\tau _{i+1}] \end{aligned}$$

holds for all $i=0,\dots ,n-1$ and implies $\Vert (\hat{\xi }-\xi )_\tau (\tau )\Vert \le \sqrt{2}\overline{\sigma }/ n \le \sqrt{2}\overline{\sigma }l / L$ for all $\tau$. Moreover, (15) implies

$$\begin{aligned} (\hat{\xi }-\xi _R)_\tau (\tau ) = -n \big ( \xi (\tau _{\lfloor n\tau \rfloor })- x_{\lfloor n\tau \rfloor } \big ) + n \big ( \xi (\tau _{\lceil n\tau \rceil })-x_{\lceil n\tau \rceil } \big ) \end{aligned}$$

and therefore $\Vert (\hat{\xi }-\xi _R)_\tau \Vert \le 2nh \le 2h(L/l +1)$ and yields the angular error bound (13) by triangle inequality. $\square$

Of course, if $l < {\hat{l}}$, then the (h, l)-dense digraph G is a subgraph of the $(h,{\hat{l}})$-dense digraph ${\hat{G}}$, provided their vertex sets coincide. Thus, the discretization error of a shortest path in ${\hat{G}}$ is less or equal to one in G—a fact that is not reflected by Theorem 3.2. The reason is the explicit rounding procedure, which tends to select arcs of length ${\hat{l}}$ in ${\hat{G}}$ even if shorter arcs of length l would be better. This effect can be essentially avoided if the connectivity length l is chosen sufficiently small compared to the path length. It should not be chosen too small compared to h, however, because then the angular error can dominate, as the following pathological example shows.

Fig. 3

Illustration of Example 3.1. Green: continuous trajectory $\xi$, gray: rounded path $\xi _R$

Example 3.1

Consider $\xi (\tau ) = [\tau ,0]^\textrm{T}$ and

$$\begin{aligned} V = \{[l(2i+j),2j-1]^\textrm{T} \mid i,j\in \mathbb {Z}\} \cup \{[0,0]^\textrm{T}, [1,0]^\textrm{T}\} \end{aligned}$$

with $l\ll 1$ and $h=\sqrt{1+l^2/4}\approx 1$. Rounding to the nearest vertex yields a discrete zigzag path with length at least 2h/l, as illustrated in Fig. 3. Thus, the bounds (13) and (14) are asymptotically sharp for $l\rightarrow 0$.

Hence, we select a theoretically optimal l by minimizing the error bound (14).

Theorem 3.3

Under the assumptions of Theorem 3.2, including $l\le L$, the choice

$$\begin{aligned} l = L \sqrt{\frac{h}{\bar{\sigma }}} \quad \Leftrightarrow \quad h = \bar{\sigma } \frac{l^2}{L^2} \end{aligned}$$

is optimal with respect to the error bound (14) and yields the bounds

$$\begin{aligned}&\Vert \xi _R(\xi )-\xi \Vert _{L^\infty ([0,1])} \le \frac{9\overline{\sigma }l^2}{8L^2}, \end{aligned}$$

(16)

$$\begin{aligned}&\Vert (\xi _R(\xi )-\xi )_\tau \Vert _{L^\infty ([0,1])} \le \frac{11\overline{\sigma }l}{2L}, \end{aligned}$$

(17)

$$\begin{aligned}&\Vert \xi _R(\xi )-\xi \Vert _{C^{0,1}([0,1])} \le 7\bar{\sigma } \frac{l}{L}. \end{aligned}$$

(18)

Proof

Straightforward minimization of (14) yields the given optimal choice of l. Inserting this into (12), (13), and the bound (14) and using $l\le L$ yields the claims. $\square$

The pathological Example 3.1 reveals a further limitation of the derivation of bounds by employing an explicit rounding procedure: the length of the rounded path $\xi _R$ can be much larger than the length of the discretely optimal path $\xi _G$. In the example this is $\mathcal {O}(2/l)\rightarrow \infty$ compared to $\mathcal {O}(1)$, with $\xi _G$ connecting the vertices along the horizontal line $[0,1]\times \{1\}$. We point out that this susceptibility of the bound to pathological worst cases is structurally similar to common a priori error estimates for finite element methods [8]. Nevertheless, even if the angular error responsible for the pathological behavior is ignored, the same optimal order of $h=\mathcal {O}(l^2)$ is obtained.

3.3 Computable Error Bounds

Theorem 3.4

Assume that $\xi _C\in {\hat{X}} \cap C^{1,1}([0,1])^2$ is a minimizer of (7) with bounded curvature, i.e., there is $\overline{\sigma }<\infty$ such that $\Vert \xi _\tau (a)-\xi _\tau (b)\Vert \le \overline{\sigma }|a-b|$ for all $a,b\in [0,1]$. Let $L=\Vert (\xi _C)_\tau \Vert$ denote the length of the optimal flight trajectory and $\overline{\alpha }_i:= \max _{\tau \in [0,1]} \alpha _i(\xi _C(\tau ))$ with $\alpha _i$ defined in Lemma 3.1. Then, there is a constant $r>0$, such that the $local$ bound

$$\begin{aligned} T(\xi _G)-T(\xi _C) \le \frac{4\overline{\sigma }^2 l^2}{3L^2} \left( \frac{l^2}{L^2} \overline{\alpha }_0 + 5\frac{l}{L}\overline{\alpha }_1 + 23 \overline{\alpha }_2 \right) \le \frac{92\overline{\sigma }^2 \overline{\alpha }_2}{3L^2} l^2 +\mathcal {O}(l^3) \end{aligned}$$

(19)

holds for all (h, l)-dense digraphs with $l\le \min \big \{\frac{r}{7\bar{\sigma }},1\big \}L$ and $h\le \frac{\bar{\sigma }l^2}{L^2}$.

Proof

Inserting the bounds (16) and (17) from Theorem 3.3 into the claim (11), we obtain

$$\begin{aligned} T(\xi _R) - T(\xi _C)&\le \int _0^1 \bigg ( \alpha _0 \frac{81\overline{\sigma }^2l^4}{64L^4} +\, \alpha _1 \frac{9\overline{\sigma }l^2}{8L^2} \frac{11\overline{\sigma }l}{2L} +\, \alpha _2 \frac{121\overline{\sigma }^2 l^2}{4L^2} \bigg ) \, \hbox {d}\tau \\&< \frac{4\bar{\sigma }^2 l^2}{3L^2} \int _0^1 \bigg ( \alpha _0 \frac{l^2}{L^2} +\, 5\alpha _1 \frac{l}{L} +\, 23 \alpha _2 \bigg ) \, \hbox {d}\tau , \end{aligned}$$

since $l \le \min \big \{\frac{r}{7\bar{\sigma }},1\big \}L$, where r is the neighborhood radius from Theorem 3.1 and $\alpha _i$ provided by Lemma 3.1. Inserting the upper bounds $\overline{\alpha }_i$ for $\alpha _i$ yields the claim. $\square$

Note that the bound holds in a certain neighborhood of a continuous minimizer $\xi _C$ and therefore bounds the asymptotic error behavior for $h,l\rightarrow 0$, rather than providing a globally reliable error bound.

We can go one step further and eliminate the dependence on the actual optimal path $\xi _C$ by choosing appropriate global bounds on the constants and route properties. For that, we define the global bounds

$$\begin{aligned} \overline{c}_0:= \Vert w\Vert _{L^\infty (\varOmega )}, \quad \overline{c}_1:= \Vert w_x\Vert _{L^\infty (\varOmega )}, \quad \text {and}\quad \overline{c}_2:= \Vert w_{xx}\Vert _{L^\infty (\varOmega )} \end{aligned}$$

for the wind and its derivatives.

Lemma 3.2

Let $\xi _C \in {\hat{X}}$ be a minimizer of (7). Further, let $\Vert w(p)\Vert \le \overline{c}_0$ and $\Vert w_x(p)\Vert \le \overline{c}_1 ~\forall p\in \varOmega$. Then, it is twice continuously differentiable and its second derivative is bounded by

$$\begin{aligned} \Vert (\xi _C)_{\tau \tau }\Vert \le \overline{\sigma }:= \frac{2 \overline{v}\overline{c}_1 L^2}{(\overline{v}-\overline{c}_0)^3} \left( (1+\sqrt{2})\overline{v}+ \overline{c}_0\right) . \end{aligned}$$

(20)

For $\overline{c}_0 \le \overline{v}/ \sqrt{5}$ this simplifies to $\overline{\sigma }\le 17\frac{\overline{c}_1L^2}{\overline{v}}$.

Again, the proof of this Lemma is rather long and can be found in the “Appendix”.

Lemma 3.3

Assume that $\xi _C$ is a global minimizer of (7) with path length L and that $\overline{c}_0 \le \overline{v}/ \sqrt{5}$. Then

$$\begin{aligned} \Vert x_D-x_O\Vert \le L \le \frac{\overline{v}+ \overline{c}_0}{\overline{v}- \overline{c}_0}\Vert x_D-x_O\Vert < \frac{8}{3}\Vert x_D-x_O\Vert . \end{aligned}$$

(21)

Proof

The lower bound is clear, since the trajectory cannot be shorter than the straight connection. The flight time $T_s$ on the straight line is at most $\frac{\Vert x_D-x_O\Vert }{\overline{v}-\overline{c}_0}$. Since $\xi _G$ is optimal, we obtain

$$\begin{aligned} T_s \ge T(\xi _C) \ge \frac{L}{\overline{v}+\overline{c}_0}, \end{aligned}$$

which yields the upper bound for L. $\square$

We can now completely eliminate the need for a posteriori information about $\xi _C$ and derive an a priori error bound.

Theorem 3.5

Let $\xi _C\in {\hat{X}} \cap C^{1,1}([0,1])^2$ be a global continuous minimizer and $\xi _G$ be a shortest path in the (h, l)-dense graph G. Moreover, let $\tilde{L}:= \Vert x_D - x_O\Vert$ and assume $\overline{c}_0 \le \overline{v}/ \sqrt{5}$. Then, with $\overline{\sigma }$ from Lemma 3.2, the following a priori error bound holds for sufficiently dense graphs:

$$\begin{aligned} T(\xi _G)-T(\xi _C)&\le \frac{4\overline{\sigma }^2}{3\tilde{L}^3 \overline{v}} \left( 14 l^2 \left( \frac{7}{2}\overline{c}_1^2 + \overline{v}\overline{c}_2\right) + \frac{51\overline{c}_1l}{\overline{v}} + 52 \right) \; l^2 \end{aligned}$$

(22)

$$\begin{aligned}&\le 1.5\cdot 10^5 \frac{\overline{c}_1^2 \tilde{L}}{\overline{v}^3} l^2 + \mathcal {O}(l^3). \end{aligned}$$

(23)

Proof

For $\underline{v}(p) = \sqrt{\overline{v}^2 - c_0^2(p)}$ we obtain $8\overline{v}/9 < \underline{v} \le \overline{v}$. Lemma 3.3 together with $\alpha _i$ from Lemma 3.1 now yields the global bounds

$$\begin{aligned} \alpha _0(p)&\le \frac{8\tilde{L}}{3\underline{v}(p)^3} (12\overline{c}_1^2 + 4\underline{v}(p)\overline{c}_2) \le \frac{14\tilde{L}}{\overline{v}} \left( \frac{7}{2}\overline{c}_1^2 + \overline{c}_2\overline{v}\right) =: \tilde{\alpha }_0, \\ \alpha _1(p)&\le \frac{8\overline{c}_1}{\underline{v}(p)^2} \le \frac{81\overline{c}_1}{8\overline{v}^2} =: \tilde{\alpha }_1, \quad \text {and} \\ \alpha _2(p)&\le \frac{2}{\tilde{L} \underline{v}(p)} \le \frac{9}{4\tilde{L} \overline{v}} =: \tilde{\alpha }_2. \end{aligned}$$

Inserting them into (19) provides the bound

$$\begin{aligned} T(\xi _G) - T(\xi _C)&\le T(\xi _R) - T(\xi _C) \\&\le \frac{4\overline{\sigma }^2 l^2}{3\tilde{L}^2} \left( \frac{l^2}{\tilde{L}^2}\tilde{\alpha }_0 + 5\frac{l}{\tilde{L}} \tilde{\alpha }_1 + 23\tilde{\alpha }_2\right) \\&\le \frac{4\overline{\sigma }^2 l^2}{3\tilde{L}^2} \left( \frac{14l^2}{\tilde{L} \overline{v}} \left( \frac{7}{2}\overline{c}_1^2 + \overline{c}_2\overline{v}\right) + \frac{405\overline{c}_1l}{8\overline{v}^2\tilde{L}} + \frac{207}{4\tilde{L} \overline{v}} \right) \\&\le \frac{4\overline{\sigma }^2 l^2}{3\tilde{L}^3 \overline{v}} \left( 14l^2 \left( \frac{7}{2}\overline{c}_1^2 + \overline{c}_2\overline{v}\right) + \frac{51\overline{c}_1l}{\overline{v}} + 52 \right) , \end{aligned}$$

which completes the proof. $\square$

4 Numerical Examples

In Sect. 3 we derived three error bounds: (i) the a priori bound (22), (ii) the local bound (19), and (iii) the computationally in general unavailable a posteriori bound (11). Now we validate these bounds with the four test problems from [7]. For this comparison and for evaluating the a posteriori bound we compute the optimal continuous trajectory $\xi _C$ numerically using a direct collocation approach to high accuracy.

4.1 Test Instances

The goal in all four test instances is to find a time-optimal trajectory from $x_O=[0,0]^\textrm{T}$ to $x_D=[1,0]^\textrm{T}$ through wind fields of varying spatial frequency, see Fig. 4. The wind speed is always bounded by $\Vert w\Vert \le {\bar{w}} = 0.5 \overline{v}$. All values are chosen dimensionless, i.e., $\overline{v}=1$. For the first test problem (a) we define the laminar shear flow

$$\begin{aligned} w(p) = \begin{bmatrix} {\bar{w}}\min (\max ( 2\frac{p_2}{H}{-}1,{-}1),1) \\ 0 \end{bmatrix}, \end{aligned}$$

with $H=0.5$, and $p_2$ denoting the y-component of p, see Fig. 4a. In problems (b)–(d), the wind w is the sum of an increasing number of non-overlapping vortices $w_i$, each of which is described by

$$\begin{aligned} w_i(p) = s_i \tilde{w}_i(r_i) \begin{bmatrix} - \sin (\alpha _i) \\ \cos (\alpha _i) \end{bmatrix}, \end{aligned}$$

where $s_i$ is the spin of the vortex ($s_i{=}{+}1$: counter-clockwise, $s_i{=}{-}1$: clockwise), $r_i = \Vert p-z_i\Vert _2$ is the distance from the vortex center $z_i$, $\alpha _i$ is the angle between p, $z_i$, and the positive x-axis with $\tan (\alpha _i) = \frac{(p-z_i)_2}{(p-z_i)_1}$ and the absolute vortex wind speed $\tilde{w}_i$ is a function of r and the vortex radius $R_i$:

$$\begin{aligned} \tilde{w}_i(r) = \left[ \begin{array}{cl} {\bar{w}} \exp \left( \frac{(r/R_i)^2}{(r/R_i)^2 - 1}\right) &{} \text {if}~r<R_i\\ 0 &{} \text {otherwise} \end{array}\right] \,. \end{aligned}$$

More precisely, problems (b)–(d) involve 1, 15, and 50 regularly aligned vortices with $R{=}1/2$, 1/8, and 1/16, respectively, see Fig. 4b–d. Vortices with positive spin (counter-clockwise) are colored green, vortices with negative spin (clockwise) are colored red.

Note that at least problem (d) is clearly an exaggeration, as no commercial plane would ever try to traverse a wind field like this. We use this instance to provide evidence that our claims hold even under the most adverse conditions.

Fig. 4

Test problems (a–d), $x_O=[0,0]^\textrm{T}$, $x_D=[1,0]^\textrm{T}$, ${\bar{w}}=\frac{1}{2} \overline{v}$. Blue dots: (h, l)-dense rectangular network-graph of some exemplary density and connectivity, satisfying $h=\overline{\sigma }l^2/L^2$, red: shortest path on the graph, green: continuous optimal trajectory. Note that in every case the straight trajectory is particularly unfavorable

4.2 Results

The three error bounds (i) a priori, (ii) local, and (iii) a posteriori involve an increasing amount of information about the optimal trajectory $\xi _C$. Hence it is not surprising that the first is far from being sharp and overestimates the actual error by several orders of magnitude. This is mainly due to the worst case estimates that must be considered in several steps finding globally valid constants. Most importantly, the wind speed and its derivatives need to be bounded globally, even though the conditions experienced along the actually flown path are usually much easier, as can be seen in Fig. 4. On the other hand, the a posteriori requires information about the optimal trajectory, that is not available at the time of calculation. However it helps us evaluate the more interesting question: how much can be gained by incorporating a posteriori information?

In order to answer this, we evaluate the error between a graph-based shortest path and the continuous optimum for various graph densities, respecting the optimal combination of node density and connectivity $h=\overline{\sigma }l^2/L^2$ as stated in Theorem 3.3. The results are shown in the double-logarithmic plots of Fig. 5 as blue dots together with a linear regression.

The three bounds are illustrated in (i) red, (ii) purple and (iii) gray. In the last case, the a posteriori bound (11), we depict the three individual parts of the integral by colored areas, from bottom to top: $\int \alpha _0 \Vert \delta \xi \Vert ^2 \hbox {d}\tau$, $\int \alpha _1 \Vert \delta \xi \Vert \Vert \delta \xi _\tau \Vert \hbox {d}\tau$, $\int \alpha _2 \Vert \delta \xi _\tau \Vert ^2 \hbox {d}\tau$. The linearly scaled depiction in Fig. 6 makes it easier to see that the relative distribution of these parts is more or less stable over a wide range of graph densities. Each part is relevant and none is vanishing even for dense graphs. This suggests that the theoretically optimal choice of (h, l) balances the error terms against each other evenly.

It needs to be mentioned, that, because we have only $\Vert w\Vert /\overline{v}\le 0.5$ here (not $\le 1/\sqrt{5}$), we cannot use $\alpha _i$ as stated in Lemma 3.1, but must revert to the results from Theorem A.1 in the “Appendix”. Since the purpose of that Lemma is solely to provide a more compact notation, however, this should not be of any concern. For the same reason, the coefficients in (22) also need to be adjusted accordingly.

Fig. 5

Evaluation of the derived error bounds (i) a priori (top, red), (ii) local (middle, purple), and (iii) a posteriori (bottom, gray), comprising three terms visualized by the colored areas. From top to bottom: $\int \alpha _0 \Vert \delta \xi \Vert ^2 \hbox {d}\tau$, $\int \alpha _1 \Vert \delta \xi \Vert \Vert \delta \xi _\tau \Vert \hbox {d}\tau$, $\int \alpha _2 \Vert \delta \xi _\tau \Vert ^2 \hbox {d}\tau$. The blue dots represent results from numerical experiments together with a linear regression line. The subfigures a–d refer to the corresponding test instances

We show linear trend lines for the a posteriori bound (iii) and the experimental results, excluding the data of the 10% sparsest graphs (rightmost data points). Results in that region are dominated by effects of local minima (e.g., the continuous optimum goes left, while the discrete path goes right, which results in a large distance error) which vanish quickly and do not contribute to the asymptotical error behavior. As our error bounds were developed to hold in a certain neighborhood of $\xi _C$, we ignore these effects. It is, however, interesting to see that the bounds hold anyway. This, together with the accumulation of data points on the left side, leads to a noticeable visual bias in the trend lines.

Remember that the (i) a priori (22) and (ii) the local error bound (19) both have the form

$$\begin{aligned} T(\xi _G) - T(\xi _C) \le c l^2 + \mathcal {O}(l^3), \end{aligned}$$

differing only in the constant c. As a first important result we point out that the quadratic order of these bounds matches the numerical results satisfyingly well. The exponents obtained from fitting a regression line to the numerical data and the evaluated a posteriori bound (iii) according to (11) are listed in Table 1.

Further, starting from the a priori bound (i), we note that the bound can be tightened significantly by incorporating a posteriori knowledge. With the local approach, the bound can already be improved by roughly four orders of magnitude, but taking all the a posteriori information into account clearly makes the biggest difference. In doing so, the bound comes close to the numerical data up to a factor of 6–11 (see Table 2) and can even resolve the aliasing artifacts.

Let us briefly discuss the visible oscillations in the actual errors. We consider the case (d), as the effect is most prominent here. The optimal solution is to quickly switch to a mostly horizontal trajectory in the middle between the first and second row of vortices and to switch back very late, using the spin of both the very first and the very last vortex. Since the horizontal part of the trajectory amounts to the majority of the travel time, it is crucial to hit the right level between the two rows.

Table 1 Exponent p for the fitted trend lines $l^p$ of the a posteriori bound (iii) and the numerical results for all four test instances (a)–(d) as depicted in Fig. 5

Table 2 Average ratio between the bounds on $T(\xi _G)-T(\xi _C)$ and the actual differences for all four test cases (a)–(d)

Fig. 6

Shares of the three parts of the a posteriori error bound; from top to bottom: $\int \alpha _0 \Vert \delta \xi \Vert ^2 \hbox {d}\tau$, $\int \alpha _1 \Vert \delta \xi \Vert \Vert \delta \xi _\tau \Vert \hbox {d}\tau$, $\int \alpha _2 \Vert \delta \xi _\tau \Vert ^2 \hbox {d}\tau$. The subfigures a–d refer to the corresponding test instances. Note that on very sparse graphs (large l, right) local minima lead to a dominating distance error term $\int \alpha _0 \Vert \delta \xi \Vert ^2 \hbox {d}\tau$ (dark, top); this is irrelevant for asymptotic considerations

Graph-based shortest paths, which, unsurprisingly, tend to mimic this strategy, are, however, restricted to certain discrete levels. Consequently, the error is sensitive to the exact node positions. If the optimal level is matched by a row of nodes, the error will attain a minimum. On the other hand, if the nodes are positioned such that the optimal trajectory lies exactly between two rows of nodes, we see a maximum error. Obviously, these are nothing more than local deviations from a clear trend.

Finally, it is interesting to notice that in all four test cases the angular error term of the a posteriori bound $\int \alpha _2 \Vert \delta \xi _\tau \Vert ^2 \hbox {d}\tau$ (light blue in Fig. 5) would have been enough to bound the numerical data (blue dots) alone, which lets us conclude that, even though the bound is sharp in the worst case, in particular the average angular error is not perfectly captured.

5 Conclusion

Discretizing the Zermelo navigation problem with a graph-based approach for computing global optima inevitably leads to approximation errors depending on the graph as well as the continuous optimal path. For a certain class of locally densely connected graphs, we have derived three bounds on the excess flight duration in terms of graph and wind properties.

While the local bound improves on the a priori bound by four orders of magnitude, stressing the importance of using localized quantities if possible, it still is far from sharp in numerical examples. The—computationally in general unavailable—a posteriori bound, in contrast, is quite sharp, and thus indicates that the use of a posteriori error estimators providing rough approximations of the actual path error $\delta \xi$ can be expected to improve the bounds further. The observed convergence rates, however, agree well with the computational bounds in both cases.

The error bounds derived here can on the one hand guide the choice of optimal graph structures—the dependence of vertex density h to connectivity length l as presented here is one example—, and on the other hand help identifying switchover points in hybrid discrete-continuous optimization algorithms [7].

Appendix

Recall from (3) the derivative

$$\begin{aligned} f(\xi ,\xi _\tau ) = \frac{-\xi _\tau ^\textrm{T}w(\xi ) + \sqrt{(\xi _\tau ^\textrm{T}w(\xi ))^2+(\overline{v}^2 - \Vert w(\xi )\Vert ^2)\xi _\tau ^2}}{\overline{v}^2 - \Vert w(\xi )\Vert ^2} \end{aligned}$$

of the time parametrization $t(\tau )$. Here, we will compute and bound its second derivative with respect to $\xi$ and $\xi _\tau$ in terms of the wind w and its derivatives.

Theorem A.1

Let $c_0(\xi )=\Vert w(\xi )\Vert <\overline{v}$, $c_1(\xi )=\Vert w_x(\xi )\Vert$, and $c_2(\xi )=\Vert w_{xx}(\xi )\Vert$. Moreover, let $L=\xi _\tau >0$. Then, the second total directional derivative of f is bounded by

$$\begin{aligned} f(\xi ,\xi _\tau )'' [\delta \xi ,\delta \xi _\tau ]^2 \le \alpha _0(\xi ) \Vert \delta \xi \Vert ^2 + \alpha _1(\xi )\Vert \delta \xi \Vert \Vert \delta \xi _\tau \Vert + \alpha _2(\xi ) \Vert \delta \xi _\tau \Vert ^2 \end{aligned}$$

with

$$\begin{aligned} {v}^2&= \overline{v}^2 - c_0^2, \\ \alpha _0&= L \Bigg [ \frac{c_1^2}{\underline{v}^3} \left( 1 +6\frac{c_0}{\underline{v}} +2\frac{\sqrt{\overline{v}^2 + c_0^2}}{\underline{v}} +6\frac{c_0^2}{\underline{v}^2} +8\frac{c_0^3}{\underline{v}^3} +8\frac{c_0^2 \sqrt{\overline{v}^2 + c_0^2}}{\underline{v}^3} \right) \\&\hspace{1cm} + \frac{c_2}{\underline{v}^2} \left( 1 +2\frac{c_0}{\underline{v}} +2\frac{c_0^2}{\underline{v}^2} +2\frac{c_0 \sqrt{\overline{v}^2 + c_0^2}}{\underline{v}^2} \right) \Bigg ], \\ \alpha _1&= \frac{c_1}{\underline{v}^2} \left[ 2 +8\frac{c_0}{\underline{v}} +4\frac{c_0^2}{\underline{v}^2} +8\frac{c_0^3}{\underline{v}^3} \right] , \\ \alpha _2&= \frac{1}{\underline{v}L} \left[ 1 + 3\frac{c_0^2}{\underline{v}^2} \right] . \end{aligned}$$

Proof

The derivative $f=t_\tau$ of parametrized time consists of two terms, the tailwind term

$$\begin{aligned} f_1 = -\frac{\xi _\tau ^\textrm{T} w}{g}, \quad g = \overline{v}^2 - w^\textrm{T} w, \end{aligned}$$

and the length term

$$\begin{aligned} f_2 = g^{-1}\left( (\xi _\tau ^\textrm{T} w)^2 + g (\xi _\tau ^\textrm{T}\xi _\tau )\right) ^{1/2}. \end{aligned}$$

At each time $\tau$, we obtain

$$\begin{aligned} \underline{v}^2:=\overline{v}^2 - c_0^2 \le g \le \overline{v}^2. \end{aligned}$$

The directional derivatives of g in direction $(\delta \xi ,\delta \xi _\tau )$ read

$$\begin{aligned} g' \delta \xi&= -2w^\textrm{T} w_x \delta \xi \quad \Rightarrow \quad \Vert g'\Vert \le 2c_0 c_1 \end{aligned}$$

and, as we are only interested in second order directional derivatives,

$$\begin{aligned} \delta \xi ^\textrm{T} g'' \delta \xi&= - 2\delta \xi ^\textrm{T}( w_x^\textrm{T} w_x+w^\textrm{T}w_{xx}) \delta \xi \quad \Rightarrow \quad \Vert g''\Vert \le 2(c_1^2+c_0c_2). \end{aligned}$$

For the tailwind term, we consider

$$\begin{aligned} f_1' \delta \xi = - g^{-2} \left( (\delta \xi _\tau ^\textrm{T} w + \xi _\tau ^\textrm{T}w_x\delta \xi )g - \xi _\tau ^\textrm{T}w g'\delta \xi \right) . \end{aligned}$$

Again, we are only interested in second directional derivatives and thus consider

Now we turn to $f_2$, first considering the term $F:= (\xi _\tau ^\textrm{T} w)^2 + g (\xi _\tau ^\textrm{T}\xi _\tau )$ with

$$\begin{aligned} \underline{v}^{2} L^2 \le F \le L^2 (\overline{v}^2 + c_0^2). \end{aligned}$$

Then,

$$\begin{aligned} F' \delta \xi = 2\xi _\tau ^\textrm{T}w(\delta \xi _\tau ^\textrm{T}w + \xi _\tau ^\textrm{T}w_x\delta \xi ) + g'\delta \xi \xi _\tau ^\textrm{T}\xi _\tau + 2g\xi _\tau ^\textrm{T} \delta \xi _\tau \end{aligned}$$

and

$$\begin{aligned} F''[\delta \xi ,\delta \xi ]&= 2(\delta \xi _\tau ^\textrm{T}w + \xi _\tau ^\textrm{T}w_x\delta \xi )^2 \\&\quad + 2\xi _\tau ^\textrm{T}w(\delta \xi _\tau ^\textrm{T} w_x \delta \xi + \delta \xi _\tau ^\textrm{T} w_x \delta \xi + \delta \xi ^\textrm{T}(\xi _\tau ^\textrm{T} w_{xx})\delta \xi ) \\&\quad + \delta \xi ^\textrm{T} g'' \delta \xi \xi _\tau ^\textrm{T}\xi _\tau + 2g'\delta \xi \xi _\tau ^\textrm{T}\delta \xi _\tau + 2g'\delta \xi \xi _\tau ^\textrm{T}\delta \xi _\tau + 2g\delta \xi _\tau ^\textrm{T}\delta \xi _\tau \\&= \delta \xi ^\textrm{T} \bigg [2w_x^\textrm{T}\xi _\tau \xi _\tau ^\textrm{T}w_x + 2\xi _\tau ^\textrm{T}w (\xi _\tau ^\textrm{T} w_{xx}) + \xi _\tau ^\textrm{T}\xi _\tau g'' \bigg ] \delta \xi \\&\quad + \delta \xi _\tau ^\textrm{T} \bigg [4w\xi _\tau ^\textrm{T}w_x + 4\xi _\tau ^\textrm{T}w w_x + 4\xi _\tau g' \bigg ]\delta \xi \\&\quad + \delta \xi _\tau ^\textrm{T}\left[ 2ww^\textrm{T} + 2g\right] \delta \xi _\tau . \end{aligned}$$

For $f_2 = g^{-1} \sqrt{F}$, we thus obtain

$$\begin{aligned} f_2' \delta \xi&= -g^{-2} g'\delta \xi F^{1/2} + \frac{1}{2}g^{-1} F^{-1/2} F'\delta \xi . \end{aligned}$$

The second directional derivative is

$$\begin{aligned} f_2''[\delta \xi ,\delta \xi ]&= 2g^{-3} (g'\delta \xi )^{2} F^{1/2} -g^{-2} \delta \xi ^\textrm{T} g''\delta \xi F^{1/2} -g^{-2}g'\delta \xi \frac{1}{2}F^{-1/2}F'\delta \xi \\&\quad - \frac{1}{2}g^{-2}g'\delta \xi F^{-1/2}F'\delta \xi -\frac{1}{4}g^{-1}F^{-3/2}(F'\delta \xi )^2 \\&\quad + \frac{1}{2}g^{-1}F^{-1/2} \delta \xi ^\textrm{T}F''\delta \xi \\&= 2g^{-3} (g'\delta \xi )^{2} F^{1/2} -g^{-2} \delta \xi ^\textrm{T} g''\delta \xi F^{1/2} -g^{-2}g'\delta \xi F^{-1/2}F'\delta \xi \\&\quad -\frac{1}{4}g^{-1}F^{-3/2}(F'\delta \xi )^2 + \frac{1}{2}g^{-1}F^{-1/2} \delta \xi ^\textrm{T}F''\delta \xi . \end{aligned}$$

Adding $f_1''$ and $f_2''$, we finally obtain

$$\begin{aligned} f''(\xi ,\xi _\tau )[\delta \xi ,\delta \xi _\tau ]^2&= (f_1'' + f_2'')[\delta \xi ,\delta \xi _\tau ]^2 \\&= - 2g^{-3} (g'\delta \xi )^2 (\xi _\tau ^\textrm{T}w){} & {} \quad + g^{-2} (\delta \xi ^\textrm{T}g''\delta \xi ) (\xi _\tau ^\textrm{T}w) \\&\quad + 2g^{-2} (g'\delta \xi ) (\xi _\tau ^\textrm{T}w_x\delta \xi ){} & {} \quad - g^{-1} w_{xx}[\xi _\tau ,\delta \xi ,\delta \xi ] \\&\quad + 2g^{-2} (g'\delta \xi ) (\delta \xi _\tau ^\textrm{T} w){} & {} \quad - 2g^{-1} (\delta \xi _\tau ^\textrm{T} w_x \delta \xi ), \\&\quad +2g^{-3} (g'\delta \xi )^2 F^{1/2}{} & {} \quad - g^{-2} (\delta \xi ^\textrm{T} g''\delta \xi ) F^{1/2} \\&\quad - g^{-2} (g'\delta \xi ) F^{-1/2} F'[\delta \xi ,\delta \xi _\tau ]{} & {} \quad + \frac{1}{2} g^{-1} F^{-1/2} F''[\delta \xi ,\delta \xi _\tau ]^2 \\&\quad - \frac{1}{4} g^{-1} F^{-3/2} (F'[\delta \xi ,\delta \xi _\tau ])^2, \end{aligned}$$

which is bounded by

$$\begin{aligned} |f''(\xi ,\xi _\tau )[\delta \xi ,\delta \xi _\tau ]^2|&\le L \Bigg [ \frac{c_1^2}{\underline{v}^3} \left( 1 +6\frac{c_0}{\underline{v}} +2\frac{\sqrt{\overline{v}^2 + c_0^2}}{\underline{v}} +6\frac{c_0^2}{\underline{v}^2} +8\frac{c_0^3}{\underline{v}^3} +8\frac{c_0^2 \sqrt{\overline{v}^2 + c_0^2}}{\underline{v}^3} \right) \\&\hspace{1cm} + \frac{c_2}{\underline{v}^2} \left( 1 +2\frac{c_0}{\underline{v}} +2\frac{c_0^2}{\underline{v}^2} +2\frac{c_0 \sqrt{\overline{v}^2 + c_0^2}}{\underline{v}^2} \right) \Bigg ] \Vert \delta \xi \Vert ^2 \\&\quad +\frac{c_1}{\underline{v}^2} \left[ 2 +8\frac{c_0}{\underline{v}} +4\frac{c_0^2}{\underline{v}^2} +8\frac{c_0^3}{\underline{v}^3} \right] \Vert \delta \xi \Vert \Vert \delta \xi _\tau \Vert \\&\quad +\frac{1}{\underline{v}L} \left[ 1 + 3\frac{c_0^2}{\underline{v}^2} \right] \Vert \delta \xi _\tau \Vert ^2. \end{aligned}$$

$\square$

Since the claim of Theorem A.1 is rather unwieldy, we simplify it, finally proving Lemma 3.1.

Lemma A.1

For any $p\in \varOmega$ let $c_0(p)=\Vert w(p)\Vert$, $c_1(p)=\Vert w_x(p)\Vert$, and $c_2(p)=\Vert w_{xx}(p)\Vert$, and assume $c_0 \le \overline{v}/ \sqrt{5}$. Moreover, let $\xi \in {\hat{X}}$, $L:=\xi _\tau >0$ and ${v}^2(p):= \overline{v}^2 - c_0^2(p)$. Then

$$\begin{aligned} \alpha _0(p)&\le \frac{L}{\underline{v}^3(p)} \left( 12 c_1^2(p) + 4 c_2(p) \underline{v}(p) \right) , \\ \alpha _1(p)&\le \frac{8c_1(p)}{\underline{v}^2(p)},\\ \alpha _2(p)&\le \frac{2}{L\underline{v}(p)} \end{aligned}$$

hold in Theorem A.1.

Proof

Let $s:=c_0/\overline{v}$ be the relative wind speed. Then

$$\begin{aligned}{} & {} \frac{c_0}{\underline{v}} = \frac{s\overline{v}}{\sqrt{\overline{v}^2-s^2\overline{v}^2}} = \frac{s}{\sqrt{1-s^2}} \le \frac{1}{2}, \\{} & {} 1 \le \frac{\sqrt{c_0^2+\overline{v}^2}}{\underline{v}} \le \sqrt{\frac{3}{2}}, \end{aligned}$$

and

$$\begin{aligned} \frac{\overline{v}}{\underline{v}} = \frac{\overline{v}}{\sqrt{\overline{v}^2 - s^2\overline{v}^2}} = \frac{1}{\sqrt{1-s^2}} \le \frac{\sqrt{5}}{2}, \end{aligned}$$

which allows to bound

$$\begin{aligned} \alpha _0&:= L \Bigg [ \frac{c_1^2}{\underline{v}^3} \left( 1 +6\frac{c_0}{\underline{v}} +2\frac{\sqrt{\overline{v}^2 + c_0^2}}{\underline{v}} +6\frac{c_0^2}{\underline{v}^2} +8\frac{c_0^3}{\underline{v}^3} +8\frac{c_0^2 \sqrt{\overline{v}^2 + c_0^2}}{\underline{v}^3} \right) \\&\hspace{1cm} + \frac{c_2}{\underline{v}^2} \left( 1 +2\frac{c_0}{\underline{v}} +2\frac{c_0^2}{\underline{v}^2} +2\frac{c_0 \sqrt{\overline{v}^2 + c_0^2}}{\underline{v}^2} \right) \Bigg ] \\&\quad \le \left[ \frac{c_1^2}{\underline{v}^3} \left( \frac{13}{2} +4\sqrt{\frac{3}{2}}\right) + \frac{c_2}{\underline{v}^2} \left( \frac{5}{2} +\sqrt{\frac{3}{2}}\right) \right] \\&\quad \le L \left( 12 \frac{c_1^2}{\underline{v}^3} + 4 \frac{c_2}{\underline{v}^2} \right) \end{aligned}$$

as well as

$$\begin{aligned} \alpha _1&:= \frac{c_1}{\underline{v}^2} \left[ 2 +8\frac{c_0}{\underline{v}} +4\frac{c_0^2}{\underline{v}^2} +8\frac{c_0^3}{\underline{v}^3} \right] \le \frac{8c_1}{\underline{v}^2} \end{aligned}$$

and

$$\begin{aligned} \alpha _2&:= \frac{1}{\underline{v}L} \left[ 1 + 3\frac{c_0^2}{\underline{v}^2} \right] \le \frac{7}{4\underline{v}L} \le \frac{2}{L\underline{v}}. \end{aligned}$$

$\square$

Lemma A.2

Let $\xi _C \in {\hat{X}}$ be a minimizer of (7). Further, let $\Vert w(p)\Vert \le \overline{c}_0$ and $\Vert w_x(p)\Vert \le \overline{c}_1 ~\forall p\in \varOmega$. Then, it is twice continuously differentiable and its second derivative is bounded by

$$\begin{aligned} \Vert (\xi _C)_{\tau \tau }\Vert \le \overline{\sigma }:= \frac{2 \overline{v}\overline{c}_1 L^2}{(\overline{v}-\overline{c}_0)^3} \left( (1+\sqrt{2})\overline{v}+ \overline{c}_0\right) . \end{aligned}$$

(24)

For $\overline{c}_0 \le \overline{v}/ \sqrt{5}$ this simplifies to $\overline{\sigma }\le 17\frac{\overline{c}_1L^2}{\overline{v}}$.

Proof

The optimal control problem (2) has originally been formulated by Zermelo [24] in terms of the heading angle $\varphi$ in unscaled time t instead of the airspeed v in scaled time $\tau$, which are related by

$$\begin{aligned} v(t) = \overline{v}\begin{bmatrix} \cos \varphi (t) \\ \sin \varphi (t) \end{bmatrix}. \end{aligned}$$

(25)

The Hamiltonian formalism yields an expression for the heading angle rate of an optimal trajectory,

$$\begin{aligned} \varphi _t = w_x: B, \quad B:=\begin{bmatrix} \cos \varphi \sin \varphi &{}\quad - \cos ^2 \varphi \\ \sin ^2\varphi &{}\quad - \cos \varphi \sin \varphi \end{bmatrix}, \end{aligned}$$

with “ : ” denoting tensor contraction, and confirms the regularity of $\xi$. Note that $\Vert B\Vert _F = \sqrt{2}$, where $\Vert \cdot \Vert _F$ denotes the Frobenius norm. By the chain rule, (25) yields

$$\begin{aligned} v_t = \overline{v}\begin{bmatrix} -\sin \varphi \\ \cos \varphi \end{bmatrix} \varphi _t \end{aligned}$$

and the bound

$$\begin{aligned} \Vert v_t\Vert = \overline{v}|\varphi _t| \le \overline{v}\Vert w_x\Vert _F \Vert B\Vert _F \le \sqrt{2}\,\overline{v}\, \overline{c}_1. \end{aligned}$$

(26)

For the ground speed $x_t(t) \underset{(1)}{=} v(t) + w(x(t))$ we thus obtain

$$\begin{aligned} \overline{v}- \overline{c}_0 \le \Vert x_t\Vert \le \overline{v}+ \overline{c}_0 \end{aligned}$$

(27)

and

$$\begin{aligned} \Vert x_{tt}\Vert&\le \Vert v_t + w_x x_t\Vert \underset{\begin{array}{c} (26)\\ (27) \end{array}}{\le } \overline{c}_1 \left( (1+\sqrt{2})\overline{v}+ \overline{c}_0\right) . \end{aligned}$$

(28)

The flight path $\xi _C$ (omitting the subscript C in the following) with constant ground speed $\Vert \xi _\tau \Vert = L$ is related to the actual flight path x by $\xi (\tau ) = x(t(\tau ))$ with $t:[0,T]\rightarrow [0,1]$ being a monotone bijection. Therefore,

$$\begin{aligned} \xi _\tau&= x_t(t) t_\tau \end{aligned}$$

(29)

yields

$$\begin{aligned} L = \Vert \xi _\tau \Vert = \Vert x_t\Vert \, t_\tau \underset{(27)}{\ge } (\overline{v}- \overline{c}_0) t_\tau \quad \Rightarrow \quad t_\tau \le \frac{L}{\overline{v}- \overline{c}_0}. \end{aligned}$$

(30)

For the second derivative, we note that

$$\begin{aligned} 0 = (L^2)_\tau = (\xi _\tau ^\textrm{T} \xi _\tau )_\tau = 2 \xi _{\tau \tau }^\textrm{T} \xi _\tau , \end{aligned}$$

which means that the curvature vector $\xi _{\tau \tau }$ is orthogonal to the path and the ground velocity $\xi _\tau$ (and $x_t$). Consequently, we obtain

$$\begin{aligned} 0&= (\xi _\tau ^\textrm{T} \xi _\tau )_\tau \underset{(29)}{=} (\xi _\tau ^\textrm{T} x_t t_\tau )_\tau = \underbrace{\xi _{\tau \tau }^\textrm{T} x_\tau }_{=0} t_\tau + \xi _\tau ^\textrm{T} x_ {\tau \tau } t_\tau ^2 + \xi _\tau ^\textrm{T} x_\tau t_{\tau \tau } \underset{(29)}{=} \xi _\tau ^\textrm{T}x_{\tau \tau } t_\tau ^2 + \frac{L^2}{t_\tau } t_{\tau \tau } \end{aligned}$$

and therefore

$$\begin{aligned} |t_{\tau \tau }|&\le \frac{t_\tau }{L^2} \Vert \xi _\tau \Vert \, \Vert x_{\tau \tau }\Vert \, t_\tau ^2 \underset{(28)}{\le } \frac{\overline{c}_1 t_\tau ^3}{L} \left( (1+\sqrt{2})\overline{v}+ \overline{c}_0\right) . \end{aligned}$$

(31)

Now we can bound

$$\begin{aligned} \Vert \xi _{\tau \tau }\Vert&\underset{(29)}{=} \Vert x_{tt} t_\tau ^2 + x_t t_{\tau \tau } \Vert \\&\le \Vert x_{tt} \Vert \, t_\tau ^2 + \Vert x_t\Vert \, |t_{\tau \tau }| \\&\underset{\begin{array}{c} (27)\\ (28)\\ (31) \end{array}}{\le } \overline{c}_1 \left( (1+\sqrt{2})\overline{v}+ \overline{c}_0\right) t_\tau ^2 + (\overline{v}+ \overline{c}_0) \frac{\overline{c}_1 t_\tau ^3}{L} \left( (1+\sqrt{2})\overline{v}+ \overline{c}_0\right) \\&= \overline{c}_1 t_\tau ^2 \left( (1+\sqrt{2})\overline{v}+ \overline{c}_0\right) \left( 1 + (\overline{v}+ \overline{c}_0) \frac{t_\tau }{L} \right) \\&\underset{(30)}{\le } \frac{\overline{c}_1 L^2}{(\overline{v}-\overline{c}_0)^2} \left( (1+\sqrt{2})\overline{v}+ \overline{c}_0\right) \left( 1 + \frac{\overline{v}+ \overline{c}_0}{\overline{v}- \overline{c}_0} \right) \\&= \frac{2 \overline{v}\overline{c}_1 L^2}{(\overline{v}-\overline{c}_0)^3} \left( (1+\sqrt{2})\overline{v}+ \overline{c}_0\right) , \end{aligned}$$

which completes the proof. $\square$