In the coastal waters of Taiwan, freshwater discharge from rivers can be an important source of uncertainty in regional ocean simulations. This effect becomes especially acute during extreme storm events, such as typhoons. In particular, record-breaking discharge caused by Typhoon Morakot (August 6-10, 2009) was observed to significantly affect near-shore temperature and salinity during the Intensive Observation Period-09 (IOP09) of the Quantifying, Predicting and Exploiting Uncertainty (QPE) research initiative. In this report, a river discharge model is developed to account for the sudden large influx of freshwater during and after the typhoon. The discharge model is then evaluated by comparison with the discharge time series for the Zhuoshu and Gaoping Rivers and by its utilization as forcing in ocean simulations. The parameters of the discharge and river forcing models and their effects on ocean simulations are discussed. The reanalysis ocean simulations with river forcing are shown to capture several of the independently observed features in the evolution of the coastal salinity field as well as the magnitude of the freshening of the ocean caused by runoff from Typhoon Morakot.

In this thesis, we explore the use of stochastic Navier-Stokes equations through the Dynamically Orthogonal (DO) methodology developed at MIT in the Multidisciplinary Simulation, Estimation, and Assimilation Systems Group. Specifically, we examine the effects of the Reynolds number on stochastic fluid flows behind a square cylinder and evaluate computational schemes to do so. We review existing literature, examine our simulation results and validate the numerical solution. The thesis uses a novel open boundary condition formulation for DO stochastic Navier-Stokes equations, which allows the modeling of a wide range of random inlet boundary conditions with a single DO simulation of low stochastic dimensions, reducing computational costs by orders of magnitude. We first test the numerical convergence and validating the numerics. We then study the sensitivity of the results to several parameters, focusing for the dynamics on the sensitivity to the Reynolds number. For the method, we focus on the sensitivity to the: resolution of in the stochastic subspace, resolution in the physical space and number of open boundary conditions DO modes. Finally, we evaluate and study how key dynamical characteristics of the flow such as the recirculation length and the vortex shedding period vary with the Reynolds number.

The quantification of uncertainties is critical when systems are nonlinear and have uncertain terms in their governing equations or are constrained by limited knowledge of initial and boundary conditions. Such situations are common in multi-scale, intermittent and non-homogeneous fluid and ocean flows, and other non-linear dynamical systems. The Dynamically Orthogonal (DO) field equations provide an efficient time- dependent adaptive methodology to predict the probability density functions of such dynamics. The present work derives efficient computational schemes for the DO methodology applied to unsteady stochastic Navier- Stokes and Boussinesq equations, and illustrates and studies the numerical aspects of these schemes. Semi- implicit projection methods are developed for the mean and for the orthonormal modes that define a basis for the evolving DO subspace, and time-marching schemes of first to fourth order are used for the stochastic coefficients. Conservative second-order nite-volumes are employed in physical space with new advection schemes based on Total Variation Diminishing methods. Other results specific to the DO equations include: (i) the definition of pseudo-stochastic pressures to obtain a number of pressure equations that is linear in the subspace size instead of quadratic; (ii) symmetric advection schemes for the stochastic velocities; (iii) the use of generalized inversion to deal with singular subspace covariances or deterministic modes; and (iv) schemes to maintain orthonormal modes at the numerical level. While (i) and (ii) are specific to fluid flows, (iii) and (iv) are important for any system of equations discretized using the DO methodology. To verify the correctness of our implementation and study the properties of our schemes and their variations, a set of stochastic flow benchmarks are defined including asymmetric Dirac and symmetric lock-exchange flows, lid-driven cavity flows, and flows past objects in a confined channel. Different Reynolds number and Grashof number regimes are employed to illustrate robustness. Optimal convergence under both time and space refinements is shown as well as the convergence of the probability density functions with the number of stochastic realizations.

Data assimilation, as presented in this thesis, is the statistical merging of sparse observational data with computational models so as to optimally improve the probabilistic description of the field of interest, thereby reducing uncertainties. The centerpiece of this thesis is the introduction of a novel such scheme that overcomes prior shortcomings observed within the community. Adopting techniques prevalent in Machine Learning and Pattern Recognition, and building on the foundations of classical assimilation schemes, we introduce the GMM-DO filter: Data Assimilation with Gaussian mixture models using the Dynamically Orthogonal field equations.

We combine the use of Gaussian mixture models, the EM algorithm and the Bayesian Information Criterion to accurately approximate distributions based on Monte Carlo data in a framework that allows for efficient Bayesian inference. We give detailed descriptions of each of these techniques, supporting their application by recent literature. One novelty of the GMM-DO filter lies in coupling these concepts with an efficient representation of the evolving probabilistic description of the uncertain dynamical field: the Dynamically Orthogonal field equations. By limiting our attention to a dominant evolving stochastic subspace of the total state space, we bridge an important gap previously identified in the literature caused by the dimensionality of the state space.

We successfully apply the GMM-DO filter to two test cases: (1) the Double Well Diffusion Experiment and (2) the Sudden Expansion fluid flow. With the former, we prove the validity of utilizing Gaussian mixture models, the EM algorithm and the Bayesian Information Criterion in a dynamical systems setting. With the application of the GMM-DO filter to the two-dimensional Sudden Expansion fluid flow, we further show its applicability to realistic test cases of non-trivial dimensionality. The GMM-DO filter is shown to consistently capture and retain the far-from-Gaussian statistics that arise, both prior and posterior to the assimilation of data, resulting in its superior performance over contemporary filters. We present the GMM-DO filter as an efficient, data-driven assimilation scheme, focused on a dominant evolving stochastic subspace of the total state space, that respects nonlinear dynamics and captures non-Gaussian statistics, obviating the use of heuristic arguments.

From naval operations to ocean science missions, the importance of autonomous vehicles is increasing with the advances in underwater robotics technology. Due to the dynamic and intermittent underwater environment and the physical limitations of autonomous underwater vehicles, feasible and optimal path planning is crucial for autonomous underwater operations. The objective of this thesis is to develop and demonstrate an efficient underwater path planning algorithm based on the level set method. Specifically, the goal is to compute the paths of autonomous vehicles which minimize travel time in the presence of ocean currents. The approach is to either utilize or avoid any type of ocean flows, while allowing for currents that are much larger than the nominal vehicle speed and for three-dimensional currents which vary with time. Existing path planning methods for the fields of ocean science and robotics are first reviewed, and the advantages and disadvantages of each are discussed. The underpinnings of the level set and fast marching methods are then reviewed, including their new extension and application to underwater path planning. Finally, a new feasible and optimal time-dependent underwater path planning algorithm is derived and presented. In order to demonstrate the capabilities of the algorithm, a set of idealized test-cases of increasing complexity are first presented and discussed. A real three-dimensional path planning example, involving strong current conditions, is also illustrated. This example utilizes four-dimensional ocean flows from a realistic ocean prediction system which simulate the ocean response to the passage of a tropical storm in the Middle Atlantic Bight region.

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.

Understanding the complex dynamics of coastal upwelling is essential for coastal ocean dynamics, phytoplankton blooms, and pollution transport. Atmospheric-driven coastal upwelling often occurs when strong alongshore winds and the Coriolis force combine to displace warmer surface waters offshore, leading to upward motions of deeper cooler, nutrient-dense waters to replace these surface waters. Using the models of the MIT Multidisciplinary Simulation, Estimation, and Assimilation System (MSEAS) group, we conduct a large set of simulation sensitivity studies to determine which variables are dominant controls for upwelling events in the Monterey Bay region. Our motivations include determining the dominant atmospheric fluxes and the causes of high-frequency fluctuations found in ocean thermal balances. We focus on the first upwelling event from August 1-5, 2006 in Monterey Bay that occurred during the Monterey Bay 06 (MB06) at-sea experiment, for which MSEAS data-assimilative baseline simulations already existed.

Using the thermal energy (temperature), salinity and momentum (velocity) conservation equations, full ocean fields in the region as well as both control volume (flux) balances and local differential term-by-term balances for the upwelling event events were computed. The studies of ocean fields concentrate on specific depths: surface-0m, thermocline-30m and undercurrent-150m. Effects of differing atmospheric forcing contributions (wind stress, surface heating/cooling, and evaporation-precipitation) on these full fields and on the volume and term-by-term balances are analyzed. Tidal effects are quantified utilizing pairs of simulations in which tides are either included or not. Effects of data assimilation are also examined.

We find that the wind stress forcing is the most important dynamical parameter in explaining the extent and shape of the upwelling event. This is verified using our large set of sensitivity studies and examining the heat flux balances. The assimilation of data has also an impact because this first upwelling event occurs during the initialization. Tidal forcing and, to a lesser extent, the daily atmospheric and data assimilation cycles explain the higher frequency fluctuations found in the volume averaged time rate of change of thermal energy.

A new generation of efficient parallel, multi-scale, and interdisciplinary ocean models is required for better understanding and accurate predictions. The purpose of this thesis is to quantitatively identify promising numerical methods that are suitable to such predictions. In order to fulfill this purpose, current efforts towards creating new ocean models are reviewed, an understanding of the most promising methods used by other researchers is developed, the most promising existing methods are studied and applied to idealized cases, new methods are incubated and evaluated by solving test problems, and important numerical issues related to efficiency are examined.

The results of other research groups towards developing the second generation of ocean models are first reviewed. Next, the Discontinuous Galerkin (DG) method for solving advection-diffusion problems is described, including a discussion on schemes for solving higher order derivatives. The discrete formulation for advection-diffusion problems is detailed and implementation issues are discussed. The Hybrid Discontinuous Galerkin (HDG) Finite Element Method (FEM) is identified as a promising new numerical scheme for ocean simulations. For the first time, a DG FEM scheme is used to solve ocean biogeochemical advection-diffusion-reaction equations on a two- dimensional idealized domain, and p-adaptivity across constituents is examined. Each aspect of the numerical solution is examined separately, and p-adaptive strategies are explored. Finally, numerous solver-preconditioner combinations are benchmarked to identify an efficient solution method for inverting matrices, which is necessary for implicit time integration schemes. From our quantitative incubation of numerical schemes, a number of recommendations on the tools necessary to solve dynamical equations for multiscale ocean predictions are provided.

A fundamental requirement in realistic computational geophysical fluid dynamics is the optimal estimation of gridded fields and of spatial-temporal scales directly from the spatially irregular and multivariate data sets that are collected by varied instruments and sampling schemes. In this work, we derive and utilize new schemes for the mapping and dynamical inference of ocean fields in complex multiply-connected domains, study the computational properties of our new mapping schemes, and derive and investigate new schemes for adaptive estimation of spatial and temporal scales. Objective Analysis (OA) is the statistical estimation of fields using the Bayesian- based Gauss-Markov theorem, i.e. the update step of the Kalman Filter. The existing multi-scale OA approach of the Multidisciplinary Simulation, Estimation and Assimilation System consists of the successive utilization of Kalman update steps, one for each scale and for each correlation across scales. In the present work, the approach is extended to field mapping in complex, multiply-connected, coastal regions and archipelagos. A reasonably accurate correlation function often requires an estimate of the distance between data and model points, without going across complex land- forms. New methods for OA based on estimating the length of optimal shortest sea paths using the Level Set Method (LSM) and Fast Marching Method (FMM) are derived, implemented and utilized in general idealized and realistic ocean cases. Our new methodologies could improve widely-used gridded databases such as the climatological gridded fields of the World Ocean Atlas (WOA) since these oceanic maps were computed without accounting for coastline constraints. A new FMM-based methodology for the estimation of absolute velocity under geostrophic balance in complicated domains is also outlined. Our new schemes are compared with other approaches, including the use of stochastically forced differential equations (SDE). We find that our FMM-based scheme for complex, multiply-connected, coastal regions is more efficient and accurate than the SDE approach. We also show that the field maps obtained using our FMM-based scheme do not require postprocessing (smoothing) of fields. The computational properties of the new mapping schemes are studied in detail. We find that higher-order schemes improve the accuracy of distance estimates. We also show that the covariance matrices we estimate are not necessarily positive definite because the Weiner Khinchin and Bochner relationships for positive definiteness are only valid for convex simply-connected domains. Several approaches to overcome this issue are discussed and qualitatively evaluated. The solutions we propose include introducing a small process noise or reducing the covariance matrix based on the dominant singular value decomposition. We have also developed and utilized novel methodologies for the adaptive estimation of spatial-temporal scales from irregularly spaced ocean data. The three novel methodologies are based on the use of structure functions, short term Fourier transform and second generation wavelets. To our knowledge, this is the first time that adaptive methodologies for the spatial-temporal scale estimation are proposed. The ultimate goal of all these methods would be to create maps of spatial and temporal scales that evolve as new ocean data are fed to the scheme. This would potentially be a significant advance to the ocean community for better understanding and sampling of ocean processes.

In this thesis, we conduct research toward understanding coupled physics-biology processes in ocean straits. Our focus is on new analytical studies and higher-order simulations of idealized dynamics that are relevant to generic biological processes. The details of coupled physics-biology models are reviewed and an in-depth global equilibrium and local stability analysis of a Nutrient-Phytoplankton-Zooplankton (NPZ) model is performed. This analysis includes parameter studies and methods to evaluate parameter sensitivity, especially in the case where some system parameters are unknown. As an initial step toward investigating the interaction between physics and biology in ocean straits, we develop and verify a new coupled physics-biology model for two-dimensional idealized physical processes including tides and apply it to the San Bernardino Strait in the Philippine Archipelago. This two-dimensional numerical model is created on a structured grid using operator splitting and masking. This model is able to accurately represent biology for various physical flows, including advection-dominated flows over discontinuities, by using the Weighted Essentially Non-Oscillatory (WENO) scheme. The numerical model is verified against a Discontinuous-Galerkin (DG) numerical scheme on an unstructured grid. Several simulations of tidal flow are completed using bathymetry and flow magnitudes com- parable to those found in the San Bernardino Strait with different sets of parameters, tidal periods, and levels of diffusion. Results are discussed and compared to those of a three-dimensional modeling system. New results include: new methods for analyzing stability, the robust two-dimensional model designed to best represent advection-dominant flows with minimal numerical diffusion and computational time, and a novel technique to initialize three-dimensional biology fields using satellite data. Additionally, application of the two-dimensional model with tidal forcing to the San Bernardino Strait reveals that flow frequencies have strong influence on biology, as very fast oscillations act to stabilize biology in the water column, while slower frequencies provide sufficient transport for increased biological activity.

In this work, we explore the different methods for parameter estimation in straightforward diffusion problems and develop ideas and distributed computational schemes for the automated evaluation of physical and numerical parameters of ocean models. This is one step of “adaptive modeling”. Adaptive modeling consists of the automated adjustment of self-evaluating models in order to best represent an observed system. In the case of dynamic parameterizations, self-modifying schemes are used to learn the correct model for a particular regime as the physics change and evolve in time.

The parameter estimation methods are tested and evaluated on one-dimensional tracer diffusion problems. Existing state estimation methods and new filters, such as the unscented transform Kalman filter, are utilized in carrying out parameter estimation. These include the popular Extended Kalman Filter (EKF), the Ensemble Kalman Filter (EnKF) and other ensemble methods such as Error Subspace Statistical Estimation (ESSE) and Ensemble Adjustment Kalman Filter (EAKF), and the Unscented Kalman Filter (UKF). Among the aforementioned recursive state estimation methods, the so-called “adjoint method” is also applied to this simple study.

Finally, real data is examined for the applicability of such schemes in real-time fore- casting using the MIT Multidisciplinary Simulation, Estimation, and Assimilation System (MSEAS). The MSEAS model currently contains the free surface hydrostatic primitive equation model from the Harvard Ocean Prediction System (HOPS), a barotropic tidal prediction scheme, and an objective analysis scheme, among other models and developing routines. The experiment chosen for this study is one which involved the Monterey Bay region off the coast of California in 2006 (MB06). Accurate vertical mixing parameterizations are essential in this well known upwelling region of the Pacific. In this realistic case, parallel computing will be utilized by scripting code runs in C-shell. The performance of the simulations with different parameters is evaluated quantitatively using Pattern Correlation Coefficient, Root Mean Squared error, and bias error. Comparisons quantitatively determined the most adequate model setup.