A Hamilton-Jacobi Formulation for Optimal Coordination of Heterogeneous Multiple Vehicle Systems

03/12/2020 ∙ by Matthew R. Kirchner, et al. ∙ The Regents of the University of California U.S. Navy 0

We present a method for optimal coordination of multiple vehicle teams when multiple endpoint configurations are equally desirable, such as seen in the autonomous assembly of formation flight. The individual vehicles' positions in the formation are not assigned a priori and a key challenge is to find the optimal configuration assignment along with the optimal control and trajectory. Commonly, assignment and trajectory planning problems are solved separately. We introduce a new multi-vehicle coordination paradigm, where the optimal goal assignment and optimal vehicle trajectories are found simultaneously from a viscosity solution of a single Hamilton-Jacobi (HJ) partial differential equation (PDE), which provides a necessary and sufficient condition for global optimality. Intrinsic in this approach is that individual vehicle dynamic models need not be the same, and therefore can be applied to heterogeneous systems. Numerical methods to solve the HJ equation have historically relied on a discrete grid of the solution space and exhibits exponential scaling with system dimension, preventing their applicability to multiple vehicle systems. By utilizing a generalization of the Hopf formula, we avoid the use of grids and present a method that exhibits polynomial scaling in the number of vehicles.



There are no comments yet.


page 1

This week in AI

Get the week's most popular data science and artificial intelligence research sent straight to your inbox every Saturday.

I Introduction

Multi-robot motion planning presents several challenges, and key among them are trajectory planning and formation control, as well as the assignment of vehicles to goal states. When the assignment of vehicles to goal states is made beforehand, many methods have been proposed. These include control theory approaches [40, 11], graph-based techniques [14], and optimization methods [27, 2]. Approximations to assist in dimensionality reductions have also been proposed, such as sequential methods [8] and hybrid approaches [10, 37].

Fig. 1: An illustration of the vehicle coordination problem. If a vehicle’s position in the formation is not assigned a priori, then many different configurations are possible to achieve the desired formation, each with a differing system objective cost. Two such possibilities are shown, one in red and an alternate assignment in blue.

When the goal states of multi-vehicle teams are not set a priori, it is necessary to determine which vehicle to allocate to which terminal goal. There is a large body of literature developing allocation and assignment algorithms [5], and applying these methods to engineering problems has appeared in diverse forms including sensor coverage [26], weapon-target assignment [1], or network routing [24]. Dynamics are not considered in this body of work, which commonly study worst-case system performance under equilibrium conditions using a game-theoretical framework [38, 25, 21]. Without consideration of dynamics, these methods have limited applicability to vehicle coordination problems and there exists a technical gap relating vehicle-goal allocation and motion planning, as most prior work assumes each is done independently. The critical dependency between assignment and trajectory planning is illustrated in Figure 1.

There have been recent attempts to close this gap as in [17]. But this work only evaluates discrete spatial waypoints on a rectangular grid without considering vehicle dynamics, and therefore the paths generated are not feasible for most physical systems. Morgan et al. [31]

generalizes this concept to linear, discrete-time dynamics and constructs a non-linear programming problem that is solved with sequential convex programming in conjunction with an auction algorithm

[5]. However, both [31] and [17] restrict the cost function to one of additive individual vehicle weights (similar to the cost function presented in [5]), which excludes many vehicle optimization problems of interest, including assembling the desired formation shape in minimal time.

We present how the optimal assignment and trajectory can be found simultaneously from the viscosity solution to a single Hamilton–Jacobi PDE, a necessary and sufficient condition for optimality. We assume the vehicles can each have unique linear dynamics. Remarkably, we show that the global viscosity solution can be found as the solution to a linear bottleneck assignment problem (LBAP) [15]. This fact can be utilized to construct a level set method that has polynomial computational scaling with the number of vehicles and still ensures a global optimum.

Our formulation has a close relation to reachability analysis which can be useful for the design and analysis of heterogeneous systems. When the dynamics of some vehicles differ greatly from the rest, some formations may not be feasible (e.g. a slow-moving vehicle attempting to join a formation with a much faster vehicle). The zero sub-level set of the viscosity solution of a related HJ PDE defines an implicit surface representation of the backwards reachable set [28], and from this one determines whether the formation can be achieved given the unique collection of vehicle dynamics.

Traditionally, numerical solutions to HJ equations require a dense, discrete grid of the solution space [33, 30]

. Computing the elements of this grid scales poorly with dimension and has limited use for problems with dimension greater than four. The exponential dimensional scaling in optimization is sometimes referred to as the “curse of dimensionality

[4, 3]. A new result in [9] discovered a numerical solution based on the Hopf formula [18] that does not require a grid and can be used to efficiently compute solutions of a certain class of Hamilton–Jacobi PDEs. However, that only applies to systems with time-independent Hamiltonians of the form , and has limited use for general linear control problems. Recently, the classes of systems were expanded upon and generalizations of the Hopf formula are used to solve optimal linear control problems in high-dimensions [19] and differential games as applied to multi-vehicle, collaborative pursuit-evasion problems [20].

The key contributions of this paper are to formulate the multi-vehicle coordination problem as a single Hamilton–Jacobi equation, show that the solution of which can be found from an LBAP, and use the generalized Hopf formula to simultaneously solve for optimal trajectory, vehicle control, and goal assignment. By utilizing the generalized Hopf formula, we create a level-set method to find the viscosity solution of the HJ in a computationally efficient manner.

We introduce a multi-vehicle dynamical system and formulate the coordination problem in Section II. Section III formulates the optimal vehicle coordination problem as the viscosity solution to a single Hamilton–Jacobi PDE. Section IV constructs a level-set method for the multi-vehicle coordination problem and utilizes the generalized Hopf formula as a computational tool to find the viscosity solution to the HJ PDE of the joint, mutli-vehicle system. Section V shows the method as applied to several examples including a case of planar motion planning with 4 vehicles. This example shows how assignment and path planning cannot be decoupled and even if a single vehicle has a different initial condition, it may result in a different assignment for optimal coordination.

Ii Problem Formulation

We consider a system that consists of vehicles and each vehicle, , has linear dynamics


for , where is the system state and is the control input, constrained to a convex admissible control set

. We add the assumption that no eigenvalue of matrix

has a strictly positive real part. We let denote a state trajectory for vehicle that evolves in time with measurable control sequence , according to starting from initial state at . The trajectory is a solution of in that it satisfies almost everywhere:

Ii-a Multi-Vehicle Model

For the set of

vehicles, we construct a joint state space with state vector

and control which is written as follows


We reiterate the above definition of the joint state space with the fact that when quantities , , , , and appear without subscripts, they refer to the joint system in , and when subscripts are used, they refer to a specific individual vehicle as defined in .

Ii-B Vehicle Coordination

We assume there exists a set of goals, and for each vehicle we associate closed convex sets with the understanding that means the vehicle is assigned to goal . Our goal is to make sure that we have one vehicle at each goal, but it does not matter which vehicle is at each goal. This goal can be expressed by the requirements that the multi-vehicle state, , belongs to the following goal set


where denotes the set of all permutations of . We represent implicitly with the function such that


and use it to construct a cost functional for the system trajectory , given terminal time as


where the function represents the rate that cost is accrued over time. The value function is defined as the minimum cost, , among all admissible controls for a given initial state with


Iii Hamilton–Jacobi Equations For Optimal Coordination

The value function in satisfies the dynamic programming principle [6, 12] and also satisfies the following initial value Hamilton–Jacobi (HJ) equation with being the viscosity solution of


for , where the Hamiltonian is defined by


where in denotes the costate. We denote by the costate trajectory that can be shown to satisfy almost everywhere:

with initial costate denoted by . With a slight abuse of notation, we will hereafter use to denote , when the initial state and control sequence can be inferred through context with the corresponding state trajectory, .

Iii-a System Hamiltonian

We consider a time-optimal formulation with


is the characteristic function of the set of admissible controls and is defined by

In this case, the integral term in disappears (as long as remains in and is simply the value of at the terminal time, . For this problem

means that there exists a feasible trajectory for the vehicles the ends at the state ; whereas means that such trajectories do not exist under the system dynamics and initial conditions. Since each vehicle has independent dynamics, the Hamiltonian in is of the form


where each vehicle’s Hamiltonian is given by


Iii-B Linear Bottleneck Assignment

Our goal is solving the Hamilton–Jacobi PDE in , and any that satisfies is, in general, non-convex and presents numerical challenges. We show that we can overcome this challenge by alternatively solving for the global value function with the following linear bottleneck assignment problem [15]:


where is the viscosity solution to


with Hamiltonian defined by . The function is an implicit surface representation of such that


The solution to can be found from using the appropriate linear bottleneck assignment algorithms (for example [7, Section 6.2]), requiring only evaluations of lower dimensional viscosity solutions and avoiding computation involving all permutations of vehicle-goal pairs. Also of note is that if each HJ equation has convex initial data, it enables the use of computationally efficient techniques that can guarantee convergence to the appropriate viscosity solution. Additionally, since each of the solutions are independent, the solutions can be computed in parallel.

We introduce a set of mild regularity assumptions that guarantee Hamilton–Jacobi equations have a unique viscosity solution [41, Chapter 7, p. 63]:

  1. Each Hamiltonian

    is continuous.

  2. There exists a constant such that for all and for all , the following inequalities hold


    with .

  3. For any compact set there exists a constant such that for all and for all the inequality holds,

    with .

  4. Each terminal cost function

    is continuous.

Lemma 1

If each vehicle Hamiltonian, , meets assumptions 1-3, then the Hamiltonian defined in also meets assumptions 1-3.

The proof is given in the appendix.

Theorem 2

Under assumptions 1-4, is a unique viscosity solution to for all , and there exists a that satisfies such that with the Hamiltonian given by , is a viscosity solution to .

We will prove the theorem constructively by proposing a particular given as




where is the implicit representation of , defined in , with . We see that for any when which implies . We also see if there exists an such that , then and implies . Therefore the proposed in satisfies and is an implicit surface representation of the set . Furthermore, since by assumption each is continuous, then is convex as the max of a finite number of continuous functions is also continuous. This implies is continuous as the minimum of a finite number of coninuous functions is continuous. Therefore satisfies assumption 4.

From Lemma 1, the system Hamiltonian defined in meets assumptions 1-3 and implies that has a unique viscosity solution, denoted as , when the initial cost function is given by [41, Theorem 8.1, p. 70]. The uniqueness of the solution implies that the viscosity solution is equivalent to the value function in . It follows that


since when . Denoting by as the control that optimizes , we have

Substituting and we have


where we use the notation to represent the -th block of the vector . Likewise has a unique viscosity solution,, for each with initial cost given , and as such


for each and with denoting the control that optimizes . Note that if it would contradict that is the optimal control. Also if , then it would contradict that is the optimal control of the entire system. Therefore, and our value function in becomes

and we arrive at our result.

Iv A Level Set Method with the Generalized Hopf Formula

First, we introduce the Fenchel–Legendre transform of a convex, proper, lower semicontinuous function , denoted as , and is defined as [16]


We propose to find the viscosity solutions of using the generalized Hopf formula.

Theorem 3

Each viscosity solution of can be found from the formula




Proceeding similar to [20], we apply a change of variables to with


which results in the following system


By utilizing this change of variables, we construct an equivalent HJ equation


From Lemma 5 (given in the appendix), meets assumption 1-3, and by composition rule is continuous [39, Theorem 4.7] and meets assumption 4. Therefore has a unique viscosity solution that is equivalent to the cost functional

where is a solution trajectory that satisfies almost everywhere. Since

we have

noting that since by the transform implies at . Thus, and we can find by finding the viscosity solution to .

Since is assumed closed and convex, this implies is convex and by assumption 4 is continuous in . By assumption 1, is continuous in , therefore gives an exact, point-wise viscosity solution for [23, Section 5.3.2, p. 215]. Formula shows that we can compute a viscosity solution to by solving a finite dimensional optimization problem. This avoids constructing a discrete spatial grid, and is numerically efficient to compute even when the dimension of the state space is large. Additionally, no spatial derivative approximations are needed with Hopf formula-based methods, and this eliminates the numeric dissipation introduced with the Lax–Friedrichs scheme [34], which is necessary to maintain numeric stability in grid-based methods.

Iv-a Numerical Optimization of the Hopf Formula

We transcribe the Hopf formula into a non-linear programming problem by approximating the integral in with an -point quadrature rule sampled on the time grid

with and corresponding weights . Additionally, we make a simple, invertible change of variable and substituting into results in the following unconstrained optimization problem that solves :


with now being defined as


The variable transformation is done so that when the matrix has at least one eigenvalue with strictly negative real part, we avoid computation of which would have exponentially unstable poles since we are evaluating the matrix exponential of . This divergence would cause sensitivity in the evaluation of the Hopf formula in with respect to small changes of . By utilizing the variable transformation and optimizing , we avoid this sensitivity.

Remark 4

Often the expression in is known in closed form and we present one such common case. Recall that denotes the Fenchel–Legendre transform of a function defined in , and suppose is the closed convex set defined as

where is any norm. Then defines a norm, which we denote with , which is the dual norm to [16]. From this we write as


Iv-B Time-Optimal Control to a Goal Set

The task of determining the control that drives the system into in minimal time can be determined by finding the smallest such that


When the system is constrained controllable, then the set of times such that is reachable with respect to contains the open interval [13]. This insures that if is outside the set at time , then and there exists a time such that for all . This implies standard root-finding algorithms can be used to find . As noted in [19], we solve for the minimum time to reach the set by constructing a newton iteration, starting from an initial guess, , with


where is the solution to at time . The value function must satisfy the HJ equation

where and each is the argument of the minimizer in . We iterate until convergence at the optimal time to reach, which we denote as . This process is summarized in Algorithm 1.

The optimal control and trajectory for each vehicle is found directly from the necessary conditions of optimality established by Pontryagin’s principal [35] by noting the optimal trajectory must satisfy

where is the optimal costate trajectory and is given by

This implies our optimal control is


for all time , provided the gradient exists.

4:      , ,
5:while  do
6:     for all  do
8:     end for
12:end while
14:      ,
Algorithm 1 Algorithm to compute the viscosity solution at the minimum time to reach of the set .

V Results

We present results for two examples. For the first we choose a toy example of two vehicles, each with differing single integrator dynamics. The purpose of this example is that it allows easy visualization of the level set propagation and the solution can be verified by comparing with methods such as [29], which cannot be applied to the higher dimensional second example to follow. The LBAP was solved using [7, Algorithm 6.1] and the optimization used sequential quadratic programming [32, Chapter 18] with initial conditions for each vehicle. To prevent a singular control condition, a slight smoothing was applied to the Hamiltonian in and using [36] with parameter .

V-a Toy Problem

Fig. 2: The two vehicle example of Section V-A. The zero level set evolution solved at various times shown in green.

We present a two vehicle problem where the dynamics are for the first vehicle and for the second vehicle. The states are the linear position of each respective vehicle. The control is bounded by . The goal states are for and for . The level set contours of the joint space are shown in Figure V-A, for ten time samples equally spaced on . Notice this seemingly simple heterogeneous system can give rise to counterexamples for the assignment formulations given in [31] where the sum of the vehicle distances are used as the assignment metric and [17] where the time of arrival of vehicles was used. Take the initial state , which is shown on Figure 2. For the time metric, the the assignment gives , while the assignment gives , indicating is the clear choice. However, is the global optimum in this example, as both vehicles reach their desired goal states in a time of , as the zero level set of the global value function intersects the point at that time, giving the minimum time-to-reach. For the distance metric, the misassignment is more pronounced, with the assignment giving as compared to the assignment gives .

V-B Planar Motion

(a) Optimal trajectories at (10% of path).
(b) Complete optimal trajectories.
(c) Complete optimal trajectories for an alternate initial condition.
Fig. 3: The example presented in Section V-B with four vehicles guiding to four possible goal sets. The optimal formation is achieved in . The trajectories are shown for two different times as they traverse the optimal path in the left two figures. The figure on the right is the trajectories when the initial condition for vehicle 4 is different.

We choose for dynamics with state , where is spatial position of a robot and is the velocity and

for each vehicle . The control is constrained to lie in the set . The robots are tasked with reaching the goal formation and coming to rest, in minimum time. The goal sets each have radius of and the centers of the goals are located spatially at , , , and . Since the 2-norm is self-dual, the Hamiltonian for each vehicle is

Since contains eigenvalues with negative real part, we optimize using the variable transformation of Section IV-A. The initial conditions of the vehicles were

Only iterations of were needed to solve for the minimum time to reach the goal formation, which was found to be . Figures (a)a and (b)b show the optimal paths found from . Figure (c)c shows the optimal paths when the initial condition for vehicle 4 was changed to and a different optimal assigment results.

Vi Conclusions and Future Work

We presented how to formulate vehicle coordination problems with unknown goal assignments as the viscosity solution to a single Hamilton–Jacobi equation. We show the solution of this single HJ PDE is equivalent to decomposing the problem and performing a linear bottleneck assignment using the viscosity solutions of independent single-vehicle problems. This allows quadratic computational scaling in the number of vehicles. Finally, a level set method based on the Hopf formula was presented for efficient computation of the vehicle value functions, in which each can be computed in parallel. The Hamilton–Jacobi formulation presented has other advantages for multi-robot systems, such as the compensation of time delays which can be induced in several ways including computation, sensing, and inter-robot communication. See for example [22].

Future work includes to expand the class of allowable vehicle dynamics from general linear models to certain types of nonlinear dynamics and to allow for dependencies between vehicles both in the dynamics models and in the cost functional.


Lemma 5

Given that satisfies assumptions 1-3, then also satisfies assumption 1-3.

From assumption 1, is continuous by composition rule [39, Theorem 4.7] and since does not depend on state, it therefore trivially meets assumption 3. It remains to show that meets assumption 2. We know assumption 2 holds for all and for all , therefore the following holds

When we have

for any . And since for all , it follows that assumption 2 is met.

Proof of Lemma 1. The sum of continuous functions is also continuous [39, Theorem 4.9], therefore assumption 1 is met. To prove assumption 2, we first write, for all ,