BZPEERPapers, people, evidence and research paths
Menu

Discover

SearchFind papers, people, topics, institutions and opportunities from one place.GraphSee how papers, authors, topics and evidence connect.

Research tools

HomeContinue from the papers, people and topics that matter now.CollectionsTurn saved papers into reading lists, evidence maps and shareable context.ResearchReconstruct the path that led you from search to paper to question.CollaborationMove from a paper, topic or expert request into a focused thread.PublishCreate full-text research articles with metadata, formulas and versions.

Network

ExpertsFind specialists whose work and reviews explain why they can help.CommunitiesJoin research groups organized around papers, topics and shared evidence.

Opportunities

ListingsBrowse roles, calls and collaborations connected to real research context.

Account

Log inReturn to your saved papers, graph views and active threads.Create accountStart saving papers, building a research trail and joining teams.
Log inCreate account

Source author record

Adrian Sandu

Adrian Sandu appears in the imported research catalog. Authorship, coauthor and topic links are available while profile ownership is still unclaimed.

ResearcherUnclaimed source record
Numerical Analysismath.NAComputational Engineering, Finance, and ScienceComputationmath.OCmath.STStatistics TheoryMachine LearningApplicationsSystems and ControlMethodologynlin.CD

Catalog footprint

What is connected

49works
12topics
4close collaborators

Actions

Connect this record

Log in to claim
Open graphBrowse works

Research graph

See the researcher in context

Open full explorer

Inspect adjacent papers, topics, institutions and collaborators without losing the researcher page.

Building this map preview

BZPEER is loading the nearby papers, people, topics and institutions for this page.

Published work

49 published item(s)

preprint2022arXiv

A Class of Multirate Infinitesimal GARK Methods

Differential equations arising in many practical applications are characterized by multiple time scales. Multirate time integration seeks to solve them efficiently by discretizing each scale with a different, appropriate time step, while ensuring the overall accuracy and stability of the numerical solution. In a seminal paper Knoth and Wolke (APNUM, 1998) proposed a hybrid solution approach: discretize the slow component with an explicit Runge-Kutta method, and advance the fast component via a modified fast differential equation. The idea led to the development of multirate infinitesimal step (MIS) methods by Wensch et al. (BIT, 2009.)Günther and Sandu (BIT, 2016) explained MIS schemes as a particular case of multirate General-structure Additive Runge-Kutta (MR-GARK) methods. The hybrid approach offers extreme flexibility in the choice of the numerical solution process for the fast component. This work constructs a family of multirate infinitesimal GARK schemes (MRI-GARK) that extends the hybrid dynamics approachin multiple ways. Order conditions theory and stability analyses are developed, and practical explicit and implicit methods of up to order four are constructed. Numerical results confirm the theoretical findings. We expect the new MRI-GARK family to be most useful for systems of equations with widely disparate time scales, where the fast process is dispersive, and where the influence of the fast component on the slow dynamics is weak.

0 citations0 reviewsNo code linkedOpen access
preprint2022arXiv

A Meta-learning Formulation of the Autoencoder Problem for Non-linear Dimensionality Reduction

A rapidly growing area of research is the use of machine learning approaches such as autoencoders for dimensionality reduction of data and models in scientific applications. We show that the canonical formulation of autoencoders suffers from several deficiencies that can hinder their performance. Using a meta-learning approach, we reformulate the autoencoder problem as a bi-level optimization procedure that explicitly solves the dimensionality reduction task. We prove that the new formulation corrects the identified deficiencies with canonical autoencoders, provide a practical way to solve it, and showcase the strength of this formulation with a simple numerical illustration.

0 citations0 reviewsNo code linkedOpen access
preprint2022arXiv

A Stochastic Covariance Shrinkage Approach in Ensemble Transform Kalman Filtering

The Ensemble Kalman Filters (EnKF) employ a Monte-Carlo approach to represent covariance information, and are affected by sampling errors in operational settings where the number of model realizations is much smaller than the model state dimension. To alleviate the effects of these errors EnKF relies on model-specific heuristics such as covariance localization, which takes advantage of the spatial locality of correlations among the model variables. This work proposes an approach to alleviate sampling errors that utilizes a locally averaged-in-time dynamics of the model, described in terms of a climatological covariance of the dynamical system. We use this covariance as the target matrix in covariance shrinkage methods, and develop a stochastic covariance shrinkage approach where synthetic ensemble members are drawn to enrich both the ensemble subspace and the ensemble transformation. We additionally provide for a way in which this methodology can be localized similar to the state-of-the-art LETKF method, and that for a certain model setup, our methodology significantly outperforms it.

0 citations0 reviewsNo code linkedOpen access
preprint2022arXiv

Eliminating Order Reduction on Linear, Time-Dependent ODEs with GARK Methods

When applied to stiff, linear differential equations with time-dependent forcing, Runge-Kutta methods can exhibit convergence rates lower than predicted by the classical order condition theory. Commonly, this order reduction phenomenon is addressed by using an expensive, fully implicit Runge-Kutta method with high stage order or a specialized scheme satisfying additional order conditions. This work develops a flexible approach of augmenting an arbitrary Runge-Kutta method with a fully implicit method used to treat the forcing such as to maintain the classical order of the base scheme. Our methods and analyses are based on the general-structure additive Runge-Kutta framework. Numerical experiments using diagonally implicit, fully implicit, and even explicit Runge-Kutta methods confirm that the new approach eliminates order reduction for the class of problems under consideration, and the base methods achieve their theoretical orders of convergence.

0 citations0 reviewsNo code linkedOpen access
preprint2022arXiv

Physics-informed neural networks for PDE-constrained optimization and control

A fundamental problem in science and engineering is designing optimal control policies that steer a given system towards a desired outcome. This work proposes Control Physics-Informed Neural Networks (Control PINNs) that simultaneously solve for a given system state, and for the optimal control signal, in a one-stage framework that conforms to the underlying physical laws. Prior approaches use a two-stage framework that first models and then controls a system in sequential order. In contrast, a Control PINN incorporates the required optimality conditions in its architecture and in its loss function. The success of Control PINNs is demonstrated by solving the following open-loop optimal control problems: (i) an analytical problem, (ii) a one-dimensional heat equation, and (iii) a two-dimensional predator-prey problem.

0 citations0 reviewsNo code linkedOpen access
preprint2021arXiv

A Stochastic Covariance Shrinkage Approach to Particle Rejuvenation in the Ensemble Transform Particle Filter

Rejuvenation in particle filters is necessary to prevent the collapse of the weights when the number of particles is insufficient to sample the high probability regions of the state space. Rejuvenation is often implemented in a heuristic manner by the addition of stochastic samples that widen the support of the ensemble. This work aims at improving canonical rejuvenation methodology by the introduction of additional prior information obtained from climatological samples; the dynamical particles used for importance sampling are augmented with samples obtained from stochastic covariance shrinkage. The ensemble transport particle filter, and its second order variant, are extended with the proposed rejuvenation approach. Numerical experiments show that modified filters significantly improve the analyses for low dynamical ensemble sizes.

0 citations0 reviewsNo code linkedOpen access
preprint2021arXiv

Multirate Linearly-Implicit GARK Schemes

Many complex applications require the solution of initial-value problems where some components change fast, while others vary slowly. Multirate schemes apply different step sizes to resolve different components of the system, according to their dynamics, in order to achieve increased computational efficiency. The stiff components of the system, fast or slow, are best discretized with implicit base methods in order to ensure numerical stability. To this end, linearly implicit methods are particularly attractive as they solve only linear systems of equations at each step. This paper develops the Multirate GARK-ROS/ROW (MR-GARK-ROS/ROW) framework for linearly-implicit multirate time integration. The order conditions theory considers both exact and approximative Jacobians. The effectiveness of implicit multirate methods depends on the coupling between the slow and fast computations; an array of efficient coupling strategies and the resulting numerical schemes are analyzed. Multirate infinitesimal step linearly-implicit methods, that allow arbitrarily small micro-steps and offer extreme computational flexibility, are constructed. The new unifying framework includes existing multirate Rosenbrock(-W) methods as particular cases, and opens the possibility to develop new classes of highly effective linearly implicit multirate integrators.

0 citations0 reviewsNo code linkedOpen access
preprint2020arXiv

A Multifidelity Ensemble Kalman Filter with Reduced Order Control Variates

This work develops a new multifidelity ensemble Kalman filter (MFEnKF) algorithm based on linear control variate framework. The approach allows for rigorous multifidelity extensions of the EnKF, where the uncertainty in coarser fidelities in the hierarchy of models represent control variates for the uncertainty in finer fidelities. Small ensembles of high fidelity model runs are complemented by larger ensembles of cheaper, lower fidelity runs, to obtain much improved analyses at only small additional computational costs. We investigate the use of reduced order models as coarse fidelity control variates in the MFEnKF, and provide analyses to quantify the improvements over the traditional ensemble Kalman filters. We apply these ideas to perform data assimilation with a quasi-geostrophic test problem, using direct numerical simulation and a corresponding POD-Galerkin reduced order model. Numerical results show that the two-fidelity MFEnKF provides better analyses than existing EnKF algorithms at comparable or reduced computational costs.

0 citations0 reviewsNo code linkedOpen access
preprint2020arXiv

An Explicit Probabilistic Derivation of Inflation in a Scalar Ensemble Kalman Filter for Finite Step, Finite Ensemble Convergence

This paper uses a probabilistic approach to analyze the converge of an ensemble Kalman filter solution to an exact Kalman filter solution in the simplest possible setting, the scalar case, as it allows us to build upon a rich literature of scalar probability distributions and non-elementary functions. To this end we introduce the bare-bones Scalar Pedagogical Ensemble Kalman Filter (SPEnKF). We show that in the asymptotic case of ensemble size, the expected value of both the analysis mean and variance estimate of the SPEnKF converges to that of the true Kalman filter, and that the variances of both tend towards zero, at each time moment. We also show that the ensemble converges in probability in the complementary case, when the ensemble is finite, and time is taken to infinity. Moreover, we show that in the finite-ensemble, finite-time case, variance inflation and mean correction can be leveraged to coerce the SPEnKF converge to its scalar Kalman filter counterpart. We then apply this framework to analyze perturbed observations and explain why perturbed observations ensemble Kalman filters underperform their deterministic counterparts.

0 citations0 reviewsNo code linkedOpen access
preprint2020arXiv

Convergence Results for Implicit--Explicit General Linear Methods

This paper studies fixed-step convergence of implicit-explicit general linear methods. We focus on a subclass of schemes that is internally consistent, has high stage order, and favorable stability properties. Classical, index-1 differential algebraic equation, and singular perturbation convergence analyses results are given. For all these problems IMEX GLMs from the class of interest converge with the full theoretical orders under general assumptions. The convergence results require the time steps to be sufficiently small, with upper bounds that are independent on the stiffness of the problem.

0 citations0 reviewsNo code linkedOpen access
preprint2020arXiv

Coupled Multirate Infinitesimal GARK Schemes for Stiff Systems with Multiple Time Scales

Traditional time discretization methods use a single timestep for the entire system of interest and can perform poorly when the dynamics of the system exhibits a wide range of time scales. Multirate infinitesimal step (MIS) methods (Knoth and Wolke, 1998) offer an elegant and flexible approach to efficiently integrate such systems. The slow components are discretized by a Runge-Kutta method, and the fast components are resolved by solving modified fast differential equations. Sandu (2018) developed the Multirate Infinitesimal General-structure Additive Runge-Kutta (MRI-GARK) family of methods that includes traditional MIS schemes as a subset. The MRI-GARK framework allowed the construction of the first fourth order MIS schemes. This framework also enabled the introduction of implicit methods, which are decoupled in the sense that any implicitness lies entirely within the fast or slow integrations. It was shown by Sandu that the stability of decoupled implicit MRI-GARK methods has limitations when both the fast and slow components are stiff and interact strongly. This work extends the MRI-GARK framework by introducing coupled implicit methods to solve stiff multiscale systems. The coupled approach has the potential to considerably improve the overall stability of the scheme, at the price of requiring implicit stage calculations over the entire system. Two coupling strategies are considered. The first computes coupled Runge-Kutta stages before solving a single differential equation to refine the fast solution. The second alternates between computing coupled Runge-Kutta stages and solving fast differential equations. We derive order conditions and perform the stability analysis for both strategies. The new coupled methods offer improved stability compared to the decoupled MRI-GARK schemes. The theoretical properties of the new methods are validated with numerical experiments.

0 citations0 reviewsNo code linkedOpen access
preprint2020arXiv

Goal-oriented a posteriori estimation of numerical errors in the solution of multiphysics systems

This paper develops a general methodology for a posteriori error estimation in time-dependent multiphysics numerical simulations. The methodology builds upon the generalized-structure additive Runge--Kutta (GARK) approach to time integration. GARK provides a unified formulation of multimethods that simulate complex systems by applying different discretization formulas and/or different time steps to individual components of the system. We derive discrete GARK adjoints and analyze their time accuracy. Based on the adjoint method, we establish computable a posteriori identities for the impacts of both temporal and spatial discretization errors on a given goal function. Numerical examples with reaction-diffusion systems illustrate the accuracy of the derived error measures. Local error decompositions are used to illustrate the power of this framework in adaptive refinements of both temporal and spatial meshes.

0 citations0 reviewsNo code linkedOpen access
preprint2020arXiv

Implicit multirate GARK methods

This work considers multirate generalized-structure additively partitioned Runge-Kutta (MrGARK) methods for solving stiff systems of ordinary differential equations (ODEs) with multiple time scales. These methods treat different partitions of the system with different timesteps for a more targeted and efficient solution compared to monolithic single rate approaches. With implicit methods used across all partitions, methods must find a balance between stability and the cost of solving nonlinear equations for the stages. In order to characterize this important trade-off, we explore multirate coupling strategies, problems for assessing linear stability, and techniques to efficiently implement Newton iterations for stage equations. Unlike much of the existing multirate stability analysis which is limited in scope to particular methods, we present general statements on stability and describe fundamental limitations for certain types of multirate schemes. New implicit multirate methods up to fourth order are derived, and their accuracy and efficiency properties are verified with numerical tests.

0 citations0 reviewsNo code linkedOpen access
preprint2020arXiv

Linearly implicit GARK schemes

Systems driven by multiple physical processes are central to many areas of science and engineering. Time discretization of multiphysics systems is challenging, since different processes have different levels of stiffness and characteristic time scales. The multimethod approach discretizes each physical process with an appropriate numerical method; the methods are coupled appropriately such that the overall solution has the desired accuracy and stability properties. The authors developed the general-structure additive Runge-Kutta (GARK) framework, which constructs multimethods based on Runge-Kutta schemes. This paper constructs the new GARK-ROS/GARK-ROW families of multimethods based on linearly implicit Rosenbrock/Rosenbrock-W schemes. For ordinary differential equation models, we develop a general order condition theory for linearly implicit methods with any number of partitions, using exact or approximate Jacobians. We generalize the order condition theory to two-way partitioned index-1 differential-algebraic equations. Applications of the framework include decoupled linearly implicit, linearly implicit/explicit, and linearly implicit/implicit methods. Practical GARK-ROS and GARK-ROW schemes of order up to four are constructed.

0 citations0 reviewsNo code linkedOpen access
preprint2020arXiv

Parallel implicit-explicit general linear methods

High-order discretizations of partial differential equations (PDEs) necessitate high-order time integration schemes capable of handling both stiff and nonstiff operators in an efficient manner. Implicit-explicit (IMEX) integration based on general linear methods (GLMs) offers an attractive solution due to their high stage and method order, as well as excellent stability properties. The IMEX characteristic allows stiff terms to be treated implicitly and nonstiff terms to be efficiently integrated explicitly. This work develops two systematic approaches for the development of IMEX GLMs of arbitrary order with stages that can be solved in parallel. The first approach is based on diagonally implicit multistage integration methods (DIMSIMs) of types 3 and 4. The second is a parallel generalization of IMEX Euler and has the interesting feature that the linear stability is independent of the order of accuracy. Numerical experiments confirm the theoretical rates of convergence and reveal that the new schemes are more efficient than serial IMEX GLMs and IMEX Runge-Kutta methods.

0 citations0 reviewsNo code linkedOpen access
preprint2017arXiv

Solving Parameter Estimation Problems with Discrete Adjoint Exponential Integrators

The solution of inverse problems in a variational setting finds best estimates of the model parameters by minimizing a cost function that penalizes the mismatch between model outputs and observations. The gradients required by the numerical optimization process are computed using adjoint models. Exponential integrators are a promising family of time discretizations for evolutionary partial differential equations. In order to allow the use of these discretizations in the context of inverse problems adjoints of exponential integrators are required. This work derives the discrete adjoint formulae for a W-type exponential propagation iterative methods of Runge-Kutta type (EPIRK-W). These methods allow arbitrary approximations of the Jacobian while maintaining the overall accuracy of the forward integration. The use of Jacobian approximation matrices that do not depend on the model state avoids the complex calculation of Hessians in the discrete adjoint formulae, and allows efficient adjoint code generation via algorithmic differentiation. We use the discrete EPIRK-W adjoints to solve inverse problems with the Lorenz-96 model and a computational magnetics benchmark test. Numerical results validate our theoretical derivations.

0 citations0 reviewsNo code linkedOpen access
preprint2016arXiv

A Parallel Implementation of the Ensemble Kalman Filter Based on Modified Cholesky Decomposition

This paper discusses an efficient parallel implementation of the ensemble Kalman filter based on the modified Cholesky decomposition. The proposed implementation starts with decomposing the domain into sub-domains. In each sub-domain a sparse estimation of the inverse background error covariance matrix is computed via a modified Cholesky decomposition; the estimates are computed concurrently on separate processors. The sparsity of this estimator is dictated by the conditional independence of model components for some radius of influence. Then, the assimilation step is carried out in parallel without the need of inter-processor communication. Once the local analysis states are computed, the analysis sub-domains are mapped back onto the global domain to obtain the analysis ensemble. Computational experiments are performed using the Atmospheric General Circulation Model (SPEEDY) with the T-63 resolution on the Blueridge cluster at Virginia Tech. The number of processors used in the experiments ranges from 96 to 2,048. The proposed implementation outperforms in terms of accuracy the well-known local ensemble transform Kalman filter (LETKF) for all the model variables. The computational time of the proposed implementation is similar to that of the parallel LETKF method (where no covariance estimation is performed). Finally, for the largest number of processors, the proposed parallel implementation is 400 times faster than the serial version of the proposed method.

0 citations0 reviewsNo code linkedOpen access
preprint2016arXiv

An Ensemble Kalman Filter Implementation Based on Modified Cholesky Decomposition for Inverse Covariance Matrix Estimation

This paper develops an efficient implementation of the ensemble Kalman filter based on a modified Cholesky decomposition for inverse covariance matrix estimation. This implementation is named EnKF-MC. Background errors corresponding to distant model components with respect to some radius of influence are assumed to be conditionally independent. This allows to obtain sparse estimators of the inverse background error covariance matrix. The computational effort of the proposed method is discussed and different formulations based on various matrix identities are provided. Furthermore, an asymptotic proof of convergence with regard to the ensemble size is presented. In order to assess the performance and the accuracy of the proposed method, experiments are performed making use of the Atmospheric General Circulation Model SPEEDY. The results are compared against those obtained using the local ensemble transform Kalman filter (LETKF). Tests are performed for dense observations ($100\%$ and $50\%$ of the model components are observed) as well as for sparse observations (only $12\%$, $6\%$, and $4\%$ of model components are observed). The results reveal that the use of modified Cholesky for inverse covariance matrix estimation can reduce the impact of spurious correlations during the assimilation cycle, i.e., the results of the proposed method are of better quality than those obtained via the LETKF in terms of root mean square error.

0 citations0 reviewsNo code linkedOpen access
preprint2016arXiv

Cluster Sampling Filters for Non-Gaussian Data Assimilation

This paper presents a fully non-Gaussian version of the Hamiltonian Monte Carlo (HMC) sampling filter. The Gaussian prior assumption in the original HMC filter is relaxed. Specifically, a clustering step is introduced after the forecast phase of the filter, and the prior density function is estimated by fitting a Gaussian Mixture Model (GMM) to the prior ensemble. Using the data likelihood function, the posterior density is then formulated as a mixture density, and is sampled using a HMC approach (or any other scheme capable of sampling multimodal densities in high-dimensional subspaces). The main filter developed herein is named "cluster HMC sampling filter" (ClHMC). A multi-chain version of the ClHMC filter, namely MC-ClHMC is also proposed to guarantee that samples are taken from the vicinities of all probability modes of the formulated posterior. The new methodologies are tested using a quasi-geostrophic (QG) model with double-gyre wind forcing and bi-harmonic friction. Numerical results demonstrate the usefulness of using GMMs to relax the Gaussian prior assumption in the HMC filtering paradigm.

0 citations0 reviewsNo code linkedOpen access
preprint2016arXiv

LIRK-W: Linearly-implicit Runge-Kutta methods with approximate matrix factorization

This paper develops a new class of linearly implicit time integration schemes called Linearly-Implicit Runge-Kutta-W (LIRK-W) methods. These schemes are based on an implicit-explicit approach which does not require a splitting of the right hand side and allow for arbitrary, time dependent, and stage varying approximations of the linear systems appearing in the method. Several formulations of LIRK-W schemes, each designed for specific approximation types, and their associated order condition theories are presented.

0 citations0 reviewsNo code linkedOpen access
preprint2016arXiv

The Reduced-Order Hybrid Monte Carlo Sampling Smoother

Hybrid Monte-Carlo (HMC) sampling smoother is a fully non-Gaussian four-dimensional data assimilation algorithm that works by directly sampling the posterior distribution formulated in the Bayesian framework. The smoother in its original formulation is computationally expensive due to the intrinsic requirement of running the forward and adjoint models repeatedly. Here we present computationally efficient versions of the HMC sampling smoother based on reduced-order approximations of the underlying model dynamics. The schemes developed herein are tested numerically using the shallow-water equations model on Cartesian coordinates. The results reveal that the reduced-order versions of the smoother are capable of accurately capturing the posterior probability density, while being significantly faster than the original full order formulation.

0 citations0 reviewsNo code linkedOpen access
preprint2015arXiv

A Derivative-Free Trust Region Framework for Variational Data Assimilation

This study develops a hybrid ensemble-variational approach for solving data assimilation problems. The method, called TR-4D-EnKF, is based on a trust region framework and consists of three computational steps. First an ensemble of model runs is propagated forward in time and snapshots of the state are stored. Next, a sequence of basis vectors is built and a low-dimensional representation of the data assimilation system is obtained by projecting the model state onto the space spanned by the ensemble deviations from the mean. Finally, the low-dimensional optimization problem is solved in the reduced-space using a trust region approach; the size of the trust region is updated according to the relative decrease of the reduced order surrogate cost function. The analysis state is projected back onto the full space, and the process is repeated with the current analysis serving as a new background. A heuristic approach based on the trust region size is proposed in order to adjust the background error statistics from one iteration to the next. Experimental simulations are carried out using the Lorenz and the quasi-geostrophic models. The results show that TR-4D-EnKF is an efficient computational approach, and is more accurate than the current state of the art 4D-EnKF implementations such as the POD-4D-EnKF and the Iterative Subspace Minimization methods.

0 citations0 reviewsNo code linkedOpen access
preprint2015arXiv

A Hybrid Monte-Carlo Sampling Smoother for Four Dimensional Data Assimilation

This paper constructs an ensemble-based sampling smoother for four-dimensional data assimilation using a Hybrid/Hamiltonian Monte-Carlo approach. The smoother samples efficiently from the posterior probability density of the solution at the initial time. Unlike the well-known ensemble Kalman smoother, which is optimal only in the linear Gaussian case, the proposed methodology naturally accommodates non-Gaussian errors and non-linear model dynamics and observation operators. Unlike the four-dimensional variational met\-hod, which only finds a mode of the posterior distribution, the smoother provides an estimate of the posterior uncertainty. One can use the ensemble mean as the minimum variance estimate of the state, or can use the ensemble in conjunction with the variational approach to estimate the background errors for subsequent assimilation windows. Numerical results demonstrate the advantages of the proposed method compared to the traditional variational and ensemble-based smoothing methods.

0 citations0 reviewsNo code linkedOpen access
preprint2015arXiv

A Time-parallel Approach to Strong-constraint Four-dimensional Variational Data Assimilation

A parallel-in-time algorithm based on an augmented Lagrangian approach is proposed to solve four-dimensional variational (4D-Var) data assimilation problems. The assimilation window is divided into multiple sub-intervals that allows to parallelize cost function and gradient computations. Solution continuity equations across interval boundaries are added as constraints. The augmented Lagrangian approach leads to a different formulation of the variational data assimilation problem than weakly constrained 4D-Var. A combination of serial and parallel 4D-Vars to increase performance is also explored. The methodology is illustrated on data assimilation problems with Lorenz-96 and the shallow water models.

0 citations0 reviewsNo code linkedOpen access
preprint2015arXiv

A-posteriori error estimates for inverse problems

Inverse problems use physical measurements along with a computational model to estimate the parameters or state of a system of interest. Errors in measurements and uncertainties in the computational model lead to inaccurate estimates. This work develops a methodology to estimate the impact of different errors on the variational solutions of inverse problems. The focus is on time evolving systems described by ordinary differential equations, and on a particular class of inverse problems, namely, data assimilation. The computational algorithm uses first-order and second-order adjoint models. In a deterministic setting the methodology provides a posteriori error estimates for the inverse solution. In a probabilistic setting it provides an a posteriori quantification of uncertainty in the inverse solution, given the uncertainties in the model and data. Numerical experiments with the shallow water equations in spherical coordinates illustrate the use of the proposed error estimation machinery in both deterministic and probabilistic settings.

0 citations0 reviewsNo code linkedOpen access
preprint2015arXiv

An Efficient Implementation of the Ensemble Kalman Filter Based on an Iterative Sherman-Morrison Formula

We present a practical implementation of the ensemble Kalman (EnKF) filter based on an iterative Sherman-Morrison formula. The new direct method exploits the special structure of the ensemble-estimated error covariance matrices in order to efficiently solve the linear systems involved in the analysis step of the EnKF. The computational complexity of the proposed implementation is equivalent to that of the best EnKF implementations available in the literature when the number of observations is much larger than the number of ensemble members. Even when this conditions is not fulfilled, the proposed method is expected to perform well since it does not employ matrix decompositions. Computational experiments using the Lorenz 96 and the oceanic quasi-geostrophic models are performed in order to compare the proposed algorithm with EnKF implementations that use matrix decompositions. In terms of accuracy, the results of all implementations are similar. The proposed method is considerably faster than other EnKF variants, even when the number of observations is large relative to the number of ensemble members.

0 citations0 reviewsNo code linkedOpen access
preprint2015arXiv

Efficient approximation of sparse Jacobians for time-implicit reduced order models

This paper introduces a sparse matrix discrete interpolation method to effectively compute matrix approximations in the reduced order modeling framework. The sparse algorithm developed herein relies on the discrete empirical interpolation method and uses only samples of the nonzero entries of the matrix series. The proposed approach can approximate very large matrices, unlike the current matrix discrete empirical interpolation method which is limited by its large computational memory requirements. The empirical interpolation indexes obtained by the sparse algorithm slightly differ from the ones computed by the matrix discrete empirical interpolation method as a consequence of the singular vectors round-off errors introduced by the economy or full singular value decomposition (SVD) algorithms when applied to the full matrix snapshots. When appropriately padded with zeros the economy SVD factorization of the nonzero elements of the snapshots matrix is a valid economy SVD for the full snapshots matrix. Numerical experiments are performed with the 1D Burgers and 2D Shallow Water Equations test problems where the quadratic reduced nonlinearities are computed via tensorial calculus. The sparse matrix approximation strategy is compared against five existing methods for computing reduced Jacobians: a) matrix discrete empirical interpolation method, b) discrete empirical interpolation method, c) tensorial calculus, d) full Jacobian projection onto the reduced basis subspace, and e) directional derivatives of the model along the reduced basis functions. The sparse matrix method outperforms all other algorithms. The use of traditional matrix discrete empirical interpolation method is not possible for very large instances due to its excessive memory requirements.

0 citations0 reviewsNo code linkedOpen access
preprint2015arXiv

Efficient Construction of Local Parametric Reduced Order Models Using Machine Learning Techniques

Reduced order models are computationally inexpensive approximations that capture the important dynamical characteristics of large, high-fidelity computer models of physical systems. This paper applies machine learning techniques to improve the design of parametric reduced order models. Specifically, machine learning is used to develop feasible regions in the parameter space where the admissible target accuracy is achieved with a predefined reduced order basis, to construct parametric maps, to chose the best two already existing bases for a new parameter configuration from accuracy point of view and to pre-select the optimal dimension of the reduced basis such as to meet the desired accuracy. By combining available information using bases concatenation and interpolation as well as high-fidelity solutions interpolation we are able to build accurate reduced order models associated with new parameter settings. Promising numerical results with a viscous Burgers model illustrate the potential of machine learning approaches to help design better reduced order models.

0 citations0 reviewsNo code linkedOpen access
preprint2015arXiv

Ensemble Kalman Filter Implementations Based on Covariance Matrix Estimation

This paper develops efficient ensemble Kalman filter (EnKF) implementations based on shrinkage covariance estimation. The forecast ensemble members at each step are used to estimate the background error covariance matrix via the Rao-Blackwell Ledoit and Wolf estimator, which has been developed specifically developed to approximate high-dimensional covariance matrices using a small number of samples. Additional samples are taken from the normal distribution described by the background ensemble mean and the estimated background covariance matrix in order to increase the size of the ensemble and reduce the sampling error of the filter. This increase in the size of the ensemble is obtained without running the forward model. After the assimilation step, the additional samples are discarded and only the initial members are propagated. Two implementations are considered. In the EnKF Full-Space (EnKF-FS) approach the assimilation process is performed in the model space, while the EnKF Reduce-Space (EnKF-RS) formulation performs the analysis in the subspace spanned by the ensemble members. Numerical experiments carried out with a quasi-geostrophic model show that the proposed implementations outperform current methods such as the traditional EnKF formulation, square root filters, and inflation-free EnKF implementations. The proposed implementations provide good results with small ensemble sizes ($\sim 10$) and small percentages of observed components from the vector state. These results are similar (and in some cases better) to traditional methods using large ensemble sizes ($\sim 80$) and large percentages of observed components. The computational times of the new implementations remain reasonably low.

0 citations0 reviewsNo code linkedOpen access
preprint2015arXiv

Exponential-Krylov methods for ordinary differential equations

This paper develops a new class of exponential-type integrators where all the matrix exponentiations are performed in a single Krylov space of low dimension. The new family, called Lightly Implicit Krylov-Exponential (LIKE), is well suited for solving large scale systems of ODEs or semi-discrete PDEs. The time discretization and the Krylov space approximation are treated as a single computational process, and the Krylov space properties are an integral part of the new LIKE order condition theory developed herein. Consequently, LIKE methods require a small number of basis vectors determined solely by the temporal order of accuracy. The subspace size is independent of the ODE under consideration, and there is no need to monitor the errors in linear system solutions at each stage. Numerical results illustrate the favorable properties of new family of methods.

0 citations0 reviewsNo code linkedOpen access
preprint2015arXiv

POD/DEIM Reduced-Order Strategies for Efficient Four Dimensional Variational Data Assimilation

This work studies reduced order modeling (ROM) approaches to speed up the solution of variational data assimilation problems with large scale nonlinear dynamical models. It is shown that a key requirement for a successful reduced order solution is that reduced order Karush-Kuhn-Tucker conditions accurately represent their full order counterparts. In particular, accurate reduced order approximations are needed for the forward and adjoint dynamical models, as well as for the reduced gradient. New strategies to construct reduced order based are developed for Proper Orthogonal Decomposition (POD) ROM data assimilation using both Galerkin and Petrov-Galerkin projections. For the first time POD, tensorial POD, and discrete empirical interpolation method (DEIM) are employed to develop reduced data assimilation systems for a geophysical flow model, namely, the two dimensional shallow water equations. Numerical experiments confirm the theoretical framework for Galerkin projection. In the case of Petrov-Galerkin projection, stabilization strategies must be considered for the reduced order models. The new reduced order shallow water data assimilation system provides analyses similar to those produced by the full resolution data assimilation system in one tenth of the computational time.

0 citations0 reviewsNo code linkedOpen access
preprint2015arXiv

Robust data assimilation using $L_1$ and Huber norms

Data assimilation is the process to fuse information from priors, observations of nature, and numerical models, in order to obtain best estimates of the parameters or state of a physical system of interest. Presence of large errors in some observational data, e.g., data collected from a faulty instrument, negatively affect the quality of the overall assimilation results. This work develops a systematic framework for robust data assimilation. The new algorithms continue to produce good analyses in the presence of observation outliers. The approach is based on replacing the traditional $Ł_2$ norm formulation of data assimilation problems with formulations based on $Ł_1$ and Huber norms. Numerical experiments using the Lorenz-96 and the shallow water on the sphere models illustrate how the new algorithms outperform traditional data assimilation approaches in the presence of data outliers.

0 citations0 reviewsNo code linkedOpen access
preprint2015arXiv

Rosenbrock-Krylov Methods for Large Systems of Differential Equations

This paper develops a new class of Rosenbrock-type integrators based on a Krylov space solution of the linear systems. The new family, called Rosenbrock-Krylov (Rosenbrock-K), is well suited for solving large scale systems of ODEs or semi-discrete PDEs. The time discretization and the Krylov space approximation are treated as a single computational process, and the Krylov space properties are an integral part of the new Rosenbrock-K order condition theory developed herein. Consequently, Rosenbrock-K methods require a small number of basis vectors determined solely by the temporal order of accuracy. The subspace size is independent of the ODE under consideration, and there is no need to monitor the errors in linear system solutions at each stage. Numerical results show favorable properties of Rosenbrock-K methods when compared to current Rosenbrock and Rosenbrock-W schemes.

0 citations0 reviewsNo code linkedOpen access
preprint2014arXiv

A Sampling Filter for Non-Gaussian Data Assimilation

Data assimilation combines information from models, measurements, and priors to estimate the state of a dynamical system such as the atmosphere. The Ensemble Kalman filter (EnKF) is a family of ensemble-based data assimilation approaches that has gained wide popularity due its simple formulation, ease of implementation, and good practical results. Most EnKF algorithms assume that the underlying probability distributions are Gaussian. Although this assumption is well accepted, it is too restrictive when applied to large nonlinear models, nonlinear observation operators, and large levels of uncertainty. Several approaches have been proposed in order to avoid the Gaussianity assumption. One of the most successful strategies is the maximum likelihood ensemble filter (MLEF) which computes a maximum a posteriori estimate of the state assuming the posterior distribution is Gaussian. MLEF is designed to work with nonlinear and even non-differentiable observation operators, and shows good practical performance. However, there are limits to the degree of nonlinearity that MLEF can handle. This paper proposes a new ensemble-based data assimilation method, named the "sampling filter", which obtains the analysis by sampling directly from the posterior distribution. The sampling strategy is based on a Hybrid Monte Carlo (HMC) approach that can handle non-Gaussian probability distributions. Numerical experiments are carried out using the Lorenz-96 model and observation operators with different levels of non-linearity and differentiability. The proposed filter is also tested with shallow water model on a sphere with linear observation operator. The results show that the sampling filter can perform well even in highly nonlinear situations were EnKF and MLEF filters diverge.

0 citations0 reviewsNo code linkedOpen access
preprint2014arXiv

Application of approximate matrix factorization to high order linearly implicit Runge-Kutta methods

Linearly implicit Runge-Kutta methods with approximate matrix factorization can solve efficiently large systems of differential equations that have a stiff linear part, e.g. reaction-diffusion systems. However, the use of approximate factorization usually leads to loss of accuracy, which makes it attractive only for low order time integration schemes. This paper discusses the application of approximate matrix factorization with high order methods; an inexpensive correction procedure applied to each stage allows to retain the high order of the underlying linearly implicit Runge-Kutta scheme. The accuracy and stability of the methods are studied. Numerical experiments on reaction-diffusion type problems of different sizes and with different degrees of stiffness illustrate the efficiency of the proposed approach.

0 citations0 reviewsNo code linkedOpen access
preprint2014arXiv

Approximate Exponential Algorithms to Solve the Chemical Master Equation

This paper discusses new simulation algorithms for stochastic chemical kinetics that exploit the linearity of the chemical master equation and its matrix exponential exact solution. These algorithms make use of various approximations of the matrix exponential to evolve probability densities in time. A sampling of the approximate solutions of the chemical master equation is used to derive accelerated stochastic simulation algorithms. Numerical experiments compare the new methods with the established stochastic simulation algorithm and the tau-leaping method.

0 citations0 reviewsNo code linkedOpen access
preprint2014arXiv

Comparison of POD reduced order strategies for the nonlinear 2D Shallow Water Equations

This paper introduces tensorial calculus techniques in the framework of Proper Orthogonal Decomposition (POD) to reduce the computational complexity of the reduced nonlinear terms. The resulting method, named tensorial POD, can be applied to polynomial nonlinearities of any degree $p$. Such nonlinear terms have an on-line complexity of $\mathcal{O}(k^{p+1})$, where $k$ is the dimension of POD basis, and therefore is independent of full space dimension. However it is efficient only for quadratic nonlinear terms since for higher nonlinearities standard POD proves to be less time consuming once the POD basis dimension $k$ is increased. Numerical experiments are carried out with a two dimensional shallow water equation (SWE) test problem to compare the performance of tensorial POD, standard POD, and POD/Discrete Empirical Interpolation Method (DEIM). Numerical results show that tensorial POD decreases by $76\times$ times the computational cost of the on-line stage of standard POD for configurations using more than $300,000$ model variables. The tensorial POD SWE model was only $2-8\times$ slower than the POD/DEIM SWE model but the implementation effort is considerably increased. Tensorial calculus was again employed to construct a new algorithm allowing POD/DEIM shallow water equation model to compute its off-line stage faster than the standard and tensorial POD approaches.

0 citations0 reviewsNo code linkedOpen access
preprint2014arXiv

Dynamic Response Optimization of Complex Multibody Systems in a Penalty Formulation using Adjoint Sensitivity

Multibody dynamics simulations are currently widely accepted as valuable means for dynamic performance analysis of mechanical systems. The evolution of theoretical and computational aspects of the multibody dynamics discipline make it conducive these days for other types of applications, in addition to pure simulations. One very important such application is design optimization. A very important first step towards design optimization is sensitivity analysis of multibody system dynamics. Dynamic sensitivities are often calculated by means of finite differences. Depending of the number of parameters involved, this procedure can be computationally expensive. Moreover, in many cases, the results suffer from low accuracy when real perturbations are used. The main contribution to the state-of-the-art brought by this study is the development of the adjoint sensitivity approach of multibody systems in the context of the penalty formulation. The theory developed is demonstrated on one academic case study, a five-bar mechanism, and on one real-life system, a 14-DOF vehicle model. The five-bar mechanism is used to illustrate the sensitivity approach derived in this paper. The full vehicle model is used to demonstrate the capability of the new approach developed to perform sensitivity analysis and gradient-based optimization for large and complex multibody systems with respect to multiple design parameters.

0 citations0 reviewsNo code linkedOpen access
preprint2014arXiv

High Order Implicit-Explicit General Linear Methods with Optimized Stability Regions

In the numerical solution of partial differential equations using a method-of-lines approach, the availability of high order spatial discretization schemes motivates the development of sophisticated high order time integration methods. For multiphysics problems with both stiff and non-stiff terms implicit-explicit (IMEX) time stepping methods attempt to combine the lower cost advantage of explicit schemes with the favorable stability properties of implicit schemes. Existing high order IMEX Runge Kutta or linear multistep methods, however, suffer from accuracy or stability reduction. This work shows that IMEX general linear methods (GLMs) are competitive alternatives to classic IMEX schemes for large problems arising in practice. High order IMEX-GLMs are constructed in the framework developed by the authors [34]. The stability regions of the new schemes are optimized numerically. The resulting IMEX-GLMs have similar stability properties as IMEX Runge-Kutta methods, but they do not suffer from order reduction, and are superior in terms of accuracy and efficiency. Numerical experiments with two and three dimensional test problems illustrate the potential of the new schemes to speed up complex applications.

0 citations0 reviewsNo code linkedOpen access
preprint2014arXiv

Optimization of Vehicle Dynamics based on Multibody Models using Adjoint Sensitivity Analysis

Multibody dynamics simulations have become widely used tools for vehicle systems analysis and design. As this approach evolves, it becomes able to provide additional information for various types of analyses. One very important direction is the optimization of multibody systems. Sensitivity analysis of multibody system dynamics is essential for design optimization. Dynamic sensitivities, when needed, are often calculated by means of finite differences. However, depending of the number of parameters involved, this procedure can be computationally expensive. Moreover, in many cases the results suffer from low accuracy when real perturbations are used. This paper develops the adjoint sensitivity analysis of multibody systems in the context of penalty formulations. The resulting sensitivities are applied to perform dynamical optimization of a full vehicle system.

0 citations0 reviewsNo code linkedOpen access
preprint2014arXiv

Solving stochastic chemical kinetics by Metropolis Hastings sampling

This study considers using Metropolis-Hastings algorithm for stochastic simulation of chemical reactions. The proposed method uses SSA (Stochastic Simulation Algorithm) distribution which is a standard method for solving well-stirred chemically reacting systems as a desired distribution. A new numerical solvers based on exponential form of exact and approximate solutions of CME (Chemical Master Equation) is employed for obtaining target and proposal distributions in Metropolis-Hastings algorithm to accelerate the accuracy of the tau-leap method. Samples generated by this technique have the same distribution as SSA and the histogram of samples show it's convergence to SSA.

0 citations0 reviewsNo code linkedOpen access
preprint2013arXiv

A class of generalized additive Runge-Kutta methods

This work generalizes the additively partitioned Runge-Kutta methods by allowing for different stage values as arguments of different components of the right hand side. An order conditions theory is developed for the new family of generalized additive methods, and stability and monotonicity investigations are carried out. The paper discusses the construction and properties of implicit-explicit and implicit-implicit,methods in the new framework. The new family, named GARK, introduces additional flexibility when compared to traditional partitioned Runge-Kutta methods, and therefore offers additional opportunities for the development of flexible solvers for systems with multiple scales, or driven by multiple physical processes.

0 citations0 reviewsNo code linkedOpen access
preprint2013arXiv

An Optimization Framework to Improve 4D-Var Data Assimilation System Performance

This paper develops a computational framework for optimizing the parameters of data assimilation systems in order to improve their performance. The approach formulates a continuous meta-optimization problem for parameters; the meta-optimization is constrained by the original data assimilation problem. The numerical solution process employs adjoint models and iterative solvers. The proposed framework is applied to optimize observation values, data weighting coefficients, and the location of sensors for a test problem. The ability to optimize a distributed measurement network is crucial for cutting down operating costs and detecting malfunctions.

0 citations0 reviewsNo code linkedOpen access
preprint2013arXiv

Efficient methods for computing observation impact in 4D-Var data assimilation

This paper presents a practical computational approach to quantify the effect of individual observations in estimating the state of a system. Such an analysis can be used for pruning redundant measurements, and for designing future sensor networks. The mathematical approach is based on computing the sensitivity of the reanalysis (unconstrained optimization solution) with respect to the data. The computational cost is dominated by the solution of a linear system, whose matrix is the Hessian of the cost function, and is only available in operator form. The right hand side is the gradient of a scalar cost function that quantifies the forecast error of the numerical model. The use of adjoint models to obtain the necessary first and second order derivatives is discussed. We study various strategies to accelerate the computation, including matrix-free iterative solvers, preconditioners, and an in-house multigrid solver. Experiments are conducted on both a small-size shallow-water equations model, and on a large-scale numerical weather prediction model, in order to illustrate the capabilities of the new methodology.

0 citations0 reviewsNo code linkedOpen access
preprint2013arXiv

Extrapolation-based implicit-explicit general linear methods

For many systems of differential equations modeling problems in science and engineering, there are natural splittings of the right hand side into two parts, one non-stiff or mildly stiff, and the other one stiff. For such systems implicit-explicit (IMEX) integration combines an explicit scheme for the non-stiff part with an implicit scheme for the stiff part. In a recent series of papers two of the authors (Sandu and Zhang) have developed IMEX GLMs, a family of implicit-explicit schemes based on general linear methods. It has been shown that, due to their high stage order, IMEX GLMs require no additional coupling order conditions, and are not marred by order reduction. This work develops a new extrapolation-based approach to construct practical IMEX GLM pairs of high order. We look for methods with large absolute stability region, assuming that the implicit part of the method is A- or L-stable. We provide examples of IMEX GLMs with optimal stability properties. Their application to a two dimensional test problem confirms the theoretical findings.

0 citations0 reviewsNo code linkedOpen access
preprint2013arXiv

Implicit Simulation Methods for Stochastic Chemical Kinetics

In biochemical systems some of the chemical species are present with only small numbers of molecules. In this situation discrete and stochastic simulation approaches are more relevant than continuous and deterministic ones. The fundamental Gillespie's stochastic simulation algorithm (SSA) accounts for every reaction event, which occurs with a probability determined by the configuration of the system. This approach requires a considerable computational effort for models with many reaction channels and chemical species. In order to improve efficiency, tau-leaping methods represent multiple firings of each reaction during a simulation step by Poisson random variables. For stiff systems the mean of this variable is treated implicitly in order to ensure numerical stability. This paper develops fully implicit tau-leaping-like algorithms that treat implicitly both the mean and the variance of the Poisson variables. The construction is based on adapting weakly convergent discretizations of stochastic differential equations to stochastic chemical kinetic systems. Theoretical analyses of accuracy and stability of the new methods are performed on a standard test problem. Numerical results demonstrate the performance of the proposed tau-leaping methods.

0 citations0 reviewsNo code linkedOpen access
preprint2013arXiv

Low-rank Approximations for Computing Observation Impact in 4D-Var Data Assimilation

We present an efficient computational framework to quantify the impact of individual observations in four dimensional variational data assimilation. The proposed methodology uses first and second order adjoint sensitivity analysis, together with matrix-free algorithms to obtain low-rank approximations of ob- servation impact matrix. We illustrate the application of this methodology to important applications such as data pruning and the identification of faulty sensors for a two dimensional shallow water test system.

0 citations0 reviewsNo code linkedOpen access
preprint2013arXiv

Multirate generalized additive Runge Kutta methods

This work constructs a new class of multirate schemes based on the recently developed generalized additive Runge-Kutta (GARK) methods (Sandu and Guenther, 2013). Multirate schemes use different step sizes for different components and for different partitions of the right-hand side based on the local activity levels. We show that the new multirate GARK family includes many well-known multirate schemes as special cases. The order conditions theory follows directly from the GARK accuracy theory. Nonlinear stability and monotonicity investigations show that these properties are inherited from the base schemes provided that additional coupling conditions hold.

0 citations0 reviewsNo code linkedOpen access
preprint2013arXiv

Partitioned and implicit-explicit general linear methods for ordinary differential equations

Implicit-explicit (IMEX) time stepping methods can efficiently solve differential equa- tions with both stiff and nonstiff components. IMEX Runge-Kutta methods and IMEX linear multistep methods have been studied in the literature. In this pa- per we study new implicit-explicit methods of general linear type (IMEX-GLMs). We develop an order conditions theory for high stage order partitioned GLMs that share the same abscissae, and show that no additional coupling order conditions are needed. Consequently, GLMs offer an excellent framework for the construction of multi-method integration algorithms. Next, we propose a family of IMEX schemes based on diagonally-implicit multi-stage integration methods and construct practical schemes of order three. Numerical results confirm the theoretical findings.

0 citations0 reviewsNo code linkedOpen access