Performance Analysis of Ice-Relative Upward-Looking Doppler Navigation of Underwater Vehicles Beneath Moving Sea Ice

Abstract

This paper addresses the problem of ice-relative underwater robotic vehicle navigation relative to moving or stationary contiguous sea ice. A review of previously-reported under-ice navigation methods is given, as well as motivation for the use of under-ice robotic vehicles with precision navigation capabilities. We then describe our proposed approach, which employs two or more satellite navigation beacons atop the sea ice along with other precision vehicle and ship mounted navigation sensors to estimate vehicle, ice, and ship states by means of an Extended Kalman Filter. A performances sensitivity analysis for a simulated 7.7 km under ice survey is reported. The number and the location of ice deployed satellite beacons, rotational and translational ice velocity, and separation of ship-based acoustic range sensors are varied, and their effects on estimate error and uncertainty are examined. Results suggest that increasing the number and/or separation of ice-deployed satellite beacons reduces estimate uncertainty, whereas increasing separation of ship-based acoustic range sensors has little impact on estimate uncertainty. Decreasing ice velocity is also correlated with reduced estimate uncertainty. Our analysis suggests that the proposed method is feasible and can offer scientifically useful navigation accuracy over a range of operating conditions.

1. Introduction

This paper addresses the problem of precision ice-relative navigation of Uninhabited Underwater Vehicles (UUVs) in the upper water-column under moving sea ice floe in the polar ice-pack—a Global Positioning System (GPS)-denied undersea environment in which conventional downward-looking bottom-lock Doppler sonar navigation is generally not possible due to excessive water depth below the vehicle.

This study seeks to evaluate quantitatively, in simulation, the performance of an underwater navigation system comprised of the following navigation sensors:

• Two or more Global Navigation Satellite System (GNSS)1 transceivers deployed on the moving ice floe to instrument ice position and orientation.
• A surface ship is equipped with a precision GNSS receiver, a true-North seeking gyrocompass, and two acoustic modems providing acoustic ranging and telemetry to the underwater vehicle(s). The advantage of two acoustic beacons is that it enables you to compute a complete position fix when the underwater vehicle is not moving relative to the ship.
• The UUVs are equipped with upward-looking Doppler sonars, precision pressure depth sensors, and true-North seeking gyro compasses.

This paper addresses the ice-relative navigation problem in a state estimation framework, where states are comprised of 28 linear and angular positions and velocities of the vehicle, ship, and ice. The paper reports a novel performance analysis of an approach to ice-relative navigation originally reported in [1]. The principal goal of this study is the sensitivity analysis to evaluate quantitatively the effects on navigation precision of the following operational parameters:

• The effect of variation in separation of ice-top GNSS beacons.
• The effect of variation in number of ice-top GNSS beacons—i.e., more than two GNSS beacons.
• The effects of variation in ice velocity (rotational and translational).
• The effects of variation in separation of ship-deployed acoustic ranging modems.

The effects of these variations on the estimates’ error and covariance are examined, with special attention given to navigation accuracy at the end of the simulated mission, when ship-to-vehicle distance is highest.

The remainder of the paper is organized as follows: Section 2 provides a brief outline of previously reported under-ice navigation methods and motivates the need for precision under-ice navigation, Section 3 details our state definition, process and observation models, Section 4 describes the simulation environment, and summarizes parameters under examination. Section 5 presents the results of the simulation studies. Finally, Section 6 concludes and summarizes.

2. Background and Literature Review: Under-Ice Navigation

For an extensive review of the scientific motivation and challenges, development, and use of underwater robotic vehicles designed for use in ice-covered waters, with special attention paid to the navigation systems employed for under-ice deployments, the reader is directed to [2]. This paper reviews scientific needs for routine access under fixed and moving ice by underwater robotic vehicles are reviewed in the contexts of geology and geophysics, biology, sea ice and climate, ice shelves, and seafloor mapping. The challenges of under-ice vehicle design and navigation are summarized. The paper then reviews all known under-ice robotic vehicles and their associated navigation systems, categorizing them by vehicle type (tethered, untethered, hybrid, and glider) and by the type of ice they were designed for (fixed glacial or sea ice and moving sea ice).

Few methods presently exist for precision navigation of UUVs for the benthic survey and sampling operations under fixed and moving sea ice. Land-fast (fixed) ice, provides an inertial reference surface against which UUVs can utilize conventional upward-looking Doppler Velocity Log (DVL)-based navigation methods. For example, the Autosub2 vehicle used upward-looking and downward-looking DVLs to provide a Dead Reckoning (DR) navigation solution while traversing >20 km beneath the Fimbul Ice Shelf in Antarctica [3]. Acoustic range-base methods including Long Baseline (LBL) [4] and Ultrashort Base Line (USBL) have been deployed through the ice for vehicles operating beneath land-fast ice [5,6].

Navigation beneath free-floating ice is more challenging because ice-floes are free to rotate and translate. Doppler sonars can provide vehicle velocity measurements when in downward-looking bottom-lock or upward-looking ice-lock range (up to ∼200 m depending on sonar frequency) with Gaussian measurement noise with a frequency-dependent single-ping standard deviation on the order of few mm/s, and update rates up to 10 Hz [7], thus enabling the development of a wide variety of Doppler-based navigation techniques—e.g., [8,9,10,11]. Navigation errors arising in Doppler-based position estimates (without external corrections) grow proportional to the square-root of time [8,12]. In [13], the authors report a system in which acoustic transducers were deployed through the ice along the vehicle’s intended trajectory and used to manually update and bound the vehicle’s dead-reckoned position and error. Cress et al. report a system in which low frequency (1376 Hz) homing beacons were suspended from the sea ice during missions beneath sea ice in the Canadian Arctic [14]. They report the homing beacons were instrumental for mission success, because accumulated navigational error and drifting of the sea ice camp could have resulted in a vehicle recovery location >30 km away from the anticipated recovery location.

Other acoustic methods of under-ice navigation have also been reported in the literature. Kimball and Rock report proof-of-concept results for a two-step terrain-relative navigation scheme for iceberg-relative navigation. In the first step, a three dimensional (3-D) map of the iceberg is generated with a multibeam sonar following the circumnavigation of the iceberg. In subsequent missions, the vehicle compares sonar range measurements to the previously generated map to localize itself relative to the iceberg [15]. Kimball et al. report a hybrid under-ice navigation system employing DVL dead-reckoning over long ranges from a “dock”, inverted USBL within 700 m of the dock, and optical stereo imaging within 75 m of the dock [16]. Sayre-McCord et al. report an approach in which, in post processing, utilizes upward-looking optical imagery in a Simultaneous Localization and Mapping (SLAM) framework to correct vehicle navigation for ice-floe mapping [17].

Webster et al. report a long-range low-frequency cooperative acoustic navigation system for multiple Seaglider AUVs beneath Arctic sea ice in which eleven ice-tethered GPS/acoustic navigation beacons broadcast telemetry and control messages to the vehicles, which also provided One Way Travel Time (OWTT) time-of-flight range measurements [18]. The Seaglider vehicles are not equipped with DVLs, and routinely dive to 1000 m depth—well beyond the range of DVL-lock to surface ice. When submerged, this approach employed signals of acoustic model ranges, signals from an attitude sensor, signals from a pressure depth sensor, and signals from the Seaglider’s buoyancy engine. This approach employs a kinematic process model in an Extended Kalman Filter (EKF) for bounded-error navigation.

Webster et al. report deployments in 2017 of two Seagliders equipped with WHOI Micromodems and acoustic Doppler current profilers (ADCPs) in the Canada Basin [19]. Low-frequency (250 Hz) OWTT acoustic transmissions from fixed tomography sources and on-board ADCP water-velocity profiles were used to estimate (in post-processing) the local current profile on a per-dive basis and the glider’s relative velocity through the water. Graupe et al. report studies of the 2017 Canada Basin Seaglider deployment data, with very large OWTT acoustic ranges of up to 480 km, in which post-processed OWTT vehicle navigation estimates with the acoustic arrival matching method, which infers ranges from comparisons with acoustic propagation models, were shown to reduce navigation error by a factor of 4 or 5 in comparison to real-time OWTT vehicle navigation estimates [20].

Navigational methodology and results from a 2014 expedition to ${83}^{°}$ N $6^{°}$ W in the Arctic with the Nereid Under-Ice (NUI) vehicle are reported in [21]. During the mission, the angular position of the ice floe was assumed constant, and ice-relative Dead Reckoning (DR) with an upward-looking DVL was used to navigate the vehicle. Images from ship-based ice radar were used in post-processing to estimate ice floe position and orientation relative to the ship, and the vehicle’s trajectory was then renavigated in the ice frame without the assumption of constant ice angular position [21].

The approach described herein is inspired by and builds upon previously reported approaches to ice-relative navigation including that employed upward-looking DVLs for ice-relative dead-reckoning [3,14,21,22], and the approach reported in [18] which modeled UUV process dynamics and sensor observations in an EKF framework, but for the case of mesoscale operations with deep Seaglider Autonomous Underwater Vehicles (AUVs) that are not equipped with DVLs. What is new about the present proposed approach is the introduction of process and observation models for ice motion, in addition to process and observation models for the ship and UUV, and the use of ice-relative upward-looking DVL velocity measurements between the UUV and moving sea ice, all in an EKF framework.

Acoustic propagation in polar waters differs from that at temperate latitudes. In temperate waters a sound velocity minimum generally occurs at depths below 500 m thus creating a convergent deep sound channel. In polar waters such as the central Arctic Ocean, sound velocity increases monotonically with depth, thus creating an upwardly refracting waveguide near the underside of the ice [23]. It is well known that extreme sound velocity profile variations with depth can ray-bending of acoustic paths that reduce the accuracy of acoustic ranging [24]. We also note that upper water column sound-speed profiles under moving sea ice that we have experienced with NUI during Arctic expeditions have been relatively benign. For example, the CTD profile for NUI dive 16 at 87${}^{°}$ N 61${}^{°}$ E on 27 September 2016 show a sound velocity profile in which the sound velocity increases monotonically with depth by 23 m/s in the top 300 m of the water column (from 1437 m/s at the surface to 1460 m/s at 300 m depth) and remains nearly constant from 300 m to 600 m.

Scientific Motivation for Under Ice Navigation

The world’s oceans cover 71% of the Earth’s surface, 12% of which is largely inaccessible to scientific research due to being covered by ice all or part of the year. In the Northern Hemisphere, sea ice coverage varies seasonally from 102% to 192% the size of the United States, and in the Southern Hemisphere sea ice coverage varies seasonally from 39% to 260% the size of Australia [25].

The importance and complexities of ice-covered waters in global climatological, geophysical, and biological processes are still being discovered. During a 2011 expedition to the Arctic, scientists observed a phytoplankton bloom in which algal concentrations indicated that previous estimates of net primary productivity by under-ice Phytoplankton may have been 10-fold too low [26]. Many ultra-slow Mid-Oceanic Ridgess (MORs) occur in geographic regions where weather windows are extremely narrow or there is ice cover (e.g., the Southern Ocean between 40${}^{°}$ S and 50${}^{°}$ S, the Australian-Antarctic Discordance, and in the Arctic—the Gakkel Ridge [27,28]). With the recent identification and first-order mapping studies of ultra-slow spreading ridges in the Arctic [28] and the Indian Ocean [29,30], ocean scientists are poised to make breakthroughs in our understanding of this important end-member of the seafloor spreading environment.

Geo-referenced UUV navigation is a highly desirable component of modern expeditionary ocean science. Geo-referenced UUV navigation enables multi-modal data sets to be co-registered across multiple dives with a single asset and across different assets. It enables photographic, bathymetric, and other geophysical maps to be related to generally lower-resolution surface-derived contextual maps. For a survey of well-established navigation methods commonly employed for underwater vehicles in temperate latitudes, which is beyond the scope of the present paper, the reader is directed to [11]. Vehicle operations in ice-covered waters present additional challenges, which resulted in the development of a relatively new class of UUVs. For a survey of under-ice UUVs, and the navigation systems used in each, the reader is referred to [31].

A motivating use-case for the navigation approach described herein is [32], in which the authors describe a survey of light transmission through sea ice with the Nereid Under-Ice (NUI) vehicle at ${83}^{°}$ N $6^{°}$ W in which the co-registration of above-ice survey measurements and under-ice survey measurements were established manually by drilling holes in the ice along a 90 m transect, and dropping markers (poles) through holes in the ice that NUI could visualize with its cameras and sonars. The goal of the methods described herein is to provide a navigation approach that could enable much larger scale ice-relative under-ice surveys without resorting to manual markers.

3. Formulation: Ice-Relative Navigation

We consider the problem of estimating the six-degree of freedom (6-DOF) position and orientation and 6-DOF linear and angular velocities of an underwater vehicle relative to moving Polar sea ice floes. Polar ice floes exhibit time-varying translational and rotational motion in response to forces arising from the wind, ocean currents, and pressure of the ice-pack itself [25].

For this preliminary study, we chose to implement a centralized EKF state estimator with kinematic process models for the vehicle, ship, and ice, inspired by the centralized EKFs reported in [18,33]. The centralized EKF formulation is appropriate for use in post-processing navigation data, where sensor information from the UUV, ship, and ice are all available globally. If such a system were to be implemented for real-time operation, where sensor information is only known locally, we would adopt a decentralized estimation approach similar to that reported in [34]. Alternative centralized and decentralized approaches such as error-state Kalman filters are also possible, but are beyond the scope of the present study. The present preliminary study addresses the performance of a loosely-coupled GNSS/INS system. The additional study of the performance benefits of a tightly-coupled GNSS/INS for the ship, in comparison to more conventional approaches, is beyond the scope of the present study.

Sensor measurements from the sensors outlined in Table 1 provide asynchronous observations of the system, while the process model is propagated with constant time steps in between asynchronous observations.

GNSS receivers directly estimate geodetic position, and only indirectly estimate geodetic velocity from the position estimates. The manufacturer’s spec sheets of the GNSS sensors that we simulated (Garmin 18x LVC [35] and Trimble SPS852 [36]) provide precise specifications for the directly-measured position noise statistics. The manufacturers do not provide noise statistics for the indirectly-estimated velocity estimates. In consequence, we chose to simulate the GPS position outputs for which we have well-documented noise statistics.

We assume the following: First, we assume that the ice floe above the vehicle is a single contiguous and rigid body that can translate horizontally and rotate in yaw (only) on the ocean surface. This assumption is reasonable for vehicle operations beneath large ice-floes, shelves, or floating glaciers, but breaks down in areas occupied by numerous smaller ice floes. Note that we neglect the roll and pitch of the ship and the ice floe because ocean swells are not present in dense pack ice. It is our experience in the Arctic that floes of several km in diameter and larger are common at higher latitudes in heavy ice cover, and that low frequency (3.5 kHz) acoustic modem communication is often possible at ranges of up to several km. Second, we assume the ability to instrument the ice with at least two GNSS nodes, labeled GNSS₁ and GNSS₂ respectively. Third, we assume an idealized Doppler Velocity Log (DVL), which measures velocity relative to a point along the instrument’s z-axis; In reality DVLs typically provide velocity information from three or more beams typically angled 30${}^{°}$ off the instrument’s z-axis with equal radial spacing.

We make several assumptions about the under-ice acoustical environment to make the simulation more tenable, namely instantaneous, straight, and fully available acoustic propagation. The complexity of the under-ice acoustic environment cannot be understated, and there is evidence to suggest that climate-induced changes to Arctic waters could significantly inhibit both long and short-range acoustic communications [37]. The authors acknowledge that a more detailed treatment of the acoustic environment could make use of the BELLHOP algorithm for estimating the availability, and true distance traveled by an acoustic ray under an assumed sound velocity profile [38], but as the simulation rarely exceeds a 1 km standoff distance from the vehicle to the ship, and the acoustic paths are largely horizontal, such effects should be minimal. This is consistent with our experience that ice-relative navigation with an upward-looking DVL necessitates that the vehicle remain within Doppler-lock range of the ice, typically with a survey depth less than 50 m, with a maximum DVL-lock depth of about 150 m in our experience, for a 300 kHz RDI Broadband DVL. This means that the slant ranges are largely horizontal for under-ice surveys in which the vehicle horizontal range might extend to 1 km.

3.1. Notation

We employ the following notation:

• ${}^a{\mathit{p}}_b\in {\mathbb{R}}^3$ is the vector to the origin of frame b, represented in frame a (i.e., the vector from the origin of frame a to the origin of frame b, represented in frame a). This may be written component-wise as ${}^a{\mathit{p}}_b={\left[{}^a{x}_b,{}^a{y}_b,{}^a{z}_b\right]}^{\top }$.
• ${}^a\mathit{p}\in {\mathbb{R}}^4$ is the homogeneous point p, represented in the a frame.
• ${}_b^a\mathit{R}\in SO\left(3\right)$ is the rotation matrix from frame b to frame a where $SO\left(3\right)=\{\mathit{R} : \mathit{R} \in {\mathbb{R}}^{3\times 3}, {\mathit{R}}^{\top }\mathit{R}=\mathit{I}, \det \left(\mathit{R}\right)=+1\}$.
• ${}^c{\mathbf{\nu }}_a^b\in {\mathbb{R}}^3$ is the relative linear velocity between frames a and b, represented in frame c.
• ${}^c{\mathbf{\omega }}_a^b\in {\mathbb{R}}^3$ is the relative angular velocity between frames a and b, represented in frame c.
• ${}_b^a\mathit{H} \in SE\left(3\right)$ is the homogeneous transformation from frame b to frame a, where $SE\left(3\right)=\{(p,\mathit{R}) : p \in {\mathbb{R}}^3, \mathit{R} \in SO\left(3\right)\}={\mathbb{R}}^3\times SO\left(3\right)$, e.g., ${}_b^a\mathit{H}=\left[\begin{array}{cc}{}_b^a\mathit{R}& {}^a{\mathit{p}}_b\\ 0_{1\times 3}& 1\end{array}\right]$.

3.2. Coordinate Frames

We chose to reference all estimator states relative to the Local World Coordinate Frame We considered using a non-inertial ice-relative reference frame for the UUV state, but abandoned it when it became clear that it was simpler to represent the ice, ship, and UUV in a single common inertial frame. We employ the following coordinate frames in the formulation of the EKF:

• Local World Coordinate Frame: This is a Euclidean inertial coordinate frame that is fixed in geodetic coordinates and is considered to be an inertial frame (we neglect Earth rotation). The location of this frame is arbitrary, but is chosen to be coincident with the Ice Coordinate Frame at $t=t_0$. This frame is denoted by w.
• Ice Coordinate Frame: This is a non-inertial Euclidean coordinate frame that is rigidly attached to the ice floe floating on the ocean surface above the UUV. The origin of this frame is coincident with GNSS${}_1$, with its x-axis extending through GNSS${}_2$, and z-axis extending down. This frame is denoted by i.
• Vehicle Coordinate Frame: this is the UUVs body-fixed body coordinate frame, and is denoted by v.
• Ship Coordinate Frame: this is the ship’s body-fixed coordinate frame, and is denoted by s.
• Earth-Centered Earth Fixed (ECEF) Coordinate Frame: this is the Euclidean coordinate frame centered at the Earth’s center, with z-axis extending through the North Pole, and x-axis extending through the intersection of the Prime Meridian and the Equator. This frame is denoted by $ECEF$.
• DVL-insonified Coordinate Frame: This frame is on the underside of the ice and is rigidly attached to the ice at the point insonified by the UUV’s DVL at the instant a DVL measurement is taken. It is aligned with the Vehicle Coordinate Frame, and denoted by p.
• Instrument Frames: We denote the DVL, Inertial Measurement Unit (IMU), pressure, and GNSS${}_n$ sensor frames by $D, I, P$ and $G_n$ respectively.

3.3. State Definition

The complete state vector for the system, $\mathit{x}\in {\mathbb{R}}^{28}$ is composed of the combined states of the vehicle, ${\mathit{x}}_v\in {\mathbb{R}}^{12}$, ship, ${\mathit{x}}_s\in {\mathbb{R}}^8$, and ice, ${\mathit{x}}_i\in {\mathbb{R}}^8$, and is given by

$$\mathit{x}={[{\mathit{x}}_v^{\top }, {\mathit{x}}_s^{\top }, {\mathit{x}}_i^{\top }]}^{\top }.$$

(1)

The UUV’s 12-DOF state vector (3 positions, 3 Euler angles, and 6 velocities), ${\mathit{x}}_v$, is given by

$${\mathit{x}}_v={[{}^w{\mathit{p}}_v^{\top },{\mathit{\phi }}^{\top },{\mathit{\nu }}^{\top },{\mathit{\omega }}^{\top }]}^{\top },$$

(2)

with constituent vectors defined as

$${}^w{\mathit{p}}_v=\left[\begin{array}{c}{}^w{x}_v\\ {}^w{y}_v\\ {}^w{z}_v\end{array}\right],\mathit{\phi }=\left[\begin{array}{c}\varphi \\ \theta \\ \psi \end{array}\right],\mathit{\nu }=\left[\begin{array}{c}u\\ v\\ w\end{array}\right],\mathit{\omega }=\left[\begin{array}{c}a\\ b\\ c\end{array}\right]$$

(3)

where ${}^w{\mathit{p}}_v\in {\mathbb{R}}^3$ and $\mathit{\phi }\in {\mathbb{R}}^3$ are the UUV’s position and attitude expressed as XYZ position and and Euler-angle roll, pitch, and yaw coordinates in the world frame, and $\mathit{\nu }\in {\mathbb{R}}^3$ and $\mathit{\omega }\in {\mathbb{R}}^3$ are the UUV’s linear and angular velocities, expressed in vehicle coordinates. While the Euler angle representation of attitude suffers from a singularity at $\theta =\frac{\pi }{2}$, which can, in turn, lead to loss of estimator consistency, we argue that the possibility of such an orientation can be reasonably ignored given the passive stability in roll and pitch of most survey-grade UUVs.

The ship’s 8-DOF state vector (3 positions, heading, and 4 velocities) is given by

$${\mathit{x}}_s={[{}^w{\mathit{p}}_s^{\top }, {\psi }_s, {}^w{\dot{\mathit{p}}}_s^{\top }, {\dot{\psi }}_s]}^{\top }$$

(4)

where ${}^w{\mathit{p}}_s\in {\mathbb{R}}^3$ are the ship’s Cartesian world coordinates and ${\psi }_s$ the ship’s yaw, and the dot operator signifies the time derivative.

The ice floe 8-DOF (three positions, heading, and four velocities) state vector is given by

$${\mathit{x}}_i={[{}^w{\mathit{p}}_i^{\top },{\psi }_i,{}^w{\dot{\mathit{p}}}_i^{\top },{\dot{\psi }}_i]}^{\top }$$

(5)

where ${}^w{\mathit{p}}_i\in {\mathbb{R}}^3$ and ${\psi }_i$ are the floe’s Cartesian coordinates and yaw as represented in the world frame, and the dot operator signifies the time derivative.

3.4. Process Models: UUV, Ship, Ice

We assume a constant-velocity kinematic process model (as in [33,34,39]) for the system, with the continuous-time non-linear equation of motion given by

$$\dot{\mathit{x}}\left(t\right)=\mathit{f}\left(\mathit{x}\left(t\right)\right)+\mathit{G}\mathit{w}\left(t\right)$$

(6)

where $\mathit{f}\left(\mathit{x}\left(t\right)\right)$ is the deterministic portion of the differential equation, and $\mathit{G}$ is a diagonal matrix which maps the process noise, $\mathit{w}\left(t\right)\sim N(0,\mathit{Q})$ to the rate derivatives $\dot{\mathbf{\nu }}, \dot{\mathbf{\omega }}, {}^w{\dot{\mathit{p}}}_s, {\dot{\psi }}_s, {}^w{\dot{\mathit{p}}}_i$ and ${\dot{\psi }}_i$, representing the inaccuracies in our kinematic model.

The term “constant velocity” for a kinematic vehicle model in a Kalman Filter (KF) or Extended Kalman Filter (EKF) refers to the behavior of the process model between observations from navigation sensors including GNSS fixes, modem OWTT ranges, DVL velocities, Attitude and Heading Reference System (AHRS) attitude measurements, and pressure depth sensor measurements. Whenever a new sensor reading (observation) is available, the KF/EKF update innovation adjusts the full estimated state—including both position and velocity. Thus the estimated plant velocity is actually continually varying in response to navigation sensor readings (observations).

3.4.1. Vehicle Process Model

We employ a constant-velocity nonlinear process model for the vehicle

$$\begin{array}{cc}\hfill {\dot{\mathit{x}}}_v& =\underset{\mathit{f}\left({\mathit{x}}_v\left(t\right)\right)}{\underbrace{\left[\begin{array}{cccc}\mathbf{0}& \mathbf{0}& {}_v^w\mathit{R}\left(\mathbf{\phi }\right)& \mathbf{0}\\ \mathbf{0}& \mathbf{0}& \mathbf{0}& \mathcal{E}\left(\mathbf{\phi }\right)\\ \mathbf{0}& \mathbf{0}& \mathbf{0}& \mathbf{0}\\ \mathbf{0}& \mathbf{0}& \mathbf{0}& \mathbf{0}\end{array}\right]{\mathit{x}}_v}}+\underset{{\mathit{G}}_v}{\underbrace{\left[\begin{array}{cc}\mathbf{0}& \mathbf{0}\\ \mathbf{0}& \mathbf{0}\\ \mathit{I}& \mathbf{0}\\ \mathbf{0}& \mathit{I}\end{array}\right]}}{\mathit{w}}_v\hfill \end{array}$$

(7)

where ${}_v^w\mathit{R}\left(\mathbf{\phi }\right)$ is the transformation from body-frame to local-level linear velocities, $\mathcal{E}\left(\mathbf{\phi }\right)$ is the transformation from body-frame angular velocities to Euler rates given by

$$\begin{array}{cc}\hfill \mathcal{E}\left(\mathbf{\phi }\right)& =\left[\begin{array}{ccc}1& sin\varphi tan\theta & cos\varphi tan\theta \\ 0& cos\varphi & -sin\varphi \\ 0& sin\varphi sec\theta & cos\varphi sec\theta \end{array}\right].\hfill \end{array}$$

(8)

and ${\mathit{w}}_v\sim \mathcal{N}(0,{\mathit{Q}}_v)$ is the independent zero-mean Gaussian process noise in the acceleration term.

3.4.2. Ship Process Model

We use a constant-velocity linear process model for the ship, which is reasonable given the high update rates of the ship’s GNSS and gyrocompass in comparison to the comparatively slow ship’s motion

$$\begin{array}{cc}\hfill {\dot{\mathit{x}}}_s& =\underset{{\mathit{F}}_s}{\underbrace{\left[\begin{array}{cc}\mathbf{0}& \mathit{I}\\ \mathbf{0}& \mathbf{0}\end{array}\right]}}{\mathit{x}}_s+\underset{{\mathit{G}}_s}{\underbrace{\left[\begin{array}{c}\mathbf{0}\\ \mathit{I}\end{array}\right]}}{\mathit{w}}_s\hfill \end{array}$$

(9)

where ${\mathit{w}}_s\sim \mathcal{N}(0,{\mathit{Q}}_s)$ is the independent zero-mean Gaussian process noise in the acceleration term.

3.4.3. Ice Process Model

Similar to the ship, we use a linear constant-velocity process model for the ice, which is reasonable given the high update rates of the ice’s GNSS beacons, given by

$$\begin{array}{cc}\hfill {\dot{\mathit{x}}}_i& =\underset{{\mathit{F}}_i}{\underbrace{\left[\begin{array}{cc}\mathbf{0}& \mathit{I}\\ \mathbf{0}& \mathbf{0}\end{array}\right]}}{\mathit{x}}_i+\underset{{\mathit{G}}_i}{\underbrace{\left[\begin{array}{c}\mathbf{0}\\ \mathit{I}\end{array}\right]}}{\mathit{w}}_i\hfill \end{array}$$

(10)

where ${\mathit{w}}_i\sim \mathcal{N}(0,{\mathit{Q}}_i)$ is the independent zero-mean Gaussian process noise in the acceleration term.

3.5. Linearization and Descretization

For a small time step, $\Delta t=\left[t_k, t_{k+1}\right)$, the continuous non-linear system in (6) is linearized and discretized about the estimate ${\mu }_k={\mathit{x}}_k$, which results in the linear, discrete-time equation

$${\mathit{x}}_{k+1}={\mathit{F}}_k{\mathit{x}}_k+{\mathit{B}}_k{\mathit{u}}_k+{\mathit{w}}_k$$

(11)

where

$${\mathit{F}}_k=e^{{\mathit{F}}_x\Delta t}$$

(12)

is the state transition matrix with the system’s Jacobian matrix, ${\mathit{F}}_x$ evaluated at ${\mu }_k$. The matrix

$${\mathit{B}}_k=e^{{\mathit{F}}_{\mathit{x}}{t}_{k+1}}{\int }_{t_k}^{t_{k+1}}{e}^{-{\mathit{F}}_x\tau }d\tau$$

(13)

is the discrete-time control gain for the constant pseudo-input control

$${\mathit{u}}_k=\mathit{f}\left({\mu }_k\right)-{\mathit{F}}_x{\mathit{x}}_k.$$

(14)

Discrete-time noise is given by

$${\mathit{w}}_k={\int }_{t_k}^{t_{k+1}}{e}^{-{\mathit{F}}_x\tau }d\tau$$

(15)

with covariance

$${\mathit{Q}}_k={\int }_{t_k}^{t_{k+1}}{e}^{{\mathit{F}}_x(t_{k+1}-\tau )}\mathit{G}\mathit{Q}{\mathit{G}}^{\top }{e}^{{\mathit{F}}_x^{\top }(t_{k+1}-\tau )}d\tau .$$

(16)

In practice, the matrix ${\mathit{B}}_k$ is evaluated numerically by means of Simpson’s composite method [40], due to the singular nature of ${\mathit{F}}_x$. The process covariance matrix, ${\mathit{Q}}_k$ can be approximated by

$${\mathit{Q}}_k\approx \mathit{G}\mathit{Q}{\mathit{G}}^{\top }\Delta t$$

(17)

for small $\Delta t$ [41].

3.6. Sensor Observation Models: Navigation Sensors

3.6.1. Vehicle DVL

We assume a simplified DVL which measures its velocity relative to the point on the ice located at frame p, which is insonified by the DVL. In order to estimate the relative velocity between frames p and D, we estimate their respective velocities as seen from frame w, and project their relative velocities into the D frame.

The DVL’s world linear velocity is given by the linear and angular components rotated into the world frame,

$${}^w{\mathbf{\nu }}_D^w={}_v^w\mathit{R}\left[\mathbf{\nu }+J^{\left(\omega \right)}{}^v{\mathit{p}}_D\right],$$

(18)

where $so\left(3\right)$ is the set of $3\times 3$ skew-symmetric matrices,

$$so\left(3\right)=\{\mathit{S}:\mathit{S}\in {\mathbb{R}}^{3\times 3},\mathit{S}={\mathit{S}}^{\top }\}$$

(19)

and $\mathcal{J}$ is the function that maps a $3\times 1$ vector into the corresponding $3\times 3$ skew-symmetric matrix, $\mathcal{J}:{\mathbb{R}}^3\to so\left(3\right)$. For $\omega \in {\mathbb{R}}^3$,

$$J\left(\mathbf{\omega }\right)=\left[\begin{array}{ccc}0& -{\omega }_3& {\omega }_2\\ {\omega }_3& 0& -{\omega }_1\\ -{\omega }_2& {\omega }_1& 0\end{array}\right].$$

(20)

To estimate the velocity ${}^w\nu {}_p^w$, where $\mathit{p}$ is the frame relative to which the DVL measures velocity, we utilize the DVL’s range measurement to construct ${}^D\mathit{p}={[0,0,r,1]}^{\top }$, assuming r is the DVL-measured range along the instrument’s z-axis. We then project this point into the ice frame

$${}^i\mathit{p}={}_i^w{\mathit{H}}^{-1}{}_v^w\mathit{H}{}_D^v\mathit{H}{}^D\mathit{p}.$$

(21)

Thus

$${}^w\mathbf{\nu }{}_p^w={}^w{\dot{\mathit{p}}}_i+{}_i^w\mathit{R}{\mathcal{J}}^{\left({[0,0,{\dot{\psi }}_i]}^{\top }\right)}{}^i{\mathit{p}}_p,$$

(22)

where ${}^i{\mathit{p}}_p$ is taken from the first three elements of ${}^i\mathit{p}$. The relative velocity between the ice and the DVL, as seen by the DVL can then be computed by

$${}^D\mathbf{\nu }{}_p^D={\left[{}_v^w\mathit{R}{}_D^v\mathit{R}\right]}^{\top }\left[{}^w\mathbf{\nu }{}_p^w{-}^w{\mathbf{\nu }}_D^w\right].$$

(23)

3.6.2. Vehicle and Ship Gyrocompass

We employed a measurement model for the gyrocompass attitude sensor in order to uniformly account for sensor noise in all sensors. We assume a 6-DOF true-North seeking gyrocompass for both the vehicle and ship which measures instrument angular rate and roll, pitch, yaw angles. When mounted to the UUV the instrument’s roll, pitch, and yaw angles are given by

$$\begin{array}{cc}\hfill {}_I^w\mathit{R}& ={}_v^w\mathit{R}{}_I^v\mathit{R}\hfill \\ \hfill \varphi & =\mathrm{atan}2({}_I^w{\mathit{R}}_{3,2},{}_I^w{\mathit{R}}_{3,3})\hfill \\ \hfill \theta & =\mathrm{atan}2(-{}_I^w{\mathit{R}}_{3,1},\sqrt{{}_I^w{\mathit{R}}_{3,2}^2+{}_I^w{\mathit{R}}_{3,3}^2})\hfill \\ \hfill \psi & =\mathrm{atan}2({}_I^w{\mathit{R}}_{2,1},{}_I^w{\mathit{R}}_{1,1}).\hfill \end{array}$$

(24)

The gyrocompass’s angular rate is given by

$${}^I\mathbf{\omega }{}_w^I={}_I^v{\mathit{R}}^{\top }\mathbf{\omega }.$$

(25)

The ship’s gyrocompass model is obtained by substituting ${}_I^s\mathit{R}$ for ${}_I^v\mathit{R}$, and ${[0,0,{\dot{\psi }}_s]}^{\top }$ for $\mathbf{\omega }$.

3.6.3. Vehicle Pressure Sensor

Assuming a constant water density, $\rho$, with appropriate units, the anticipated pressure sensor measurement, $z_p$, is given by

$$z_p=\left(\rho g\left[0,0,1,0\right]{}_v^w\mathit{H}\left[\begin{array}{c}{}^v{\mathit{p}}_P\\ 1\end{array}\right]\right)+z_{p_0},$$

(26)

where g and $z_{p_0}$ are gravitational acceleration and pressure at mean sea-level.

3.6.4. Ship–Vehicle Range Sensor

A One Way Travel Time (OWTT) range sensor measures the distance between two acoustic transceivers, by means of measuring the time of flight of a one-way signal through water of a known density profile. As with the DVL model, our range sensor model assumes instantaneous measurement, and neglects acoustic travel time through the water. With these simplifying assumptions, the anticipated range r is given by the Euclidean norm of ${}^T\mathit{p}{}_B$, where the T and B are the frames corresponding to the ship-mounted acoustic transducer and vehicle-mounted acoustic beacon, respectively, which is given by

$$\begin{array}{cc}\hfill \left[\begin{array}{c}{}^T{\mathit{p}}_B\\ 1\end{array}\right]& ={}_T^s{\mathit{H}}^{-1}{}_s^w{\mathit{H}}^{-1}{}_v^w\mathit{H}\left[\begin{array}{c}{}^v{\mathit{p}}_B\\ 1\end{array}\right]\hfill \\ \hfill r& =\parallel {}^T{\mathit{p}}_B{\parallel }_2.\hfill \end{array}$$

(27)

3.6.5. GNSS

We assume a 3-DOF GNSS measurement, consisting of ${\varphi }_g,{\theta }_g,$ and h which represent geodetic latitude, longitude and altitude (relative to the WGS84 Ellipsoid). The expected GNSS position in Earth-Centered Earth-Fixed (ECEF) coordinates is given by

$$\left[\begin{array}{c}{}^{ECEF}{\mathit{p}}_{GNSS}\\ 1\end{array}\right]{=}_{ECEF}^w{\mathit{H}}^{-1}{}_s^w\mathit{H}\left[\begin{array}{c}{}^s{\mathit{p}}_{GNSS}\\ 1\end{array}\right].$$

(28)

To calculate the expected ${\varphi }_g,{\theta }_g$ and h, we use the geodetic-ECEF conversions in given in ([42] Appendix B).

The appropriate ice-GNSS transformations, state-variables, and position vectors can be substituted into (28) to obtain the expected measurements for the GNSS beacons located on top of the ice.

4. Comparative Numerical Simulation Performance Analysis: Ice-Relative Navigation

We employed the Gazebo simulation environment [43] to generate simulated instrument measurements and ground truth data against which the EKF’s performance was evaluated. The simulation consisted of a 2.1 km × 2.1 km × 2 m contiguous ice floe, ship, and UUV performing a rectangular survey relative to the ice floe. Survey lines were 1 km long, with 100 m spacing. The UUV model, trajectory controller, and pressure sensor plug-in were provided by UUV Simulator [44], existing Gazebo plug-ins were used to generate gyrocompass and GNSS measurements, and custom Gazebo plug-ins were written to simulate DVL and OWTT range measurements.

In these simulations, the simulated underwater vehicle was commanded to follow a 1 km × 0.7 km ice-relative rectangular survey, with 100 m leg spacing at a depth of 10 m with a commanded ice-relative velocity of 1 m/s.

We assume that the ship’s telemetry downlink to the submerged vehicle (via acoustic communications or via a lightweight fiber optic tether) contains regular reports of the estimated state (position, angular position, and associated velocities) of the ice and the ship. With this information, the vehicle can use its on-board sensors of the upward-looking DVL together with its AHRS and depth-sensor to estimate its absolute position and velocity relative to the geode as well as its position and velocity relative to the ice floe, thus enabling the vehicle to perform closed-loop control of the vehicle trajectory as desired to accomplish mission goals. With the proposed navigation system it is possible for the vehicle to execute closed-loop control relative to any desired frame of reference, including the ice-frame, geodetic world frame, or other frames as required to achieve mission objectives. Because the vehicle is expected to operate in closed-loop control, under-ice water currents relative to the ice flow will not perturb the vehicle trajectories.

4.1. Sensor Observation Noise

Zero-mean Gaussian measurement (observation) noise was assumed for each instrument, with noise and update rates listed in Table 1. While a published range accuracy for the Woods Hole Oceanographic Institution (WHOI) Micromodem is unavailable, the modem measures one-way travel-times with a measurement resolution of 125 microseconds (0.1875 m assuming a sound velocity of 1500 m/s), and reports measurements with a decimal resolution of 100 microseconds [45], thus the primary source of acoustic range estimation error is uncertainty or variation of sound velocity in the water column, and, in some cases, the effect of the acoustic path between source and receiver. Based upon our experience in moderate-range OWTT navigation, we simulated OWTT range-measurement noise to be zero mean with a standard deviation of 0.2% of actual slant-range between transmitter and receiver. The data sheet for the RDI Workhorse Navigator series [7] specifies them to have excellent long-term accuracy (bias) within +/−2 mm/s for the 300 kHz units, and +/−1 mm/s for the 600 kHz and 1200 kHz units when at low velocities, thus we consider the bias terms to be negligible. The IXSEA Phins true-North-seeking Fiber-optic Gyroscope (FOG) and Inertial Measurement Unit (IMU) has internal compensation for accelerometer and angular rate sensor biases [46].

Table 1. Simulated sensor update rates and noise statistics.

Instrument	Model	Measurement	$\mathit{\sigma }$	Update Rate
Gyro-compass	IXSEA	Heading	0.1${}^{°}sec\left(lat\right)$	10 Hz
	PHINS [46]	Pitch, Roll	0.01${}^{°}$
DVL	Teledyne RDI	Velocity	3 mm/s	2 Hz
	300 kHz [7]	Range	1% Range
Pressure	Paroscientific	Pressure	0.01%	2 Hz
Sensor	Digiquartz [47]
Ice GNSS	Garmin 18x LVC [35]	Position	3 m	1 Hz
Ship GNSS	Trimble SPS852 [36]	Position	0.3 m	1 Hz
Range	WHOI	UUV-Ship	0.2% Slant	1/60 Hz
Sensor	Micromodem	Range	Range

Note: Garmin 18x position accuracy value computed from 95% CEP values

Table 1. Simulated sensor update rates and noise statistics.

Instrument	Model	Measurement	$\mathit{\sigma }$	Update Rate
Gyro-compass	IXSEA	Heading	0.1${}^{°}sec\left(lat\right)$	10 Hz
	PHINS [46]	Pitch, Roll	0.01${}^{°}$
DVL	Teledyne RDI	Velocity	3 mm/s	2 Hz
	300 kHz [7]	Range	1% Range
Pressure	Paroscientific	Pressure	0.01%	2 Hz
Sensor	Digiquartz [47]
Ice GNSS	Garmin 18x LVC [35]	Position	3 m	1 Hz
Ship GNSS	Trimble SPS852 [36]	Position	0.3 m	1 Hz
Range	WHOI	UUV-Ship	0.2% Slant	1/60 Hz
Sensor	Micromodem	Range	Range

Note: Garmin 18x position accuracy value computed from 95% CEP values

4.2. Process Noise

As is usual for kinematic vehicle plant models (e.g., [33,34,39]), we assumed zero-mean Gaussian process noise Process noise values were tuned empirically to minimize ice-relative vehicle error, with chosen process noise values shown in

Table 2. Kinematic process model noise statistics.

State Variable	Associated Process Noise (${\mathit{\sigma }}^2$)	Units
$\dot{\nu }$	${[0.05,0.005,0.05]}^{\top }$	${\left(\frac{m}{s^2}\right)}^2$
$\dot{\omega }$	${[0.02,0.02,0.02]}^{\top }$	${\left(\frac{rad}{s^2}\right)}^2$
${}^w{\dot{\mathit{p}}}_s$	${[0.005,0.005,0.005]}^{\top }$	${\left(\frac{m}{s}\right)}^2$
${\dot{\psi }}_s$	$0.002$	${\left(\frac{rad}{s}\right)}^2$
${}^w{\dot{\mathit{p}}}_i$	${[0.00005,0.00005,0.00005]}^{\top }$	${\left(\frac{m}{s}\right)}^2$
${\dot{\psi }}_i$	$0.00002$	${\left(\frac{rad}{s}\right)}^2$

.

Most of the oceanographic underwater vehicles, including Nereid Under-Ice (NUI) hybrid underwater vehicle [21], are designed to have low drag in the surge direction (for forward motion), and low drag in the heave direction (for vertical descent and ascent), and have a relatively high drag in the sway direction. These process noise statistics correspond to the relative mobility of such a vehicle based upon its directional drag differences.

4.3. Simulated Mission

A 1 km × 0.7 km ice-relative rectangular survey, with 100 m leg spacing was simulated with the UUV at a depth of 10 m, with a commanded ice-relative velocity of 1 m/s. The UUV trajectory controller utilized ground-truth ice-relative position and velocity information. An ice-relative survey path is useful, for example, for conducting oceanographic surveys of sea ice physical properties over a large area, such as that reported in [32] and making upward-looking multibeam maps of the complex underwater topography of multi-year sea ice [22].

The ship was modeled as being rigidly attached to the ice floe on its port side. Based upon the Authors’ experience, this is the most common configuration employed by icebreakers when rafting to a large ice floe because it allows the ship to maintain stable contact with the floe while simultaneously clearing a small opening in the ice on the starboard side, adjacent to the ship’s oceanographic winches and cranes, to launch and recover vehicles and instruments.

4.4. Sensitivity Analysis

To evaluate the state estimator’s performance and sensitivity to several controllable instrumentation decisions, we examined the individual effects of the following parameters: $\left(i\right)$ variation in spacing of ice-deployed GNSS beacons, $\left(ii\right)$ variation in the number of ice-deployed GNSS beacons, $\left(iii\right)$ variation in spacing between ship-deployed OWTT range transducers, and $\left(iv\right)$ variation in translational and rotational ice velocity on the accuracy of the proposed method. While ice velocities are an uncontrollable environmental factor, they nevertheless present an additional parameter that is likely to vary under realistic operating conditions. The two metrics used to investigate performance of the EKF under the simulated parameters were average ice-relative Root Mean Square (RMS) vehicle position error, and average RMS ice-relative vehicle x-y position uncertainty, as estimated by the EKF. These values were averaged over the entirety of the simulated missions.

4.4.1. GNSS Configuration

To examine the effect spacing and the number of GNSS beacons on the estimator, the simulated ice floe was instrumented with 30 distinct configurations of ice-mounted GNSS beacons. In each configuration the first beacon, GNSS${}_1$, marks the origin of the ice frame and was aligned with the stern of the ship, and 5 m off the port side. The position of the remaining GNSS units are depicted in Figure 1 with edge spacing varying from 100–1000 m. These experiments correspond to exp1, exp2, and exp3 in

Table 3. Parameters under variation in the sensitivity analysis included number (and spacing) of ice-deployed GNSS units, ice velocity, and OWTT transducer spacing.

	(x, y) Locations of GNSS in Ice Frame			Ice Velocity	OWTT Baseline, b
	GNSS${}_2$	GNSS${}_3$	GNSS${}_4$
exp1	(100n,0)	---	---	0.25 m/s 2${}^{°}$/h	40 [m]	$n=1,\dots ,10$
exp2	(100n,0)	(100n, −100n)	---	0.25 m/s 2${}^{°}$/h	40 [m]
exp3	(100n,0)	(100n, −100n)	(−100n,0)	0.25 m/s 2${}^{°}$/h	40 [m]
exp4	(500,0)	---	---	0.125n m/s 1n${}^{°}$/h	40 [m]	$n=0,\dots ,4$
exp5	(500,0)	---	---	0.25 m/s 2${}^{°}$/h	20n [m]	$n=1,\dots ,5$

.

4.4.2. Ice Velocity

Our experience is that floe velocities are highly variable depending on ice conditions, wind, and currents, but that it is common to observe ice translational velocities in the 0–0.5 m/s range and rotational velocities of 0–4 deg/h [21]. A total of five linear and angular ice floe velocities are simulated, in separate simulation runs. These experiments correspond to exp4 in Table 3.

4.4.3. One Way Travel Time (OWTT) Baseline

To examine the effect of OWTT transducer spacing, the simulated baseline was varied from 20–100 m. Previous deployments in the Arctic reported using a baseline of approximately 40 m [21]. On some icebreakers, the deck layout makes it difficult to deploy modem transducers on a baseline that is more than a fraction of the actual ship length. The range of acoustic baselines chosen for the simulations was selected to fully cover the range of baselines we have experienced in actual Arctic deployments.

Simulated range measurements were taken on an alternating 120 s cycle, with the stern modem transducer making measurements at the top of the cycle and the bow modem transducer making measurements 60 s later. These experiments correspond to exp5 in Table 3.

Because of the size of the parameter space under exploration, we limit our examination to a one-parameter variation at a time. We assume a nominal configuration of two ice-based GNSS units, spaced 500 m apart, 40 m spacing between ship-based OWTT transducers, and a constant linear and angular ice velocity of 0.25 m/s and 2 deg/h. A summary of the parameters under examination during each set of experiments is provided in Table 3.

5. Results: Ice-Relative Navigation

Five distinct simulations were conducted, corresponding to the five sets of linear and angular ice velocities outlined in Table 3. Simulated measurements from all 30 ice GNSS beacon configurations, 6 OWTT transducers, and the vehicle and ship sensors outlined in Table 1 were saved for post-processing. Measurements from the selected GNSS and OWTT units, were then used in conjunction with the process and observation models outlined in Section 3 to estimate the state of the system using the EKF according to Table 3. Process noise values used in the analysis are detailed in Table 2, and were tuned based on empirical assessment of vehicle position estimates. Position covariances were initialized as follows: UUV $x-y-z$: 10 m${}^2$, ship $x-y$: 1 m${}^2$, and ice $x-y$ covariance of 2 m². Initial position states were initialized to their true values, and velocities and accelerations were initialized at zero. Due to the limitations of the simulation environment, data could not be generated faster than real-time, and thus it was not feasible to conduct a Monte Carlo-style experiment of repeated simulation runs.

Figure 2a shows the true (black) and estimated (red) vehicle tracklines in the world frame, with 1$\sigma$ vehicle $x-y$ covariance ellipsoids (blue) estimated from exp1 using 500 m GNSS spacing. Given the estimated vehicle and ice poses in the world frame, we can project the vehicle’s estimated position and orientation into the ice frame,

$${}_v^i\mathit{H}={}_i^w{\mathit{H}}^{-1}{}_v^w\mathit{H}.$$

(29)

Using the method outlined in [48] we can then propagate the vehicle and ice linear and angular position covariance by means of first order approximation into the ice frame,

$$\mathrm{\Sigma }={\mathit{J}}_i{\mathrm{\Sigma }}_i{\mathit{J}}_i^{\top }+{\mathit{J}}_i{\mathrm{\Sigma }}_v{\mathit{J}}_i^{\top }-{\mathit{J}}_i{\mathrm{\Sigma }}_{iv}{\mathit{J}}_i^{\top }-{\mathit{J}}_i{\mathrm{\Sigma }}_{vi}{\mathit{J}}_i^{\top },$$

(30)

where ${\mathit{J}}_i\in {\mathbb{R}}^{6\times 6}$ is the Jacobian matrix corresponding to the ${}_i^w\mathit{H}$ transformation, ${\mathrm{\Sigma }}_i$ and ${\mathrm{\Sigma }}_v$$\in {\mathbb{R}}^{6\times 6}$ are the 6-DOF position and orientation covariance matrices of the ice and UUV, and ${\mathrm{\Sigma }}_{iv}$ is the off-diagonal portion of the state covariance matrix corresponding to the covariance between the 6-DOF ice and vehicle poses. Figure 2b shows the estimated and actual ice-relative trajectory, with first-order propagated covariance ellipsoids.

5.1. The Effects of Variation in GNSS Configuration

Estimator sensitivity to GNSS receiver placement on the ice is of natural interest because it is a variable which, given reasonable ice and weather conditions, is configurable by a deployment team prior to vehicle deployment. In addition to GNSS receiver spacing, estimator sensitivity to the number of deployed GNSS units is intuitively of interest as a means to provide additional independent observations of ice state.

Figure 3 shows the average estimated ice-relative vehicle $x-y$ RMS position error for the 7.7 km simulated mission. For each GNSS configuration (line, triangle, square, corresponding to exp1, exp2, exp3 respectively), we note that increased sensor separation is correlated with decreased ice-relative vehicle RMS position error.

Figure 4 shows ice-relative vehicle $x-y$ position uncertainty, computed using (30), averaged over the 7.7 km simulated mission for exp1, exp2, and exp3. Estimator uncertainty is shown component-wise in the ice frame, where the x and y components correspond to along-track and cross-track uncertainty respectively along the 1 km survey lines. The trend is similar to that of Figure 3 wherein increasing number and baseline of GNSS units is correlated with decreasing position uncertainty. Of note is that the reduction in estimated ice-relative vehicle position uncertainty by means of additional GNSS units is most pronounced at smaller spacing.

The trends observed in Figure 3 and Figure 4 clearly show a monotonic reduction in ice-relative vehicle $x-y$ position error with increasing GNSS spacing in all three experimental configurations, and are representative of the performance we observed throughout our study.

5.2. The Effects of Variation in Ice Velocity

For the purpose of examining ice-velocity on estimator performance, the ice-deployed GNSS configuration and OWTT transducer spacing were held constant at 500 m (line orientation), and 40 m respectively. Figure 5 shows the ice-relative vehicle $x-y$ position RMS error, averaged over each simulated mission. Position RMS error appears to increase slightly with increasing ice velocity, with the exception of the 0.25 m/s, 2${}^{°}$/h case which shows slightly reduced RMS error. Such an outlier could be stochastic in nature, or could be due to the fact that the filter’s process noise was tuned empirically using the using results from the 0.25 m/s, $2^{°}$/h case with 500 m GNSS spacing, and 40 m OWTT beacon spacing.

Figure 6 shows ice-relative vehicle $x-y$ position uncertainty averaged over each simulated mission. While increasing linear and angular ice velocity had little effect on RMS error, estimator uncertainty appears to grow monotonically with increased speed.

5.3. The Effects of Variation in OWTT Baseline

Spacing between ship-deployed OWTT transducers, the final parameter examined, is a variable that is typically easily adjustable prior to deployment of a vehicle. Increased transducer separation provides favorable geometry when ranging to a vehicle at smaller standoff distances, but this advantage decreases at large distances from the ship. For the final analysis, a 500 m GNSS line configuration and constant 0.25 m/s, 2 deg/h ice velocity were assumed while transducer spacing was varied from 20–100 m in 20 m increments.

Figure 7 shows average ice-relative vehicle $x-y$ position RMS over the range of simulated transducer spacing. Examining Figure 7, RMS error appears effectively insensitive to transducer spacing.

Figure 8 shows ice-relative vehicle $x-y$ position uncertainty averaged over the full vehicle trajectory for the five simulated transducer spacing. Similar to ice-relative position RMS error, the estimators uncertainty appears insensitive to transducer spacing when averaged over the full trajectory.

Figure 9 shows the magnitude of ice-relative position error averaged over the vehicle trajectory as a function of the number of GNSS receivers used (exp1–exp3). The general trend shows a correlation between the increased number of GNSS receivers, and reduced ice-relative position error, with the exception of the 100m case - possibly an artifact of the stochastic nature of the simulation.

6. Conclusions

This paper reports the results of a sensitivity analysis of an Extended Kalman Filter for use in navigation of underwater vehicles beneath moving sea ice using simulated sensor measurements. The effects on ice-relative vehicle $x-y$ navigation position RMS error and uncertainty are examined over a range of ice-deployed GNSS spacing and configurations, varying translational and rotational ice floe velocities, and ship-deployed OWTT transducer spacing. The data suggest that increasing the number and spacing of ice-deployed GNSS beacons reduces average RMS error and position uncertainty. While ice-relative vehicle position error appears relatively insensitive to increasing linear and angular ice velocity, position uncertainty (as estimated by the EKF) appears to increase monotonically with increased ice velocity. Lastly the simulations suggest that both ice-relative RMS error and uncertainty appear unaffected by varying ship-board OWTT transducer spacing (baseline) over the 20–100 m range.

These simulation results suggest that instrumenting ice floes with two or more GNSS receivers, along with other precision underwater navigation instrumentation can provide a scientifically useful means of navigation beneath moving sea ice over a range of operating conditions and vehicle-ship standoff distances.

Barring the opportunity to implement such a navigation system in the field, future studies could better model underwater acoustics, namely non-instantaneous time of flight, ray bending, and under-ice acoustical reflections.