Finite-horizon, energy-efficient trajectories in unsteady flows.

Intelligent mobile sensors, such as uninhabited aerial or underwater vehicles, are becoming prevalent in environmental sensing and monitoring applications. These active sensing platforms operate in unsteady fluid flows, including windy urban environments, hurricanes and ocean currents. Often constrained in their actuation capabilities, the dynamics of these mobile sensors depend strongly on the background flow, making their deployment and control particularly challenging. Therefore, efficient trajectory planning with partial knowledge about the background flow is essential for teams of mobile sensors to adaptively sense and monitor their environments. In this work, we investigate the use of finite-horizon model predictive control (MPC) for the energy-efficient trajectory planning of an active mobile sensor in an unsteady fluid flow field. We uncover connections between trajectories optimized over a finite-time horizon and finite-time Lyapunov exponents of the background flow, confirming that energy-efficient trajectories exploit invariant coherent structures in the flow. We demonstrate our findings on the unsteady double gyre vector field, which is a canonical model for chaotic mixing in the ocean. We present an exhaustive search through critical MPC parameters including the prediction horizon, maximum sensor actuation, and relative penalty on the accumulated state error and actuation effort. We find that even relatively short prediction horizons can often yield energy-efficient trajectories. We also explore these connections on a three-dimensional flow and ocean flow data from the Gulf of Mexico. These results are promising for the adaptive planning of energy-efficient trajectories for swarms of mobile sensors in distributed sensing and monitoring.

Introduction

The ability to generate energy-efficient trajectories that take advantage of the inherent motions of a background flow field has significant implications for monitoring large bodies of water with intelligent mobile sensors [1-3], furthering our understanding of the climate and natural ecosystems [4-6]. Developments in this area also present economic opportunities for cost reduction in industries that rely heavily on maritime transport and shipping. Self-powered mobile sensors typically have complex performance tradeoffs, limiting size, weight and power (SWAP). Furthermore, most mobile sensors will only have partial and imperfect information about the ambient flow field, resulting in a finite-horizon predictive window to make decisions about its trajectory. Improving the generation of energy-efficient trajectories that intelligently leverage the flow field to go with the flow may have significant benefits in extending the duration and reach of these mobile sensing platforms. This work provides an extensive analysis of trajectories generated through a finite-horizon model predictive control (MPC) optimization of a mobile sensor in a time-varying background flow across a wide range of system parameters. Furthermore, we establish connections between the control performance and efficiency with the alignment of these trajectories along coherent structures in the background flow.

Currently, there exists an extensive literature that has investigated various algorithms for trajectory generation for such transport applications. For example, graph search algorithms and stochastic optimization have been investigated for path planning [7-9]. Assimilating in situ observations obtained by mobile sensors in an adaptive fashion into ocean models has also been explored, for example with mixed integer programming algorithms [10,11]. Coordinated control of ocean gliders for adaptive ocean sensing has been exhaustively studied in Monterey Bay [12-15]. Algorithms inspired from computational fluid dynamics have also been used to explore coordinated control of swarms in flow fields [16-20]. More recent developments include the use of reinforcement learning algorithms to find optimal paths in turbulent flow fields [21-23]. However, there has been relatively little work in developing a deep understanding of the connection between the dynamics of the flow field and the nature of the optimal trajectories within the flow fields, with a few notable exceptions [24-27]. A key challenge in exploring this connection is the complexity of fluid flow fields, which typically involve the existence of multiple scales in space and time.

To understand the complexity of fluids, techniques from dynamical systems are often employed. Lagrangian coherent structures (LCS) have emerged as a robust and principled approach to uncover invariant manifolds that mediate the transport of material in unsteady fluid flows [28-34]. Specifically, LCS define the transport barriers in a flow field where passive drifters are attracted to or repelled by. There has been considerable work in the development of algorithms to accurately and efficiently compute these structures from data [30,33,35-42]. The finite-time Lyapunov exponent (FTLE) is a scalar value that characterizes the divergence from a trajectory over a finite time interval, and is often used to compute LCS. The FTLE method has been successfully applied to domains of bio-propulsion [43], medicine [44,45], the spread of microbes [46] and the study of aerodynamics [47,48].

The ideas from both trajectory generation and the theory of LCS have been related in the past [24-27]. A predecessor of this was the planning of space missions using invariant manifolds [49]. In the context of ocean transport, Inanc et al. [24] showed that the optimal trajectories of autonomous agents generated using a receding-horizon optimal control algorithm overlap with LCS. Moreover, Senatore & Ross [26] exploited this idea further to generate energy optimal paths by controlling the agents to track the background LCS. Recent papers have further explored the connections between optimal control and LCS [50-52] in the context of path planning in the ocean. However, there is still a need to better understand how the prediction horizon and relative cost of actuation in the autonomous agent optimization relate to the use of coherent structures in the unsteady background flow.

In this work, we investigate the explicit connections between finite-horizon energy-optimal trajectories of a mobile sensor and the underlying background flow dynamics. We specifically analyse how key parameters of the MPC-based optimization affect how the resulting autonomous agent trajectory uses unsteady fluid coherent structures for energy-efficient transport. This analysis is performed primarily on the double gyre flow field, which is a testbed to understand mixing and transport in the ocean. Subsequently, we also further verify our observations on advanced flow fields including the Arnold–Beltrami–Childress (ABC) flow, a three-dimensional incompressible analytical flow field, and unsteady flow data from the Gulf of Mexico. A summary of our methodology is shown in figure 1. The choice of MPC is particularly relevant in this work, as both the FTLE and MPC rely on finite-time horizons in their computations. To explore this connection, we perform an exhaustive search through several of the trajectory optimization parameters that are important to practitioners, including the prediction time horizon and step size for MPC, the relative cost of actuation versus state tracking error, and the maximum agent velocity. We find that there are strong correlations between the presence of background FTLE ridges and the actuation energy expenditure at the corresponding locations along the trajectory.

Figure 1

Overview of the proposed methodology for analysing the connections between finite-horizon energy optimal trajectories and the FTLE field. A self-propelling agent is controlled to transit from a starting location to a goal location through a finite-time horizon energy-optimal trajectory in a time-varying double gyre flow field. The resulting agent trajectory, along with the finite-horizon predicted trajectories at each time step, are shown and colour-coded based on instantaneous energy expenditure (top left). The trajectory history (solid) and the future forecast trajectory bundle (dashed) at an example time instant are shown (top right); the instantaneous FTLE ridges are also shown below these with blue indicating the repelling LCS and red indicating the attracting LCS. As can be observed from the snapshots taken at four particular times (bottom), the energy expenditure along the planned trajectory, given by the colour of the dashed line, and the shape of the finite-horizon trajectory depend on the evolution of the local FTLE ridges. (Online version in colour.)

The remainder of this work is organized as follows. In §2, the core methodologies of MPC and FTLE are discussed. MPC will be the primary optimization algorithm used to generate trajectories, and these will be analysed using FTLE fields. Section 3 describes models for the mobile sensor dynamics and actuation, along with the dynamics of the unsteady double gyre background flow field. The main results are presented in §4, including in-depth analysis of trajectories generated across a wide range of system parameters. In particular, the time horizon of the MPC optimization, the relative cost of actuation versus state tracking error, and the frequency of the background flow oscillation are all investigated. In §5, we highlight the use of MPC on other flow fields and demonstrate our observations on the ABC and Gulf of Mexico flow fields. Section 6 provides a summary of results and a discussion of limitations with suggestions for future work. Appendix A also provides additional plots and analysis of the data that was not presented in the main text.

Methodology

In this section, we introduce two approaches for analysing and generating fuel-efficient trajectories for a mobile sensor in an unsteady background flow: FTLE fields and MPC. First, we introduce the computation of FTLE fields [28,30,33] for passive tracer particles to extract LCS from a time-varying flow field. This method is particularly important to characterize the uncontrolled behaviour of drifters in terms of finite-time attraction and repulsion behaviours. Next, we introduce the preliminaries of finite-horizon MPC, which is an online control optimization algorithm that optimizes a cost function defined over a finite-time prediction horizon. We will use MPC for trajectory optimization of a mobile sensor in an unsteady background flow. MPC is a natural choice, since the mobile sensor will have limited actuation authority, and information about the flow field will only be approximate as it is limited to a finite-time horizon.

Finite-time Lyapunov exponents

Given a vector field , the dynamics of a passive drifter is given by Here, represents time, and is the position of the drifter, where typically or , depending on the dimension of space. The FTLE field can be used to determine the LCS of an unsteady vector field [28,30,33]. The LCS are curves or surfaces in the domain where nearby trajectories are strongly attracted to or repelled from, making them time-varying analogues of stable and unstable invariant manifolds in dynamical systems theory [53].

The FTLE algorithm is as follows. First, a grid of drifters is initialized at time and numerically integrated through the flow field for a fixed amount of time (i.e. the time horizon) , resulting in a flow map : The flow map operator takes each drifter at an initial condition and returns its new position after it is advected through the vector field for a time . Next, the Jacobian matrix of partial derivatives of the flow map, , is computed using finite differences for each drifter in the grid, represented by the grid node indices , such that where are the horizontal and vertical components of the position vector . This flow map Jacobian is used to compute the Cauchy–Green deformation tensor, given by where represents the matrix transpose, not to be confused with the duration of integration . Finally, the largest eigenvalue of for each drifter is used to compute the FTLE field: Alternatively, can be computed as the largest singular value from the singular value decomposition of . It is important to note that for unsteady flow fields, the FTLE field will also vary in time, so that at each new time step a new grid of particles must be reinitialized and advected through the flow. This procedure is typically quite expensive, although there are algorithms to eliminate redundant calculations [35,36].

LCS are often computed as ridges of the FTLE field, which requires an additional step of computing the Hessian of for ridge extraction. FTLE based on drifters integrated forward in time, , results in coherent structures that repel drifters. Similarly, FTLE based on drifters integrated backward in time, , results in coherent structures that attract drifters. These can be seen in figure 1 as red and blue curves, where the red curves are attracting and the blue are repelling. FTLE fields and the resulting LCS are related to almost invariant sets from statistical dynamical systems [54-57]. In particular, LCS act as separatrices in the flow, segmenting different regions where passive tracers remain trapped [58]. FTLE and LCS have also been used extensively to analyse ocean flows [59-61], for example to model the spread of pollution [62]. More broadly, FTLE has also been used to compute coherent structures for a wide range of other flows [63-68]. In this work, we will use FTLE fields generated from passive particles to investigate the trajectories of active mobile sensors, to understand how and when these sensors exploit structures in the flow field for energy-efficient transport.

Model predictive control

The dynamics of mobile sensors operating in real environments are often strongly nonlinear and subject to hardware constraints, time delays, non-minimum phase dynamics, instability, and restrictions on actuation capability. These limitations make the use of traditional linear control approaches challenging, motivating the powerful MPC optimization [69-72] described here. In this work, we use MPC to generate trajectories for a mobile sensor in an unsteady background flow and investigate how these trajectories vary with the optimization parameters.

In general, the dynamics of a nonlinear system with actuation can be written as where is the controlled vector field. In the context of this paper, the state can be either the position of the agent, as in the previous section, or both the position and velocity of the agent.

MPC is a powerful method for calculating the actuation by formulating an iterative optimization problem that minimizes a cost function over a finite-time horizon. The controller enacts this optimized actuation policy for a short time, often for a single time step, and then the optimization problem is recomputed and initialized at the current state. In this way, MPC is quite robust to model uncertainty and disturbances, as the optimization is continuously being reinitialized as new information is available about how the system actually responds to the actuation. Computing over a finite-time horizon might also make MPC more flexible and faster than a global optimization technique, especially for chaotic systems, which may result in stiff long-time optimizations. These benefits make MPC more versatile and widely used over other traditional trajectory generation algorithms. Finally, the FTLE and MPC computations are both performed over a finite time horizon, suggesting the potential for a connection between the outputs of the two algorithms.

Typically, the optimization cost for MPC can be formulated as subject to the system dynamics in (2.6) and control constraints imposed by physical limitations: Here, and are the minimum and maximum values the components of can take, respectively. For example, the actuators may be unable to produce thrusts beyond a certain value. The state error is given by . The finite-time horizon over which we forecast our model for the optimization is ; this term is similar to , the advection time used to calculate FTLE. is a positive definite matrix that quantifies the penalty on actuation effort, and and are positive semi-definite matrices that quantify the penalty on deviations of the state from the goal throughout the trajectory and at the final time step, respectively. For computational purposes, (2.7) is often discretized. The sampling time step is , the discretization of . In this paper, the serves as the ‘replanning frequency’, or, the time after which the a re-optimization of a trajectory is carried out. It is possible to improve the computational speed and convergence of the algorithm with a warm start, which uses the trajectory computed in a previous instance as the initial guess for the trajectory in the next instance [73].

Model problem

We now discuss the models used to simulate the agent dynamics and the unsteady flow field the mobile sensor operates within. We also provide specific parameters that are used for all numerical experiments.

Sensor dynamics

In a two-dimensional setting, a simple kinematic model for the dynamics of the mobile sensor is given by adding the velocity due to actuation, , to the background flow velocity : The state is the position vector. The key assumption in this model is that, without control, the velocity of the sensor, , matches the velocity of the background fluid flow. Thus, the uncontrolled mobile sensor can be considered as a passive Lagrangian drifter, and (3.1) degenerates to (2.1) when all times. Moreover, it assumes that the sensor can generate its own relative velocity in addition to the flow-induced velocity. It is possible to develop more sophisticated models for the mobile sensor dynamics that include inertial and rotational dynamics; in Zhang et al. [25], it was shown that trajectories based on such models also show strong correlation with the presence of background LCS.

Double gyre flow field

We will investigate the motion of the mobile sensor above in the unsteady double gyre flow field described here. The double gyre flow is an analytically defined, periodic vector field that is often used to study mixing and coherent structures related to those found in geophysical circulations. In particular, the double gyre represents a typical large-scale ocean circulation phenomenon often observed in the northern mid-latitude ocean basins. This circulation is quite dominant and is persistent, consisting of sub-polar and sub-tropical gyres. As a major type of ocean circulation, several main features of the double gyre phenomena have been identified through analysing observational data and numerical simulations [74-76].

The double gyre velocity field is derived from the stream function where the time dependency is introduced by with time-dependent coefficients This flow is defined on a non-dimensionalized domain of , where , , , , , . Here, dictates the magnitude of oscillation in the -direction, is the angular oscillation frequency, and controls the velocity magnitude. Unless stated otherwise, the parameters used for the double gyre flow field are as in Shadden et al. [30], where , and . The resulting velocity field is given by

Specific control objective

By combining the mobile sensor model and the double gyre flow field, the dynamics of the sensor are given by The objective is to move a mobile sensor from a starting location at coordinates to a goal location at . The cost function is given by where is the state tracking error. This cost function is subject to constraints on the actuation of the mobile sensor which ensure that the maximum sensor velocity is significantly smaller than the largest background flow field velocity, . This constraint is imposed to model the limited actuation available in real world scenarios. Here, and , where is the identity matrix. The discretized time step across all double gyre simulations is . Each simulation was run for steps. In the following sections, we will vary the relative cost of actuation versus state error, given by the ratio , and analyse how this impacts the mobile sensor trajectories. To solve the resulting optimization problems, we use the CasADi [77] and MPCTools [78] packages, which use an interior point filter-line search algorithm (IPOPT).

Double gyre results

In this section, we examine energy-efficient trajectories for an active mobile sensor generated using MPC across a range of hyperpameters, including the prediction horizon, penalty weights on the state error and control effort, and the double gyre oscillation frequency. The MPC results for this thorough parameter sweep are presented in §4(d) and summarized in figure 7. Our goal is to understand the sensitivity of the trajectory to parameters and to uncover performance tradeoffs, for example with the time horizon of optimization. We find a large sweet spot where effective, energy-efficient trajectories are generated. Furthermore, we establish connections between the efficient mobile sensor trajectories and the LCS of the underlying flow field.

Figure 7

Multiple simulations were carried out at each ratio spaced logarithmically, from 0 to 100, time horizon ranging from 1 to 10, and gyre frequency ranging from to . The data presented here are in the form of scatter plots for each gyre frequency with each colour representing the Pareto optimal tradeoff curve between the total energy spent along each trajectory and the sum of deviations from target along the trajectory. The trajectories shown in the bottom row correspond to the highlighted purple circles (1,2,3) in the Pareto optimal corresponding to . (Online version in colour.)

Trajectories with different relative actuation cost,

Figure 2 shows the effect of varying the ratio of control effort penalty to the state error penalty on the trajectories, for a fixed time horizon of ; similar plots for a range of time horizons from to are shown in electronic supplementary material, figures S1–S9. The ratio quantifies the relative cost of actuation, and varying this parameter is important to understand performance tradeoffs when the mobile sensor has a limited actuation budget. As is increased, corresponding to actuation being more expensive, the agent actuates less, and the state tracking error increases. This increase in state tracking error tends to correspond to larger steady-state limit cycles about the goal state. The weighted actuation cost increases with , as we fix and increase ; however, the unweighted actuation decreases with . Importantly, the trend of cost versus is not strictly monotonic, and there are discontinuous jumps corresponding to bifurcations in the orbit; the non-monotonic behaviour and bifurcations are more pronounced for other in appendix A. For small values such as and , the agent moves around the goal state in a tight orbit, and this orbit continuously expands as increases, as shown for . However, between and the trajectory undergoes a rapid qualitative change, where the radius of the orbit around the goal state jumps.

Figure 2

Dependency of the resulting sensor trajectories on the ratio between the penalties on energy expenditure () and state error (), for a fixed time horizon . The trajectories are colour-coded by instantaneous energy expenditure. There are three costs shown in the top right figure: , , and . Here, generally increases with . The trajectories undergo several qualitative changes (bifurcations) as is increased, forming different types of periodic orbits shown on the bottom of the figure. An example of this is when is changed from 25 to 26, we observe a major change in the shape of the final orbit around the goal, as opposed to the minor change from to 25, where the final eye-shaped orbit only gradually increases in size. The formation of these orbits is dependent on the background flow FTLE, as can be seen from the inset of case . We observe that as the LCS move, the stable and unstable LCS intersect at a point, whose location changes every instant, and the sensor moves in a manner in which it is right on top of this intersection point for most of the time when (top left box). (Online version in colour.)

It is interesting to note in figure 2 that the agent has an initial loop in the right basin, while the agent does not. This behaviour is counterintuitive, as the agent should expend control more freely, and thus more aggressively seek the goal state. As shown in figure 15 in appendix A, the more aggressive agent does move away from the starting state faster initially; however, it becomes trapped on the side of a repelling LCS farther away from the goal location and must make an entire orbit around the right gyre before approaching the goal state. The maximum agent velocity is smaller than the maximum gyre velocity, so even the most aggressive agents are unable to break out of the right gyre without precise timing. This type of bifurcation also occurs for fixed by varying the time horizon, as in figure 3. In this case, the behaviour is more consistent with intuition, as the longer time horizon trajectories avoid being trapped in the right gyre.

Figure 15

In figure 2, we observe that for , , the trajectory makes a loop in the right gyre before ultimately reaching the goal. This is not the case for which is more greedy in spending energy. This is somewhat counterintuitive, since we expect that spending more energy should drive the mobile sensor more rapidly to the goal. In this plot, we show the difference in how the control energy is spent comparing and . We observe that takes longer to reach the goal because good performance is very sensitive to timing in the gyre. We see that the aggressive control (grey trajectory) moves too far ahead of the non-aggressive control (purple trajectory) to make use of the gyre dynamics to move directly to the goal state. (Online version in colour.)

Figure 3

Different instances of a trajectory for different time horizons while keeping fixed. By varying the time horizon, , we see that the extra loop in figure 2 is due to a sensitivity of the planned path with respect to , where the agent becomes stuck in the right gyre for lower . We can see here with reference to figure 15 how as the time horizon increases, the agent moves less aggressively and improves its timing with the gyre oscillations to reach the goal sooner without taking an extra loop. Another interesting aspect of this plot can be seen at and where the multiple agents line up to mirror the movement of the red LCS, further confirming the strong correlation finite-horizon optimal trajectories and the FTLE ridges. (Online version in colour.)

Previous work [24] suggests that low-energy trajectories tend to coincide with the LCS of the background flow. In our example, even for and , the mobile agent can be seen aligning with and exploiting the coherent structures. For example, in the top left of figure 2, the sensor moves along on the intersection of the attracting and repelling LCS as it orbits the goal state. In the next section, we will see that the agent also precisely times its actuation before and after crossing a repelling LCS to take advantage of the background drift.

Instantaneous energy versus finite-time Lyapunov exponent ridge

Given the existing connection of low-energy trajectories and FTLE ridges, we are interested in how the energy is used along a trajectory. Figure 4 shows how the agent ‘schedules’ an increase in actuation to cross a repelling (blue) FTLE ridge. After crossing, the agent decreases its actuation, as it is naturally repelled from the blue ridge and attracted by the red ridge into the left basin. Similar timing and utilization of the FTLE ridges is observed for a wide range of time horizons and control aggressiveness. This occurs in the double gyre only across the repelling hyperbolic LCS near the middle of the domain, as opposed to the shear-driven LCS near the goal location. This can be observed across the all gyre trajectories, where the yellow (high energy) regions in the trajectory are in the middle of the domain.

Figure 4

Influence of the FTLE ridges on the energy expenditure both along prediction horizon (right top) and instantaneously (right bottom) under parameters , and . An agent’s motion along with the predicted forecast trajectory (dashed line) are shown on the left, together with the repelling (blue) LCS and the attracting (red) LCS. There is a correlation between the spike in both the instantaneous energy spent (right top) and the cost along the forecast trajectory (right bottom) with movement across an FTLE ridge. Here, unlike in figure 2, the summations and are only along the forward dashed line in the plots under ‘snapshots’ and not along the entire trajectory as in figure 2. (Online version in colour.)

Periodic orbits

We observe that controlled trajectories often form periodic orbits around the goal state, as seen in figures 2, 5 and 6. Because the background flow field is periodic, the agent would require constant actuation to stay fixed at the goal state. Instead, the agent trajectory tends to form a periodic orbit around the goal, balancing state tracking error and control expenditure. Typically, this orbit is larger for agents with a tighter energy budget (i.e. for larger ). Many past studies have focused on trajectory planning where the final state is fixed at the goal. However, given the constantly evolving background flow field and its dominant effect on mobile sensor dynamics, it is important in practice, to consider the cases where the final state cannot be fixed. Figure 6 also indicates that the shape of the final periodic orbit depends on the frequency of the double gyre oscillation, with the frequency of the agent orbit synchronizing with the gyre frequency. In the other example flows observed below, similar periodic orbits are observed, where the agent loiters around the goal state. It will be interesting to investigate these orbits in more detail, including the classes of flows they exist in, and the conditions under which they bifurcate.

Figure 5

The mobile sensor settle on periodic orbits around the goal state (a) and the magnitude of the Fourier transform of the instantaneous energy spent by the mobile sensor (b). We observe that the time series of the energy spent is periodic with frequencies at integer multiples of the double gyre oscillation frequency, which correspond to the peaks in (b). (Online version in colour.)

Figure 6

Variations of agent trajectory under different double gyre oscillation frequencies. The frequency of the periodic orbits depends on the double gyre oscillation frequency. (Online version in colour.)

Model predictive control parameter sweep

We now present an exhaustive sweep through two of the most critical parameters for MPC, the prediction horizon and the cost function penalty ratio , for different gyre oscillation frequency . The first two parameters are related to the power and prediction capability of the mobile sensor, and the third parameter characterizes the unsteadiness of the background flow. We perform a full parameter sweep for the time horizon () and the cost penalty ratio , for the double gyre frequency . For each parameter value, we compute the state tracking error and the (unweighted) actuation energy expenditure, integrated along the entire trajectory.

Figure 7 shows the results from the MPC parameter sweep. For all time horizons and gyre frequencies, we observe that the trajectories sweep out a Pareto front in control expenditure versus state tracking error as is varied logarithmically from 0 to 100. The bottom row of figure 7 shows three representative trajectories along the Pareto front. As is increased, there is often a sharp drop in control cost with a relatively small increase in state tracking error, suggesting that there are energy-efficient trajectories that achieve relatively good tracking performance. However, we observe a break point in this monotonic trend, beyond which increasing results in rapid deterioration of the state error with relatively little decrease in control cost. This break point corresponds to the scenario where the motion of the sensor is dominated by the background flow, and the chaotic nature of the flow field dominates the state and energy errors. This phenomenon is more evident for smaller time horizons.

It is observed that longer prediction horizons produce trajectories that are more energy efficient with smaller state errors. This is expected, as longer time horizons include more information about the flow field in the optimization. This trend is weaker for small-to-moderate and is more pronounced for larger . The shape of the Pareto curve also changes with the double gyre frequency. This shape change is particularly evident for moderate frequencies, suggesting a ‘resonance’ in the interaction of trajectories with the background flow. Resonance with changing gyre frequency has been explored in the context of inertial particles in the double gyre flow [34].

Sensor velocity

To gain further insight into the dependency of sensor actuation velocity on the background flow velocity, we compare their distributions along the resulting trajectory at different values. Figure 8 shows histograms of the magnitude and orientation of the sensor actuation velocity versus the background flow velocity, for a range of values. It can be observed that more aggressive agents with smaller have larger actuation velocity magnitudes and tend to move perpendicular to the background flow. Agents with larger , corresponding to more conservative actuation policies, tend to have smaller actuation velocity and align their actuation in the direction of the flow field to take advantage of the background flow. Except in the most aggressive case, the mobile sensor rarely uses the maximum control velocity. Additional plots with the - and -components of the agent velocity are presented in figure 16 in appendix A.

Figure 8

As the sensor moves in the double gyre flow field, it is constantly taking control actions , where are the x- and y-components of its actuation respectively. At each instant, the sensor is also moving over a background double gyre flow velocity vector whose components, and , are given by (3.4). The top row of histograms are of the magnitude of control actions taken (in red), against the magnitude of the background current velocity (in grey), where are the sensor coordinates at time . The second row shows the heading angle of the sensor (the orientation of ), plotted in red, against the orientation of the background flow field velocity vector, plotted in grey. The corresponding trajectories are shown on the bottom. (Online version in colour.)

Figure 16

As the sensor moves in the double gyre flow field, it is constantly taking control actions , where , are the x- and y-components of its actuation, respectively. At each instant, the sensor is also moving over a background double gyre current vector whose components , are given by (3.4). The top row is a histogram of the x-component of control actions taken (in red), against the x-component of the background current velocity (in black), where , are the sensor coordinates at time . The second row is a similar plot for the y-component. We observe that the gyre takes values beyond the actuation capacity of the sensor, which highlights the under-actuated nature of the problem. Also, at low ratios, the distribution of control actions follows a distribution with two peaks at , which corresponds to a situation similar to bang-bang control. As we increase the ratio, the distribution of control actions moves to a single peak centred around zero corresponding to the use of very little control effort when compared with the background velocity. (Online version in colour.)

Different start and end locations

Figure 9 shows the MPC optimized trajectories for six different starting locations along the right and lower boundaries of the domain. Although the paths have different initial transients, the trajectories evolve onto the same periodic loitering orbits around the goal state, indicating that they are ultimately leveraging similar flow structures. Similarly, figure 10 shows the MPC optimized trajectories for several different goal states. Some goal states are much more difficult to reach than others, because of the strong unsteady background flow field. With more aggressive control () the trajectories generally form tighter loitering orbits. However, for the case when the goal state is in the upper middle of the domain, it is clear that none of the MPC trajectories are able to find a suitable loitering pattern. This diversity of orbits highlights the importance of choosing a suitable goal location, which would likely involve a higher level of planning.

Figure 9

MPC trajectories with varying ratio, starting from different initial conditions with the same goal location. We find that eventually, the trajectories converge to similar periodic orbits thereby using the LCS in similar ways to orbit around the goal location despite having different transients. (Online version in colour.)

Figure 10

MPC trajectories to different goal locations with varying ratio for the same initial condition. We find that in the case of placing the goal near the middle to bottom region (top left and top middle plots) of double gyre, the sensor is able to form stable orbits near the attracting LCS, where the sensor moves with the base of the attracing LCS. However, placing the goal near the repelling LCS causes more difficulty for MPC in forming small stable orbits (top right). The bottom plots show that it is possible to form small periodic orbits in the corners of the double gyre flow field when the goal is placed close to them. (Online version in colour.)

Advanced test cases

From the analyses in the double gyre flow field, we observed the strong dependency of energy-efficient MPC trajectories on the LCS of unsteady periodic flow fields. In this section, we study the use of MPC on more challenging test cases: an analytical three-dimensional incompressible flow field and a real-world flow field reconstructed from ocean model datasets. The goal is to demonstrate how and when these results generalize to better understand the dependency of energy-efficient MPC trajectories on FTLE-based LCS.

Arnold–Beltrami–Childress

We now demonstrate the use of MPC on the ABC flow field [79,80], which is an incompressible model for a flow evolving in a three-dimensional periodic domain. It has been studied exhaustively in the past as a stepping stone to understanding turbulent flow fields. The three-dimensional space contains six interwoven vortices. An important feature of this flow field is that even the steady version of the flow field can give rise to chaotic trajectories. However, for our purposes, we investigate the unsteady case. Similar to (3.5) for the double gyre flow field, the equation for a mobile sensor evolving in the ABC flow is given by with parameters , where is a time-varying component that makes the flow field unsteady. We investigate trajectories in the regime where , which has been exhaustively studied numerically and analytically (with specifically in our simulations), and . All ABC simulations presented were run for time steps with a step size of . The cost function used was the same as equation (3.6). In figure 11, we show the trajectories planned by the MPC for a time horizon . Across these simulations, we observed that, similar to the double gyre, it is possible to form loitering orbits close to the goal point. These loitering orbits become larger for larger . Results across different time horizons and are further summarized in figure 12. We found that relatively short time horizon (compared to the period of oscillation ) trajectories are able to reach the goal state, and similar to the double gyre, longer time horizons reduce error with lesser energy consumption. The inflection points and breaking off of points from the curve (for example, the black dots near ) correspond to drastic changes in trajectory shape (bifurcations).

Figure 11

MPC trajectories formed by the sensor in an ABC flow field as a function of the ratio. We observe that for a time horizon of , we are able to find several cases of periodic orbits loitering close to the goal location. In these cases, the initial position is and the goal is . (Online version in colour.)

Figure 12

This scatter plot shows the performance of MPC in the ABC flow field similar to figure 7. We observe that MPC trajectories with short time horizons are able to reach the goal state and increasing the time horizon has a benefit of decreasing the actuation energy usage. We also observe that the inflection points and breaking off points correspond to bifurcations in the trajectories. (Online version in colour.)

Figure 14

This figure shows the change in MPC trajectories in the Gulf of Mexico as the ratio is varied. We see an aperiodic loitering state near the goal location. Given that the time horizon in this case h, which is relatively short compared to the total trajectory time of days, we can see that even with short time horizons the MPC trajectories are capable of making it to the goal location. (Online version in colour.)

Gulf of Mexico

For the final example, we consider the Gulf of Mexico surface velocity estimates from the HYbrid Coordinate Ocean Model (HYCOM). This data-assimilative model synthesizes remotely sensed and in situ measurements on a hybrid coordinate system. We used daily -resolution data from the HYCOM 1992–1995 experiment 19.0 (top left in figure 13) to generate a vector field. We then used linear interpolation on this vector field in space and time to generate a function that could be used for model predictive control. The parameters were chosen to be , , , and . The step size for the MPC. The full trajectories in figure 14 were computed for 1000 time steps (100 days).

Figure 13

This figure shows the use of MPC to plan trajectories in the Gulf of Mexico dataset. We have chosen the region highlighted in dashed-dotted lines in the top left plot. The bottom six plots show the trajectory generated for (full trajectory can be seen in the rightmost plot in figure 14) and in colour shading from yellow to purple to highlight how energy is spent along the path. The units on the - and -axis are longitude and latitude, respectively. The red cross shows start location, , the green cross shows the goal location, . The black dot shows the instantaneous sensor location. Viewing the six panels in sequence we observe that the sensor moves across the blue repelling LCS from to , where is the time in days. We observe a spike in the instantaneous energy spent and a slow drop as we move away from the repelling LCS (as seen in the top right plot of ). We then observe that the sensor synchronizes with the attracting LCS to move towards the goal location. (Online version in colour.)

We observe similar spiking behaviour in the energy spent when moving across a repelling LCS, which is seen in panels to of figure 13. Here, is the time in days (specifically, corresponds to the flow field on Day- of the HYCOM dataset). We observe that the sensor synchronizes with the attracting LCS to move towards the goal location. In periodic flow fields such as the double gyre, we observe the formation of periodic orbits. However, in case of the Gulf of Mexico (aperiodic), we observed the formation of an aperiodic loitering trajectory in proximity to the goal location. The time horizon here is relatively short compared to the total trajectory time. We can see that even with short time horizons the MPC trajectories are capable of making it to the goal location. Although the trajectory does evolve towards a loitering orbit near the goal state, it can be seen in figure 14 that for more aggressive control ( and ), the trajectory is able to get much closer to the goal state, with a tighter orbit.

Discussion and conclusion

In this work, we have investigated the behaviour of trajectories optimized over a finite horizon for a controlled mobile sensor in unsteady flow fields, as both the control and flow field parameters were varied. A thorough study was conducted on the double gyre flow field, and selected results were further verified on more advanced flow fields such as the ABC and the Gulf of Mexico. In particular, finite-time MPC was used to generate energy-efficient trajectories for a range of parameters, particularly the prediction horizon and the relative penalty between the state error and control effort. The double gyre oscillation frequency was varied to study its influence on the resulting trajectories. We have constrained the maximum actuation velocity to be less than the largest background flow velocity such that some degree of intelligent planning is required to efficiently traverse the flow.

Through a quantitatively exhaustive study, we have uncovered several interesting trends and established connections between the finite-horizon optimized mobile sensor trajectories and the coherent structures of the underlying flow field. By varying the relative cost of actuation and deviations in the state (i.e. ), the control cost and state error sweep out a Pareto front, and there is often a sweet spot where relatively good state tracking performance can be achieved with low actuation costs. These energy-efficient trajectories tend to align with the LCS to take advantage of the unsteady background flow. We also found that locations where the sensor spends most energy correlated with the presence of repelling hyperbolic LCS. These findings could be better understood by recalling that repelling LCS are defined as material barriers in the flow field. This explains the expenditure of more energy close to the repelling LCS, as energy needs to be spent to overcome barrier for movement. Furthermore, as the energy spent for movement goes to zero, the sensor behaves similar to a passive tracer. It is known that the short-term dynamics of passive tracers are governed by movement with the FTLE. This suggests that these connections hold for energy optimality. Importantly, we find that it is often possible to generate effective, energy-efficient mobile sensor trajectories with a relatively short prediction horizon, which is promising for the future design of trajectories with limited or partial knowledge of the background flow field.

We observe a rough trend of lower state error when control is less expensive, which agrees with the intuition that the agent is able to more directly pursue the goal state by actuating more aggressively. However, this trend is not monotonic, as there are several cases where slightly decreasing the control cost results in worse state tracking performance. These non-monotonic changes in the cost versus correspond to bifurcations in the agent trajectories, which either correspond to longer trajectories, or to discontinuous jumps in the shape and size of the periodic orbit around the goal state. These bifurcations are more common for smaller prediction horizons, which is also consistent with the intuition that smaller prediction horizons may lead the agent to get trapped by unfavourable flow structures. Similarly, for a fixed relative control cost, there are bifurcations in the optimized trajectory with variations in the time horizon. These bifurcations are relevant in the context of generating mobile sensor trajectories using MPC, as small changes in the weights can lead the drastically different trajectories. Upon closer inspection, these bifurcations correspond to the agent trajectory passing through a Lagrangian coherent structure, after which the two trajectory behaviours diverge.

It is also important to note that the energy-efficient trajectories typically result in periodic orbits around the target position, since the unsteady double gyre is periodically oscillating. Previous studies in trajectory generation have mainly focused on solving boundary value optimizations for trajectories keeping the start and end points fixed. Our results show that these assumptions can be relaxed, and moreover, it is possible to reach periodic steady states with little actuation even when the uncontrolled drifter dynamics are chaotic. These periodic orbits correspond to, often desirable, station keeping or hovering behaviour.

This work has several implications for the control of individual mobile robots and swarms of robots in geophysical flows. The ability to generate energy-efficient trajectories that take advantage of the background flow with a short prediction horizon is promising for practical applications. The ability to maintain close periodic orbits around the goal state may also enable efficient long-time monitoring. For example, fix-wing unmanned aerial vehicles must often loiter over an area for sensing and monitoring. Variations in the shape of periodic orbits and the Pareto optimal curves over different ratios with the gyre oscillation frequency have implications for ocean applications, which exhibit a wide range of spatio-temporal scales with varying oscillating frequencies. We also observed an increase in the expenditure of the sensor’s actuation energy as it approached background LCS. This result is beneficial in the context of identifying background coherent structures by observing the energy expenditure patterns of controlled agents. This is an important problem with ongoing work [81-83]. These results are also potentially useful in the design of scalable navigation algorithms for mobile sensor swarms where the objective is to maintain cohesion or connectivity between agents.

This work motivates a number of interesting future directions. Our results indicate that it is possible to design nearly optimal, energy-efficient trajectories, even with short prediction horizons for the MPC; however, it was assumed that the background flow was known perfectly for this short horizon. It will be important to further explore the robustness of these trajectory optimizations to more realistic scenarios with partial, noisy and uncertain information about the background flow. This analysis may benefit from recent works that have investigated the sensitivity of FTLE calculations to uncertain flow field data [84,85] as well as how FTLE can be used to propagate uncertainties through chaotic flow maps [86]. Because the optimization result depends strongly on how the MPC trajectories interact with LCS of the background flow field, it may also be possible to incorporate knowledge about the LCS more directly to the optimization. Even with uncertain or partially observed flow field information, often the LCS are quite persistent, and it may be possible to develop time-varying maps of the coherent structures in different geographical regions, for example off the Horn of Africa or in the Gulf of Mexico. In addition, it will be interesting to explore the use of other coherent structure and modal decomposition identification techniques [87,88]. Further study is also required to characterize the dynamics and coherent structures of the controlled vector field of the agent given a specific control policy. In addition, all the results in this paper were developed through the study of the double gyre flow field. It will be interesting to perform similar investigations for a variety of flow fields. For example, it will be important to explore how these results change when the flow exhibits a wider range of multiscale behaviour in space and time. Particularly, multiscale turbulence will affect the prediction horizon, as uncertainties will be magnified making it challenging to forecast flow structures. These multiscale structures will also impact energy efficiency, both through the forecast uncertainty, but also through making optimal paths more circuitous. Finally, extending the analysis to additional three-dimensional turbulent flows will also be critical.

Click here for additional data file.