In the ocean domain, opportunities for a paradigm shift in the science of autonomy involve fundamental theory, rigorous methods, and efficient computations for autonomous systems that collect information, learn, collaborate and make decisions under uncertainty, all in optimal integrated fashion and over long duration, persistently adapting to and sustainably utilizing the ocean environment. The ocean is a prime example of multiscale nonlinear multidisciplinary and strongly coupled dynamics where measurements of one variable can be used to infer other fields through their joint dynamic probability density functions. Integrating ocean dynamics with autonomy enables the principled exploration, sustainable utilization, and strong conservation of our oceans. The object of this thesis is to develop theory, algorithms, and computational systems for the high dimensional optimal path planning of autonomous vehicles in the physical space augmented with other dynamical fields, and for the Bayesian nonlinear assimilation of the observations gathered by these vehicles along their trajectory. The resulting high dimensional optimal path planning and generalized Lagrangian Bayesian data assimilation enable the sustained and optimal operation of autonomous vehicles over a long time duration in realistic uncertain ocean settings. With this vision, the vehicles autonomously make decisions to optimally achieve their mission targets in the augmented space of physical and collectible fields, e.g., reach the destination in minimum time while using the minimum energy, harvest the maximum wave-solar-wind energy, or farm the maximum amount of kelp. To this end, we focus on three specific theoretical and methodological goals: (i) Develop exact differential equations and accurate algorithms to efficiently predict and compute the reachable sets and optimal controls for complex high-dimensional objectives in dynamical fields and showcase these controls in realistic ocean scenarios; (ii) Develop Bayesian theory and schemes to predict Lagrangian field probability densities and rigorously assimilate Lagrangian data collected by moving vehicles or drifters, leveraging the different sensitivity and timescales of the underlying Lagrangian and Eulerian systems; and (iii) Integrate the optimal planning and assimilation to enable learning from the information gathered on-board the vehicles, from Bayesian updates of the optimal controls to the acquisition of knowledge using Bayesian learning. We showcase the theory and methods in a range of ocean applications.

In the first part, we review the theory and schemes to predict joint energy-time, harvesting-time, and energy-harvesting-time reachability fronts and optimal paths using state augmentation. We validate our energy-time algorithms in analytical and representative dynamical fields. We then derive theory to predict reachability fronts across multiple times simultaneously and obtain a closed loop control law allowing vehicles to accomplish their mission even after straying from their initial plan due to forecast errors. The theory and schemes are developed for both backward and forward reachable tubes with time-varying target and start sets. The resulting value functions elegantly capture not only the reachable tubes but also time-to-reach and time-to-leave maps as well as start time versus duration maps. We validate results with analytical solutions and demonstrate wider applications for optimal control in dynamic environments.

In the second part, we develop and implement fundamental schemes for multi-timescale Bayesian data assimilation for coupled dynamical systems, with a focus on Lagrangian-Eulerian systems. We obtain a Gaussian Mixture Model (GMM) – Dynamically Orthogonal (DO) based hybrid filter for Lagrangian and Eulerian stochastic fields and observations. We first showcase the schemes for a coupled system where we analytically validate the performance of the filter. We subsequently demonstrate results by applying the filter for a general coupled chaotic system and for a joint Lagrangian-Eulerian system with a more complex quasi-geostrophic flow.

In the third part, we integrate the schemes developed for the first two parts. We propose coupled methods that allow ocean vehicles to robustly and optimally complete their mission while continuously learning from the new information being collected, updating the Lagrangian and Eulerian fields, their joint probabilities, and the robust optimal control of their future trajectories. We showcase preliminary results using the proposed method.

We conclude by demonstrating several planning and Lagrangian algorithms in data-assimilative ocean simulations and real-time ocean experiments with real data and forecasts. This includes the characterization of residence times and connectivity in the Red Sea, the transport of plastics in the coastal ocean showcasing results for Massachusetts Bay, the subduction pathways of surface waters to intermediate depths in the Alboran Sea, and the Bayesian Eulerian-Lagrangian data assimilation of drifter data in the Balearic Sea.

Sound waves are critical for a variety of underwater applications including communication, navigation, echo-sounding, environmental monitoring, and marine biology research. However, the incomplete knowledge of the ocean environment and acoustic parameters makes reliable acoustic predictions challenging. This is due to the sparse and heterogeneous data, as well as to the complex ocean physics and acoustics dynamics, multiscale interactions, and large dimensions. There are thus several sources of uncertainty in acoustic predictions. They include the ocean current, temperature, salinity, pressure, density, and sound speed fields, the bottom topography and seabed fields, the sources and receivers properties, and finally the model equations themselves. The goals of this thesis are to address these challenges. Specifically, we: (1) obtain, solve, and verify differential equations for efficient probabilistic underwater acoustic modeling in uncertain environments; (2) develop theory and implement algorithms for the Bayesian nonlinear inference and learning of the ocean, bathymetry, seabed, and acoustic fields and parameters using sparse data; and (3) demonstrate the new methodologies in a range of underwater acoustic applications and real sea experiments, showcasing new capabilities and leading to improved understanding.

In the first part, we derive, discretize, implement, and verify stochastic differential equations that (i) capture dominant input uncertainties in the environment (e.g., ocean, bathymetry, and seabed) and in the acoustic parameters (e.g. source location, frequency, and bandwidth), and (ii) predict the acoustic pressure fields and their probability distributions, respecting the nonlinear governing equations and non-Gaussian statistics. Starting from the acoustic Parabolic Equation (PE), we develop Dynamically Orthogonal (DO) differential equations for range-optimal acoustic uncertainty quantification. Using DO expansions for the input uncertainties, we develop the reduced-order DO-PEs theory for the Narrow-Angle PE (NAPE) and Pad\’e Wide-Angle PE (WAPE) stochastic partial differential equations (PDEs). We verify the discretized DO-PEs in new stochastic range-independent and range-dependent test cases, and demonstrate their advantages over state-of-the-art methods for uncertainty quantification and wave propagation in random media. Results show that a single DO-PE simulation can accurately predict stochastic range-dependent acoustic fields and their full non-Gaussian probability distributions, with computational savings of several orders of magnitude when compared to direct Monte Carlo methods.

In the second part, we extend recent nonlinear Bayesian data assimilation (DA) to the inference and learning of ocean-bathymetry-seabed-acoustic fields and parameters using sparse acoustic and oceanographic data. We combine the acoustic DO-PEs with Gaussian mixture models (GMMs) to predict probability densities in the DO subspace, allowing for efficient non-Gaussian estimation of state variables, parameters, and model functions themselves. The joint multidisciplinary estimation is enabled by state augmentation where the ocean-acoustic-bathymetry-seabed states and parameters are fit together to GMMs within the DO subspace. The new GMM-DO ocean acoustic inference system is validated by assimilating sparse data to infer the source depth, source frequency, and acoustic and environment fields and parameters in five new high-dimensional inference test cases based on state-of-the-art oceanographic and geoacoustic benchmarks. We evaluate the convergence to inference parameters and quantify the learning skill. Results show that our PDE-based Bayesian learning successfully captures non-Gaussian statistics and acoustic ambiguities. Using Bayes’ law, it provides accurate probability distributions for the multivariate quantities and enables principled learning from noisy, sparse, and indirect data.

In the final part, we integrate our acoustic DO-PEs and GMM-DO frameworks with the MSEAS primitive equation ocean modeling system to enable unprecedented probabilistic forecasting and learning of ocean physics and acoustic pressure and transmission loss (TL) fields, accounting for uncertainties in the ocean, acoustics, bathymetry, and seabed fields. We demonstrate the use of this system for low to mid-frequency propagation with real ocean data assimilation in three regions. The first sea experiment takes place in the western Mediterranean Sea where we showcase the system’s performance in predicting ocean and acoustic probability densities, and assimilating sparse TL and sound speed data for joint ocean physics-acoustics-source depth inversion in deep ocean conditions with steep ridges. In the second application, we simulate stochastic acoustic propagation in Massachusetts Bay around Stellwagen Bank and use our GMM-DO Bayesian inference system to assimilate TL data for acoustic and source depth inversion in shallow dynamics with strong internal waves. Finally, in the third experiment in the New York Bight, we employ our system as a novel probabilistic approach for broadband acoustic modeling and inversion. Overall, our results mark significant progress toward end-to-end ocean-acoustic systems for new ocean exploration and management, risk analysis, and advanced operations.

With the expanding availability of computational power, numerical modeling plays an increasingly pivotal role in the field of oceanography, enabling scientists to explore and understand ocean processes which are otherwise inaccessible or challenging to observe directly. It provides a crucial tool for investigating a range of phenomena from large-scale circulation patterns to small-scale turbulence, shaping our understanding of marine ecosystems, global climate, and weather patterns. However, this same wide range of spatiotemporal scales presents a distinct computational challenge in capturing physical interactions extending from the diffusive scale (millimeters, seconds) to planetary length scales spanning thousands of kilometers and time scales spanning millennia. Therefore, numerical and parameterization improvements have and will continue to define the state of the art in ocean modeling, in tandem with the integration of observational data and adaptive methods. As scientists strive to better understand multiscale ocean processes, the thirst for comprehensive simulations has proceeded apace with concomitant increases in computing power, and submesoscale resolutions where nonhydrostatic effects are important are progressively becoming approachable in ocean modeling. However, few realistic ocean circulation models presently have nonhydrostatic capability, and those that do overwhelmingly use low-order finite-difference and finite-volume methods, which are plagued by dispersive errors, and are arduous to utilize in general, especially on unstructured domains and in conjunction with adaptive numerical capabilities. High-order discontinuous Galerkin (DG) finite element methods (FEMs) allow for arbitrarily high-order solutions on unstructured meshes and often out-compete low-order models with respect to accuracy per computational cost, providing significant reduction of dispersion and dissipation errors over long-time integration horizons. These properties make DG-FEMs ideal for the next generation of ocean models, and, in this thesis, we develop a novel DG-FEM ocean model with the above longer-term vision and adaptive multiscale capabilities in mind.

Using a novel hybridizable discontinuous Galerkin (HDG) spatial discretization for both the hydrostatic and nonhydrostatic ocean equations with a free surface, we develop an accurate and efficient high-order finite element ocean model. We emphasize the stability and robustness properties of our schemes within a projection method discretization. We provide detailed benchmarking and performance comparisons for the parallelized implementation, tailored to the specifics of HDG finite element methods. We demonstrate that the model achieves optimal convergence, and is capable of accurately simulating nonhydrostatic behavior. We evaluate our simulations in diverse dynamical regimes including linear gravity waves, internal solitary waves, and the formation of Rayleigh-Taylor instabilities in the mixed layer. Motivated by investigating local nonhydrostatic submesoscale dynamics using realistic ocean simulation data, we develop schemes to initialize and nest the new DG-FEM model within a comprehensive hydrostatic ocean modeling system. Nested within such data-assimilative hydrostatic simulations in the Alboran Sea, we provide a demonstration of our new model’s ability to capture both hydrostatic and nonhydrostatic dynamics that arise in the presence of wind-forced instabilities in the upper ocean layers. We show that such a model can both validate and work in tandem with larger hydrostatic modeling systems, enabling multi-dynamics simulations and enhancing the predictive fidelity of ocean forecasts.

Next, as DG-FEM methods are well-suited to adaptive refinement, we develop a method to learn new adaptive mesh refinement strategies directly from numerical simulation by formulating the adaptive mesh refinement (AMR) process as a reinforcement learning problem. Finite element discretizations of problems in computational physics can usefully rely on adaptive mesh refinement to preferentially resolve regions containing important features during simulation. However, most spatial refinement strategies are heuristic and rely on domain-specific knowledge or trial-and-error. We treat the process of adaptive mesh refinement as a local, sequential decision-making problem under incomplete information, formulating AMR as a partially observable Markov decision process. Using a deep reinforcement learning (DRL) approach, we train policy networks for AMR strategy directly from numerical simulation. The training process does not require an exact solution or a high-fidelity ground truth to the partial differential equation (PDE) at hand, nor does it require a pre-computed training dataset. The local nature of our deep reinforcement learning approach allows the policy network to be trained inexpensively on much smaller problems than those on which they are deployed, and the DRL-AMR learning process we devise is not specific to any particular PDE, problem dimension, or numerical discretization. The RL policy networks, trained on simple examples, can generalize to more complex problems and can flexibly incorporate diverse problem physics. To that end, we apply the method to a range of PDEs relevant to fluid and ocean processes, using a variety of high-order discontinuous Galerkin and hybridizable discontinuous Galerkin finite element discretizations. We show that the resultant learned policies are competitive with common AMR heuristics and strike a favorable balance between accuracy and cost such that they often lead to a higher accuracy per problem degree of freedom, and are effective across a wide class of PDEs and problems.

Advances in computational power have brought the possibility of realistically modeling our world with numerical simulations closer than ever. Nevertheless, our appetite for higher fidelity simulations and faster run times grows quickly; we will always grasp at what is just beyond our computational reach. No matter if we seek to understand biological, chemical, or physical systems, the bottleneck of scientific computing is almost always the same: high dimensionality. The goal of reduced-order modeling is to reduce the number of unknowns in a system while minimizing the loss of accuracy in the approximate solution. Ideal model-order reduction techniques are optimal compromises between computational tractability and solution fidelity. While there are plenty of such techniques to choose from, their widespread adoption remains to be seen due to several persistent challenges. Many methods are intrusive and difficult to implement, lack traditional numerical guarantees such as convergence and stability, and cannot adapt to unforeseen dynamics. We seek to promote the adoption of reduced-order models (ROMs) by creating non-intrusive, efficient, and dynamically adaptive algorithms that maintain the essential features and numerical guarantees of their full-order counterparts.

In this thesis, we derive and apply algorithms for dynamical reduced-order models. Many model-order reduction approaches project dynamical systems onto a fixed subspace obtained from either a simplification of the original equations, a set of known functions such as orthogonal polynomials, or a reduced basis of full-order simulations computed offline. However, if the true system exits the span of the prescribed subspace, such approaches quickly accumulate large errors. In contrast, dynamical ROMs adapt their subspaces as the system evolves. Geometrically, this amounts to integrating a dynamical system along a nonlinear manifold embedded in a full-order Euclidean space. We develop schemes that not only change subspaces at each discrete time step, but that change the subspace in between time steps for improved accuracy. Even further, our numerical schemes automatically detect when the dynamics depart the nonlinear manifold and may jump to a new nonlinear manifold that better captures the system state. For concreteness, we focus on a reduced-order modeling technique called the dynamical low-rank approximation (DLRA), a discrete analogue to the dynamically orthogonal (DO) differential equations. The DLRA evolves a low-rank system in time (or range) as an approximation to a full-rank system, and in contrast to many methods, the DLRA does not require an offline stage where full-order simulations are computed. It is also agnostic to the source of high dimensionality, whether it be the high resolution required, the large domain, or the stochasticity of the problem. These features make it a versatile tool suitable for a wide variety of problems. We evaluate, verify, and apply our new dynamical reduced-order models and schemes to a varied set of dynamical systems, including stochastic fluid flows and waves, videos and their dynamic compression, realistic ocean acoustics and underwater sound propagation with dynamic coordinate transforms, and stochastic reachability and time-optimal path planning.

The majority of this work is devoted to new adaptive integration schemes for the DLRA. We start by introducing perturbative retractions, which map arbitrary-rank matrices back to a manifold of fixed-rank matrices. They asymptotically approximate the truncated singular value decomposition at a greatly reduced cost while guaranteeing convergence to the best low-rank approximation in a fixed number of iterations. From these retractions, we develop the dynamically orthogonal Runge-Kutta (DORK) schemes, which change the subspace onto which the system’s dynamics are projected in between time steps. The DORK schemes are improved by using stable, optimal (so) perturbative retractions, resulting in the so-DORK schemes. They are more efficient, accurate, and stable than their predecessors. We also introduce gradient-descent (gd) retractions and the gd-DORK schemes, which tend to converge rapidly to the best low-rank approximation by recursively applying retractions. The DORK schemes may be made rank-adaptive and robust to rank overapproximation with either a pseudoinverse or by changing the gauge of the integration scheme. While the pseudoinverse technique accumulates slightly more error, it preserves mode continuity, a feature that changing the gauge lacks. Next, we derive an alternating-implicit (ai) linear low rank solver, which is used to create ai-DORK schemes. The ai-DORK schemes are a general-purpose family of implicit integration schemes that have the same algorithmic complexity as explicit schemes (provided some conditions on the dynamics), which vastly broadens the scope of problems that can be solved with the DLRA. This relieves stringent time-step restrictions and enables the DLRA to handle stiff systems. Furthermore, we develop a piecewise polynomial approximation using adaptive clustering in order to handle non-polynomial nonlinearities in reduced-order models. We thoroughly
test these numerical schemes on well-conditioned and ill-conditioned matrix differential equations; data-driven dynamical systems including videos; Schrödinger’s equation; a stochastic, viscous Burgers’ equation; a deterministic, two-dimensional, viscous Burgers’ equation; an advection-diffusion partial differential equation (PDE); a nonlinear, stochastic Fisher-KPP PDE; nonlinear, stochastic ray tracing; and a nonlinear, stochastic Hamilton-Jacobi-Bellman PDE for time-optimal path planning. We find that the reduced-order solutions may be made arbitrarily accurate using rank adaptive dynamical schemes that automatically track the true rank of the full-order simulation, and nonlinearities may be well-approximated by dynamically increasing the number of stochastic clusters, all at a greatly reduced computational cost.

In addition to DORK schemes, we create a tailor-made low-rank integration scheme for the narrow-angle parabolic wave equation called the low-rank split-step Fourier method. Acoustic simulations are often bottlenecked by the Nyquist criterion, which insists that we sample spatially at least twice per wavelength. To address this, our low-rank split-step Fourier method has an algorithmic complexity that scales sublinearly in the number of classical degrees of freedom, enabling vastly larger computational domains and higher frequencies. We demonstrate its efficacy on realistic ocean acoustics problems in Massachusetts Bay with sound speed fields obtained from our high-resolution ocean primitive equations modeling system. In comparing the low rank and full-rank simulations, we demonstrate that the dynamical low-rank method captures the full-rank features including three-dimensional acoustic energy propagation in complex ocean fields with internal waves and rapidly varying bathymetry.

Lastly, with tools from machine learning, we introduce learnable and automatically differentiable coordinate transforms. The compressibility of a system heavily depends on the choice of coordinates, and frequently a coordinate system is chosen for its simplicity rather than its efficiency. Our novel coordinate transforms are determined in a hands-off manner by minimizing a cost function that includes the environmental data expressed in terms of the non-constant coefficients and initial conditions of a PDE. Not only do we automatically obtain Jacobians and Hessians of the transforms, we also find coordinate systems that reduce the rank of solutions to PDEs. This improves the accuracy of the DLRA for the same cost as a typical low-rank simulation, and it accelerates the convergence in rank to the full-order solution. The coordinate transforms also enable low-rank domain decomposition, which is particularly useful in ocean acoustics where the water-seabed interface is discontinuous. We demonstrate this methodology on a first-order PDE with advection and a second-order PDE, the parabolic wave equation, using two examples. We first show acoustic propagation along a three-dimensional wedge and compare the accuracy of solutions computed in the original and transformed coordinate systems. We then show acoustic propagation in a realistic ocean environment over Stellwagen Bank in Massachusetts Bay with a
dynamic coordinate transform.

Complex dynamical models are used for prediction in many domains, and are useful to mitigate many of the grand challenges being faced by humanity, such as climate change, food security, and sustainability. However, because of computational costs, complexity of real-world phenomena, and limited understanding of the underlying processes involved, models are invariably approximate. The missing dynamics can manifest in the form of unresolved scales, inexact processes, or omitted variables; as the neglected and unresolved terms become important, the utility of model predictions diminishes. To address these challenges, we develop and apply novel scientific machine learning methods to learn unknown and discover missing dynamics in models of dynamical systems.

In our Bayesian approach, we develop an innovative stochastic partial differential equation (PDE) – based model learning framework for high-dimensional coupled biogeochemical-physical models. The framework only uses sparse observations to learn rigorously within and outside of the model space as well as in that of the states and parameters. It employs Dynamically Orthogonal (DO) differential equations for adaptive reduced-order stochastic evolution, and the Gaussian Mixture Model-DO (GMM-DO) filter for simultaneous nonlinear inference in the augmented space of state variables, parameters, and model equations. A first novelty is the Bayesian learning among compatible and embedded candidate models enabled by parameter estimation with special stochastic parameters. A second is the principled Bayesian discovery of new model functions empowered by stochastic piecewise polynomial approximation theory. Our new methodology not only seamlessly and rigorously discriminates between existing models, but also extrapolates out of the space of models to discover newer ones. In all cases, the results are generalizable and interpretable, and associated with probability distributions for all learned quantities. To showcase and quantify the learning performance, we complete both identical-twin and real-world data experiments in a multidisciplinary setting, for both filtering forward and smoothing backward in time. Motivated by active coastal ecosystems and fisheries, our identical twin experiments consist of lower-trophic-level marine ecosystem and fish models in a two-dimensional idealized domain with flow past a seamount representing upwelling due to a sill or strait. Experiments have varying levels of complexities due to different learning objectives and flow and ecosystem dynamics. We find that even when the advection is chaotic or stochastic from uncertain nonhydrostatic variable-density Boussinesq flows, our framework successfully discriminates among existing ecosystem candidate models and discovers new ones in the absence of prior knowledge, along with simultaneous state and parameter estimation. Our framework demonstrates interdisciplinary learning and crucially provides probability distributions for each learned quantity including the learned model functions. In the real-world data experiments, we configure a one-dimensional coupled physical-biological-carbonate model to simulate the state conditions encountered by a research cruise in the Gulf of Maine region in August, 2012. Using the observed ocean acidification data, we learn and discover a salinity based forcing term for the total alkalinity (TA) equation to account for changes in TA due to advection of water masses of different salinity caused by precipitation, riverine input, and other oceanographic processes. Simultaneously, we also estimate the multidisciplinary states and an uncertain parameter. Additionally, we develop new theory and techniques to improve uncertainty quantification using the DO methodology in multidisciplinary settings, so as to accurately handle stochastic boundary conditions, complex geometries, and the advection terms, and to augment the DO subspace as and when needed to capture the effects of the truncated modes accurately. Further, we discuss mutual-information-based observation planning to determine what, when, and where to measure to best achieve our learning objectives in resource-constrained environments.

Next, motivated by the presence of inherent delays in real-world systems and the Mori-Zwanzig formulation, we develop a novel delay-differential-equations-based deep learning framework to learn time-delayed closure parameterizations for missing dynamics. We find that our neural closure models increase the long-term predictive capabilities of existing models, and require smaller networks when using non-Markovian over Markovian closures. They efficiently represent truncated modes in reduced-order models, capture effects of subgrid-scale processes, and augment the simplification of complex physical-biogeochemical models. To empower our neural closure models framework with generalizability and interpretability, we further develop neural partial delay differential equations theory that augments low-fidelity models in their original PDE forms with both Markovian and non-Markovian closure terms parameterized with neural networks (NNs). For the first time, the melding of low-fidelity model and NNs with time-delays in the continuous spatiotemporal space followed by numerical discretization automatically provides interpretability and allows for generalizability to computational grid resolution, boundary conditions, initial conditions, and problem specific parameters. We derive the adjoint equations in the continuous form, thus, allowing implementation of our new methods across differentiable and non-differentiable computational physics codes, different machine learning frameworks, and also non-uniformly-spaced spatiotemporal training data. We also show that there exists an optimal amount of past information to incorporate, and provide methodology to learn it from data during the training process. Computational advantages associated with our frameworks are analyzed and discussed. Applications of our new neural closure modeling framework are not limited to the shown fluid and ocean experiments, but can be widely extended to other fields such as control theory, robotics, pharmacokinetic-pharmacodynamics, chemistry, economics, and biological regulatory systems.

Transport of any material quantity due to background fields, i.e. advective transport, in fluid dynamical systems has been a widely studied problem. It is of crucial importance in classical fluid mechanics, geophysical flows, micro and nanofluidics, and biological flows. Even though mathematical models that thoroughly describe such transport exist, the inherent nonlinearities and the high dimensionality of complex fluid systems make it very challenging to develop the capabilities to accurately compute and characterize advective material transport. We systematically study the problems of predicting, uncovering, and learning the principal features of advective material transport in this work. The specific objectives of this thesis are to: (i) develop and apply new numerical methodologies to compute the solutions of advective transport equations with minimal errors and theoretical guarantees, (ii) propose and theoretically investigate novel criteria to detect sets of fluid parcels that remain the most coherent / incoherent throughout an extended time interval to quantify fluid mixing, and (iii) extend and develop new machine learning methods to infer and predict the transport features, given snapshot data about passive and active material transport.

The first part of this work deals with the development of the PDE-based ‘method of flow map composition’, which is a novel methodology to compute the solutions of the partial differential equation describing classical advective and advective-diffusive-reactive transport. The method of composition yields solutions almost devoid of numerical errors, and is readily parallelizable. It can compute more accurate solutions in less time than traditional numerical methods. We also complete a comprehensive theoretical analysis and analytically obtain the value of the numerical timestep that minimizes the net error. The method of flow map composition is extensively benchmarked and its applications are demonstrated in several analytical flow fields and realistic data-assimilative ocean plume simulations.

We then utilize the method of flow map composition to analyze Lagrangian material coherence in dynamic open domains. We develop new theory and schemes to efficiently predict the sets of fluid parcels that either remain the most or the least coherent over an extended amount of time. We also prove that these material sets are the ones to maximally resist advective stretching and diffusive transport. Thus, they are of significant importance in understanding the dynamics of fluid mixing and form the skeleton of material transport in unsteady fluid systems. The developed theory and numerical methods are utilized to analyze Lagrangian coherence in analytical and realistic scenarios. We emphasize realistic marine flows with multiple time-dependent inlets and outlets, and demonstrate applications in diverse dynamical regimes and several open ocean regions.

The final part of this work investigates the machine inference and prediction of the principal transport features from snapshot data about the transport of some material quantity. Our goals include machine learning the underlying advective transport features, coherent / incoherent sets, and attracting and repelling manifolds, given the snapshots of advective and advective-diffusive material fields. We also infer and predict high resolution transport features by optimally combining coarse resolution snapshot data with localized high resolution trajectory data. To achieve these goals, we use and extend recurrent neural networks, including a combination of long short-term memory networks with hypernetworks. We develop methods that leverage our knowledge of the physical system in the design and architecture of the neural network and enforce the known constraints that the results must satisfy (e.g. mass conservation) in the training loss function. This allows us to train the networks only with partial supervision, without samples of the expected output fields, and still infer and predict physically consistent quantities. The developed theory, methods, and computational software are analyzed, validated, and applied to a variety of analytical and realistic fluid flows, including high-resolution ocean transports in the Western Mediterranean Sea.

The rapidly-growing computational power and the increasing capability of uncertainty quantification, statistical inference, and machine learning have opened up new opportunities for utilizing data to assist, identify and refine physical models. In this thesis, we focus on Bayesian learning for a particular class of models: high-dimensional nonlinear dynamical systems, which have been commonly used to predict a wide range of transient phenomena including fluid flows, heat transfer, biogeochemical dynamics, and other advection-diffusion-reaction-based transport processes. Even though such models often express the differential form of fundamental laws, they commonly contain uncertainty in their initial and boundary values, parameters, forcing and even formulation. Learning such components from sparse observation data by principled Bayesian inference is very challenging due to the systems’ high-dimensionality and nonlinearity.

We systematically study the theoretical and algorithmic properties of a Bayesian learning methodology built upon previous efforts in our group to address this challenge. Our systematic study breaks down into the three hierarchical components of the Bayesian learning and we develop new numerical schemes for each. The first component is on uncertainty quantification for stochastic dynamical systems and fluid flows. We study dynamic low-rank approximations using the dynamically orthogonal (DO) equations including accuracy and computational costs, and develop new numerical schemes for re-orthonormalization, adaptive subspace augmentation, residual-driven closure, and stochastic Navier-Stokes integration. The second part is on Bayesian data assimilation, where we study the properties of and connections among the different families of nonlinear and non-Gaussian filters. We derive an ensemble square-root filter based on minimal-correction second-moment matching that works especially well under the adversity of small ensemble size, sparse observations and chaotic dynamics. We also obtain a localization technique for filtering with high-dimensional systems that can be applied to nonlinear non-Gaussian inference with both brute force Monte Carlo (MC) and reduced subspace modeling in a unified way. Furthermore, we develop a mutual-information-based adaptive sampling strategy for filtering to identify the most informative observations with respect to the state variables and/or parameters, utilizing the sub-modularity of mutual information due to the conditional independence of observation noise. The third part is on active Bayesian model learning, where we have a discrete set of candidate dynamical models and we infer the model formulation that best explains the data using principled Bayesian learning. To predict the observations that are most useful to learn the model formulation, we further extend the above adaptive sampling strategy to identify the data that are expected to be most informative with respect to both state variables and the uncertain model identity.

To investigate and showcase the effectiveness and efficiency of our theoretical and numerical advances for uncertainty quantification, Bayesian data assimilation, and active Bayesian learning with stochastic nonlinear high-dimensional dynamical systems, we apply our dynamic data-driven reduced subspace approach to several dynamical systems and compare our results against those of brute force MC and other existing methods. Specifically, we analyze our advances using several drastically different dynamical regimes modeled by the nonlinear Lorenz-96 ordinary differential equations as well as turbulent bottom gravity current dynamics modeled by the 2-D unsteady incompressible Reynolds-averaged Navier-Stokes (RANS) partial differential equations. We compare the accuracy, efficiency, and robustness of different methodologies and algorithms. With the Lorenz-96 system, we show how the performance differs under periodic, weakly chaotic, and very chaotic dynamics and under different observation layouts. With the bottom gravity current dynamics, we show how model parameters, domain geometries, initial fields, and boundary forcing formulations can be identified and how the Bayesian methodology performs when the candidate model space does not contain the true model. The results indicate that our active Bayesian learning framework can better infer the state variables and dynamical model identity with fewer observations than many alternative approaches in the literature.

The coastal ocean is a prime example of multiscale nonlinear fluid dynamics. Ocean fields in such regions are complex, with multiple spatial and temporal scales and nonstationary heterogeneous statistics. Due to the limited measurements, there are multiple sources of uncertainties, including the initial conditions, boundary conditions, forcing, parameters, and even the model parameterizations and equations themselves. To reduce uncertainties and allow long-duration measurements, the energy consumption of ocean observing platforms need to be optimized. Predicting the distributions of reachable regions, time-optimal paths, and risk-optimal paths in uncertain, strong and dynamic flows is also essential for their optimal and safe operations. Motivated by the above needs, the objectives of this thesis are to develop and apply the theory, schemes, and computational systems for: (i) Dynamically Orthogonal ocean primitive-equations with a nonlinear free-surface, in order to quantify uncertainties and predict probabilities for four-dimensional (time and 3-d in space) coastal ocean states, respecting their nonlinear governing equations and non-Gaussian statistics; (ii) Stochastic Dynamically Orthogonal level-set optimization to rigorously incorporate realistic ocean flow forecasts and plan energy-optimal paths of autonomous agents in coastal regions; (iii) Probabilistic predictions of reachability, time-optimal paths and risk-optimal paths in uncertain, strong and dynamic flows.

For the first objective, we further develop and implement our Dynamically Orthogonal (DO) numerical schemes for idealized and realistic ocean primitive equations with a nonlinear free-surface. The theoretical extensions necessary for the free-surface are completed. DO schemes are researched and DO terms, functions, and operations are implemented, focusing on: state variable choices; DO norms; DO condition for flows with a dynamic free-surface; diagnostic DO equations for pressure, barotropic velocities and density terms; non-polynomial nonlinearities; semi-implicit time-stepping schemes; and re-orthonormalization consistent with leap-frog time marching. We apply the new DO schemes, as well as their theoretical extensions and efficient serial implementation to forecast idealized-to-realistic stochastic
coastal ocean dynamics. For the realistic simulations, probabilistic predictions for the Middle Atlantic Bight region, Northwest Atlantic, and northern Indian ocean are showcased.

For the second objective, we integrate data-driven ocean modeling with our stochastic DO level-set optimization to compute and study energy-optimal paths, speeds, and headings for ocean vehicles in the Middle Atlantic Bight region. We compute the energy-optimal paths from among exact time-optimal paths. For ocean currents, we utilize a data-assimilative multiscale re-analysis, combining observations with implicit two-way nested multi-resolution primitive-equation simulations of the tidal-to-mesoscale dynamics in the region. We solve the reduced-order stochastic DO level-set partial differential equations (PDEs) to compute the joint probability of minimum arrival-time, vehicle-speed time-series, and total energy utilized. For each arrival time, we then select the vehicle-speed time-series that minimize the total energy utilization from the marginal probability of vehicle-speed and total energy. The corresponding energy-optimal path and headings be obtained through a particle backtracking equation. For the missions considered, we analyze the effects of the regional tidal currents, strong wind events, coastal jets, shelfbreak front, and other local circulations on the energyoptimal paths.

For the third objective, we develop and apply stochastic level-set PDEs that govern the stochastic time-optimal reachability fronts and paths for vehicles in uncertain, strong, and dynamic flow fields. To solve these equations efficiently, we again employ their dynamically orthogonal reduced-order projections. We develop the theory and schemes for risk-optimal planning by combining decision theory with our stochastic time-optimal planning equations. The risk-optimal planning proceeds in three steps: (i) obtain predictions of the probability distribution of environmental flows, (ii) obtain predictions of the distribution of exact timeoptimal paths for the forecast flow distribution, and (iii) compute and minimize the risk of following these uncertain time-optimal paths. We utilize the new equations to complete stochastic reachability, time-optimal and risk-optimal path planning in varied stochastic quasi-geostrophic flows. The effects of the flow uncertainty on the reachability fronts and time-optimal paths is explained. The risks of following each exact time-optimal path is evaluated and risk-optimal paths are computed for different risk tolerance measures. Key properties of the risk-optimal planning are finally discussed.

Theoretically, the present methodologies are PDE-based and compute stochastic ocean fields, and optimal path predictions without heuristics. Computationally, they are several orders of magnitude faster than direct Monte Carlo.

Such technologies have several commercial and societal applications. Specifically, the probabilistic ocean predictions can be input to a technical decision aide for a sustainable fisheries co-management program in India, which has the potential to provide environment friendly livelihoods to millions of marginal fishermen. The risk-optimal path planning equations can be employed in real-time for efficient ship routing to reduce greenhouse gas emissions and save operational costs.

When humans or robots operate in complex dynamic environments, the planning of paths and the collection of observations are basic, indispensable problems. In the oceanic and atmospheric environments, the concurrent use of multiple mobile sensing platforms in unmanned missions is growing very rapidly. Opportunities for a paradigm shift in the science of autonomy involve the development of fundamental theories to optimally collect information, learn, collaborate and make decisions under uncertainty while persistently adapting to and utilizing the dynamic environment. To address such pressing needs, this thesis derives governing equations and develops rigorous methodologies for optimal path planning and optimal sampling using collaborative swarms of autonomous mobile platforms. The application focus is the coastal ocean where currents can be much larger than platform speeds, but the fundamental results also apply to other dynamic environments. We first undertake a theoretical synthesis of minimum-time control of vehicles operating in general dynamic flows. Using various ideas rooted in non-smooth calculus, we prove that an unsteady Hamilton-Jacobi equation governs the forward reachable sets in any type of Lipschitz-continuous flow. Next, we show that with a suitable modification to the Hamiltonian, the results can be rigorously generalized to perform time-optimal path planning with anisotropic motion constraints and with moving obstacles and unsafe ‘forbidden’ regions. We then derive a level-set methodology for distance-based coordination of swarms of vehicles operating in minimum time within strong and dynamic ocean currents. The results are illustrated for varied fluid and ocean flow simulations. Finally, the new path planning system is applied to swarms of vehicles operating in the complex geometry of the Philippine Archipelago, utilizing realistic multi-scale current predictions from a data-assimilative ocean modeling system. In the second part of the thesis, we derive a theory for adaptive sampling that exploits the governing nonlinear dynamics of the system and captures the non-Gaussian structure of the random state fields. Optimal observation locations are determined by maximizing the mutual information between the candidate observations and the variables of interest. We develop a novel Bayesian smoother for high-dimensional continuous stochastic fields governed by general nonlinear dynamics. This smoother combines the adaptive reduced-order Dynamically-Orthogonal equations with Gaussian Mixture Models, extending linearized Gaussian backward pass updates to a nonlinear, non-Gaussian setting. The Bayesian information transfer, both forward and backward in time, is efficiently carried out in the evolving dominant stochastic subspace. Building on the foundations of the smoother, we then derive an efficient technique to quantify the spatially and temporally varying mutual information field in general nonlinear dynamical systems. The globally optimal sequence of future sampling locations is rigorously determined by a novel dynamic programming approach that combines this computation of mutual information fields with the predictions of the forward reachable set. All the results are exemplified and their performance is quantitatively assessed using a variety of simulated fluid and ocean flows. The above novel theories and schemes are integrated so as to provide real-time computational intelligence for collaborative swarms of autonomous sensing vehicles. The integrated system guides groups of vehicles along predicted optimal trajectories and continuously improves field estimates as the observations predicted to be most informative are collected and assimilated. The optimal sampling locations and optimal trajectories are continuously forecast, all in an autonomous and coordinated fashion.

Accurate modeling of physical and biogeochemical dynamics in coastal ocean regions is required for multiple scientific and societal applications, covering a wide range of time and space scales. However, in light of the strong nonlinearities observed in coastal regions and in biological processes, such modeling is challenging. An important subject that has been largely overlooked is the numerical requirements for regional ocean simulation studies. Major objectives of this thesis are to address such computational questions for non-hydrostatic multiscale flows and for biogeochemical interactions, and to derive and develop numerical schemes that meet these requirements, utilizing the latest advances in computational fluid dynamics. We are interested in studying nonlinear, transient, and multiscale ocean dynamics over complex geometries with steep bathymetry and intricate coastlines, from sub-mesoscales to basin-scales. These dynamical interests, when combined with our requirements for accurate, efficient and flexible ocean modeling, led us to develop new variable resolution, higher-order and non-hydrostatic ocean modeling schemes. Specifically, we derived, developed and applied new numerical schemes based on the novel hybrid discontinuous Galerkin (HDG) method in combination with projection methods. The new numerical schemes are first derived for the Navier-Stokes equations. To ensure mass conservation, we define numerical fluxes that are consistent with the discrete divergence equation. To improve stability and accuracy, we derive a consistent HDG stability parameter for the pressure-correction equation. We also apply a new boundary condition for the pressure-corrector, and show the form and origin of the projection method’s time-splitting error for a case with implicit diffusion and explicit advection. Our scheme is implemented for arbitrary, mixed-element unstructured grids using a novel quadrature-free integration method for a nodal basis, which is consistent with the HDG method. To prevent numerical oscillations, we design a selective high-order nodal limiter. We demonstrate the correctness of our new schemes using a tracer advection benchmark, a manufactured solution for the steady diffusion and stokes equations, and the 2D lock-exchange problem. These numerical schemes are then extended for non-hydrostatic, free-surface, variable-density regional ocean dynamics. The time-splitting procedure using projection methods is derived for non-hydrostatic or hydrostatic, and nonlinear free-surface or rigid-lid, versions of the model. We also derive consistent HDG stability parameters for the free-surface and non-hydrostatic pressure-corrector equations to ensure stability and accuracy. New boundary conditions for the free-surface-corrector and pressure-corrector are also introduced. We prove that these conditions lead to consistent boundary conditions for the free-surface and pressure proper. To ensure discrete mass conservation with a moving free-surface, we use an arbitrary LagrangianEulerian (ALE) moving mesh algorithm. These schemes are again verified, this time using a tidal flow problem with analytical solutions and a 3D lock-exchange benchmark. We apply our new numerical schemes to evaluate the numerical requirements of the coupled biological-physical dynamics. We find that higher-order schemes are more accurate at the same efficiency compared to lower-order (e.g. second-order) accurate schemes when modeling a biological patch. Due to decreased numerical dissipation, the higher-order schemes are capable of modeling biological patchiness over a sustained duration, while the lower-order schemes can lose significant biomass after a few non-dimensional times and can thus solve erroneous nonlinear dynamics. Finally, inspired by Stellwagen Bank in Massachusetts Bay, we study the effect of non-hydrostatic physics on biological productivity and phytoplankton fields for tidally-driven flows over an idealized bank. We find that the non-hydrostatic pressure and flows are important for biological dynamics, especially when flows are supercritical. That is, when the slope of the topography is larger than the slope of internal wave rays at the tidal frequency. The non-hydrostatic effects increase with increasing nonlinearity, both when the internal Froude number and criticality parameter increase. Even in cases where the instantaneous biological productivity is not largely modified, we find that the total biomass, spatial variability and patchiness of phytoplankton can be significantly altered by non-hydrostatic processes. Our ultimate dynamics motivation is to allow quantitative simulation studies of fundamental nonlinear biological-physical dynamics in coastal regions with complex bathymetric features such as straits, sills, ridges and shelfbreaks. This thesis develops the necessary numerical schemes that meet the stringent accuracy requirements for these types of flows and dynamics.

In the past decades an increasing number of problems in continuum theory have been treated using stochastic dynamical theories. This is because dynamical systems governing real processes always contain some elements characterized by uncertainty or stochasticity. Uncertainties may arise in the system parameters, the boundary and initial conditions, and also in the external forcing processes. Also, many problems are treated through the stochastic framework due to the incomplete or partial understanding of the governing physical laws. In all of the above cases the existence of random perturbations, combined with the com- plex dynamical mechanisms of the system often leads to their rapid growth which causes distribution of energy to a broadband spectrum of scales both in space and time, making the system state particularly complex. Such problems are mainly described by Stochastic Partial Differential Equations and they arise in a number of areas including fluid mechanics, elasticity, and wave theory, describing phenomena such as turbulence, random vibrations, flow through porous media, and wave propagation through random media. This is but a partial listing of applications and it is clear that almost any phenomenon described by a field equation has an important subclass of problems that may profitably be treated from a stochastic point of view.

In this work, we develop a new methodology for the representation and evolution of the complete probabilistic response of infinite-dimensional, random, dynamical systems. More specifically, we derive an exact, closed set of evolution equations for general nonlinear continuous stochastic fields described by a Stochastic Partial Differential Equation. The derivation is based on a novel condition, the Dynamical Orthogonality (DO), on the representation of the solution. This condition is the “key” to overcome the redundancy issues of the full representation used while it does not restrict its generic features. Based on the DO condition we derive a system of field equations consisting of a Partial Differential Equation (PDE) for the mean field, a family of PDEs for the orthonormal basis that describe the stochastic subspace where uncertainty “lives” as well as a system of Stochastic Differential Equations that defines how the uncertainty evolves in the time varying stochastic subspace. If additional restrictions are assumed on the form of the representation, we recover both the Proper-Orthogonal-Decomposition (POD) equations and the generalized Polynomial- Chaos (PC) equations; thus the new methodology generalizes these two approaches. For the efficient treatment of the strongly transient character on the systems described above we derive adaptive criteria for the variation of the stochastic dimensionality that characterizes the system response. Those criteria follow directly from the dynamical equations describing the system.

We illustrate and validate this novel technique by solving the 2D stochastic Navier-Stokes equations in various geometries and compare with direct Monte Carlo simulations. We also apply the derived framework for the study of the statistical responses of an idealized “double gyre” model, which has elements of ocean, atmospheric and climate instability behaviors.

Finally, we use our new stochastic description for flow fields to study the motion of inertial particles in flows with uncertainties. Inertial or finite-size particles in fluid flows are commonly encountered in nature (e.g., contaminant dispersion in the ocean and atmosphere) as well as in technological applications (e.g., chemical systems involving particulate reactant mixing). As it has been observed both numerically and experimentally, their dynamics can differ markedly from infinitesimal particle dynamics. Here we use recent results from stochastic singular perturbation theory in combination with the DO representation of the random flow, in order to derive a reduced order inertial equation that will describe efficiently the stochastic dynamics of inertial particles in arbitrary random flows.

In shallow water, a large part of underwater acoustic prediction uncertainties are in- duced by sub-meso-to-small scale oceanographic variabilities. Conventional oceano- graphic measurements for capturing such ocean-acoustic environmental variabilities face the classical conflict between resolution and coverage. The Adaptive Rapid En- vironmental Assessment (AREA) project was proposed to resolve this conflict by optimizing the location of in-situ measurements in an adaptive manner. In this thesis, ideas, concepts and performance limits in AREA are clarified. Both an engineering and a mathematical model for AREA are developed. A modularized AREA simulator was developed and implemented in C++. Philosophies in AREA are discussed. Presumptions about the ocean are made to bridge the gap between the viewpoint in the oceanography community, where the ocean environment is consid- ered to be a deterministic but very complicated system, and that of the underwater acoustic community, where the ocean environment is treated as a random system. At present, how to optimally locate the in-situ measurements made by a single AUV carrying a CTD (conductivity, temperature and depth) sensor is considered in AREA. In this thesis, the AUV path planning is modeled as a Shortest Path problem. However, due to the sound velocity correlation effect, the size of this problem can be very large. A method is developed to simplify the graph for a fast solution. As a significant step, a linear approximation for acoustic Transmission Loss (TL) is investigated numerically and analytically. In addition to following a predetermined path, an AUV can also adaptively gener- ate its path on-board. This adaptive on-board AUV routing problem is modeled using Dynamic Programming (DP) in this thesis. A method based on an optimized prede- termined path is developed to reduce the size of the DP problem and approximately yet efficiently solve it using Pattern Recognition. As a special case, a thermocline- oriented AUV yoyo control and control parameter optimization methods for AREA are also developed. 2 Finally, some AUV control algorithms for capturing fronts are developed. A frame- work for real-time TL forecasts is developed. This is the first time that TL forecasts have been linked with ocean forecasts in real-time. All of the above ideas and methods developed were tested in two experiments, FAF05 in the northern Tyrrhenian Sea in 2005 and MB06 in Monterey Bay, CA in 2006. The latter MB06 sea exercise was a major field experiment sponsored by the Office of Naval Research and the thesis compiles significant findings from this effort.