@article{WormellReich2021, author = {Wormell, Caroline L. and Reich, Sebastian}, title = {Spectral convergence of diffusion maps}, series = {SIAM journal on numerical analysis / Society for Industrial and Applied Mathematics}, volume = {59}, journal = {SIAM journal on numerical analysis / Society for Industrial and Applied Mathematics}, number = {3}, publisher = {Society for Industrial and Applied Mathematics}, address = {Philadelphia}, issn = {0036-1429}, doi = {10.1137/20M1344093}, pages = {1687 -- 1734}, year = {2021}, abstract = {Diffusion maps is a manifold learning algorithm widely used for dimensionality reduction. Using a sample from a distribution, it approximates the eigenvalues and eigenfunctions of associated Laplace-Beltrami operators. Theoretical bounds on the approximation error are, however, generally much weaker than the rates that are seen in practice. This paper uses new approaches to improve the error bounds in the model case where the distribution is supported on a hypertorus. For the data sampling (variance) component of the error we make spatially localized compact embedding estimates on certain Hardy spaces; we study the deterministic (bias) component as a perturbation of the Laplace-Beltrami operator's associated PDE and apply relevant spectral stability results. Using these approaches, we match long-standing pointwise error bounds for both the spectral data and the norm convergence of the operator discretization. We also introduce an alternative normalization for diffusion maps based on Sinkhorn weights. This normalization approximates a Langevin diffusion on the sample and yields a symmetric operator approximation. We prove that it has better convergence compared with the standard normalization on flat domains, and we present a highly efficient rigorous algorithm to compute the Sinkhorn weights.}, language = {en} } @misc{vanLeeuwenKunschNergeretal.2019, author = {van Leeuwen, Peter Jan and Kunsch, Hans R. and Nerger, Lars and Potthast, Roland and Reich, Sebastian}, title = {Particle filters for high-dimensional geoscience applications: A review}, series = {Quarterly journal of the Royal Meteorological Society}, volume = {145}, journal = {Quarterly journal of the Royal Meteorological Society}, number = {723}, publisher = {Wiley}, address = {Hoboken}, issn = {0035-9009}, doi = {10.1002/qj.3551}, pages = {2335 -- 2365}, year = {2019}, abstract = {Particle filters contain the promise of fully nonlinear data assimilation. They have been applied in numerous science areas, including the geosciences, but their application to high-dimensional geoscience systems has been limited due to their inefficiency in high-dimensional systems in standard settings. However, huge progress has been made, and this limitation is disappearing fast due to recent developments in proposal densities, the use of ideas from (optimal) transportation, the use of localization and intelligent adaptive resampling strategies. Furthermore, powerful hybrids between particle filters and ensemble Kalman filters and variational methods have been developed. We present a state-of-the-art discussion of present efforts of developing particle filters for high-dimensional nonlinear geoscience state-estimation problems, with an emphasis on atmospheric and oceanic applications, including many new ideas, derivations and unifications, highlighting hidden connections, including pseudo-code, and generating a valuable tool and guide for the community. Initial experiments show that particle filters can be competitive with present-day methods for numerical weather prediction, suggesting that they will become mainstream soon.}, language = {en} } @book{VanLeeuwenChengReich2015, author = {Van Leeuwen, Peter Jan and Cheng, Yuan and Reich, Sebastian}, title = {Nonlinear data assimilation}, series = {Frontiers in applied dynamical systems: reviews and tutorials ; 2}, journal = {Frontiers in applied dynamical systems: reviews and tutorials ; 2}, publisher = {Springer}, address = {Cham}, isbn = {978-3-319-18346-6}, doi = {10.1007/978-3-319-18347-3}, pages = {xii, 118}, year = {2015}, abstract = {This book contains two review articles on nonlinear data assimilation that deal with closely related topics but were written and can be read independently. Both contributions focus on so-called particle filters. The first contribution by Jan van Leeuwen focuses on the potential of proposal densities. It discusses the issues with present-day particle filters and explorers new ideas for proposal densities to solve them, converging to particle filters that work well in systems of any dimension, closing the contribution with a high-dimensional example. The second contribution by Cheng and Reich discusses a unified framework for ensemble-transform particle filters. This allows one to bridge successful ensemble Kalman filters with fully nonlinear particle filters, and allows a proper introduction of localization in particle filters, which has been lacking up to now.}, language = {en} } @article{TaghvaeideWiljesMehtaetal.2017, author = {Taghvaei, Amirhossein and de Wiljes, Jana and Mehta, Prashant G. and Reich, Sebastian}, title = {Kalman filter and its modern extensions for the continuous-time nonlinear filtering problem}, series = {Journal of dynamic systems measurement and control}, volume = {140}, journal = {Journal of dynamic systems measurement and control}, number = {3}, publisher = {ASME}, address = {New York}, issn = {0022-0434}, doi = {10.1115/1.4037780}, pages = {11}, year = {2017}, abstract = {This paper is concerned with the filtering problem in continuous time. Three algorithmic solution approaches for this problem are reviewed: (i) the classical Kalman-Bucy filter, which provides an exact solution for the linear Gaussian problem; (ii) the ensemble Kalman-Bucy filter (EnKBF), which is an approximate filter and represents an extension of the Kalman-Bucy filter to nonlinear problems; and (iii) the feedback particle filter (FPF), which represents an extension of the EnKBF and furthermore provides for a consistent solution in the general nonlinear, non-Gaussian case. The common feature of the three algorithms is the gain times error formula to implement the update step (to account for conditioning due to the observations) in the filter. In contrast to the commonly used sequential Monte Carlo methods, the EnKBF and FPF avoid the resampling of the particles in the importance sampling update step. Moreover, the feedback control structure provides for error correction potentially leading to smaller simulation variance and improved stability properties. The paper also discusses the issue of nonuniqueness of the filter update formula and formulates a novel approximation algorithm based on ideas from optimal transport and coupling of measures. Performance of this and other algorithms is illustrated for a numerical example.}, language = {en} } @article{StaniforthWoodReich2006, author = {Staniforth, Andrew and Wood, Nigel and Reich, Sebastian}, title = {A time-staggered semi-Lagrangian discretization of the rotating shallow-water equations}, series = {Quarterly journal of the Royal Meteorological Society}, volume = {132}, journal = {Quarterly journal of the Royal Meteorological Society}, number = {621C}, publisher = {Wiley}, address = {Weinheim}, issn = {0035-9009}, doi = {10.1256/qj.06.30}, pages = {3107 -- 3116}, year = {2006}, abstract = {A time-staggered semi-Lagrangian discretization of the rotating shallow-water equations is proposed and analysed. Application of regularization to the geopotential field used in the momentum equations leads to an unconditionally stable scheme. The analysis, together with a fully nonlinear example application, suggests that this approach is a promising, efficient, and accurate alternative to traditional schemes.}, language = {en} } @article{SomogyvariReich2020, author = {Somogyv{\´a}ri, M{\´a}rk and Reich, Sebastian}, title = {Convergence tests for transdimensional Markov chains in geoscience imaging}, series = {Mathematical geosciences : the official journal of the International Association for Mathematical Geosciences}, volume = {52}, journal = {Mathematical geosciences : the official journal of the International Association for Mathematical Geosciences}, number = {5}, publisher = {Springer}, address = {Heidelberg}, issn = {1874-8961}, doi = {10.1007/s11004-019-09811-x}, pages = {651 -- 668}, year = {2020}, abstract = {Classic inversion methods adjust a model with a predefined number of parameters to the observed data. With transdimensional inversion algorithms such as the reversible-jump Markov chain Monte Carlo (rjMCMC), it is possible to vary this number during the inversion and to interpret the observations in a more flexible way. Geoscience imaging applications use this behaviour to automatically adjust model resolution to the inhomogeneities of the investigated system, while keeping the model parameters on an optimal level. The rjMCMC algorithm produces an ensemble as result, a set of model realizations, which together represent the posterior probability distribution of the investigated problem. The realizations are evolved via sequential updates from a randomly chosen initial solution and converge toward the target posterior distribution of the inverse problem. Up to a point in the chain, the realizations may be strongly biased by the initial model, and must be discarded from the final ensemble. With convergence assessment techniques, this point in the chain can be identified. Transdimensional MCMC methods produce ensembles that are not suitable for classic convergence assessment techniques because of the changes in parameter numbers. To overcome this hurdle, three solutions are introduced to convert model realizations to a common dimensionality while maintaining the statistical characteristics of the ensemble. A scalar, a vector and a matrix representation for models is presented, inferred from tomographic subsurface investigations, and three classic convergence assessment techniques are applied on them. It is shown that appropriately chosen scalar conversions of the models could retain similar statistical ensemble properties as geologic projections created by rasterization.}, language = {en} } @article{SkeelReich2011, author = {Skeel, R. D. and Reich, Sebastian}, title = {Corrected potential energy functions for constrained molecular dynamics}, series = {European physical journal special topics}, volume = {200}, journal = {European physical journal special topics}, number = {1}, publisher = {Springer}, address = {Heidelberg}, issn = {1951-6355}, doi = {10.1140/epjst/e2011-01518-8}, pages = {55 -- 72}, year = {2011}, abstract = {Atomic oscillations present in classical molecular dynamics restrict the step size that can be used. Multiple time stepping schemes offer only modest improvements, and implicit integrators are costly and inaccurate. The best approach may be to actually remove the highest frequency oscillations by constraining bond lengths and bond angles, thus permitting perhaps a 4-fold increase in the step size. However, omitting degrees of freedom produces errors in statistical averages, and rigid angles do not bend for strong excluded volume forces. These difficulties can be addressed by an enhanced treatment of holonomic constrained dynamics using ideas from papers of Fixman (1974) and Reich (1995, 1999). In particular, the 1995 paper proposes the use of "flexible" constraints, and the 1999 paper uses a modified potential energy function with rigid constraints to emulate flexible constraints. Presented here is a more direct and rigorous derivation of the latter approach, together with justification for the use of constraints in molecular modeling. With rigor comes limitations, so practical compromises are proposed: simplifications of the equations and their judicious application when assumptions are violated. Included are suggestions for new approaches.}, language = {en} } @article{ShinZoellerHolschneideretal.2011, author = {Shin, Seoleun and Z{\"o}ller, Gert and Holschneider, Matthias and Reich, Sebastian}, title = {A multigrid solver for modeling complex interseismic stress fields}, series = {Computers \& geosciences : an international journal devoted to the publication of papers on all aspects of geocomputation and to the distribution of computer programs and test data sets ; an official journal of the International Association for Mathematical Geology}, volume = {37}, journal = {Computers \& geosciences : an international journal devoted to the publication of papers on all aspects of geocomputation and to the distribution of computer programs and test data sets ; an official journal of the International Association for Mathematical Geology}, number = {8}, publisher = {Elsevier}, address = {Oxford}, issn = {0098-3004}, doi = {10.1016/j.cageo.2010.11.011}, pages = {1075 -- 1082}, year = {2011}, abstract = {We develop a multigrid, multiple time stepping scheme to reduce computational efforts for calculating complex stress interactions in a strike-slip 2D planar fault for the simulation of seismicity. The key elements of the multilevel solver are separation of length scale, grid-coarsening, and hierarchy. In this study the complex stress interactions are split into two parts: the first with a small contribution is computed on a coarse level, and the rest for strong interactions is on a fine level. This partition leads to a significant reduction of the number of computations. The reduction of complexity is even enhanced by combining the multigrid with multiple time stepping. Computational efficiency is enhanced by a factor of 10 while retaining a reasonable accuracy, compared to the original full matrix-vortex multiplication. The accuracy of solution and computational efficiency depend on a given cut-off radius that splits multiplications into the two parts. The multigrid scheme is constructed in such a way that it conserves stress in the entire half-space.}, language = {en} } @article{ShinSommerReichetal.2010, author = {Shin, Seoleun and Sommer, Matthias and Reich, Sebastian and N{\´e}vir, Peter}, title = {Evaluation of three spatial discretization schemes with the Galewsky et al. test}, issn = {1530-261X}, doi = {10.1002/Asl.279}, year = {2010}, abstract = {We evaluate the Hamiltonian particle methods (HPM) and the Nambu discretization applied to shallow-water equations on the sphere using the test suggested by Galewsky et al. (2004). Both simulations show excellent conservation of energy and are stable in long-term simulation. We repeat the test also using the ICOSWP scheme to compare with the two conservative spatial discretization schemes. The HPM simulation captures the main features of the reference solution, but wave 5 pattern is dominant in the simulations applied on the ICON grid with relatively low spatial resolutions. Nevertheless, agreement in statistics between the three schemes indicates their qualitatively similar behaviors in the long-term integration.}, language = {en} } @article{ShinReichFrank2012, author = {Shin, Seoleun and Reich, Sebastian and Frank, Jason}, title = {Hydrostatic Hamiltonian particle-mesh (HPM) methods for atmospheric modelling}, series = {Quarterly journal of the Royal Meteorological Society}, volume = {138}, journal = {Quarterly journal of the Royal Meteorological Society}, number = {666}, publisher = {Wiley-Blackwell}, address = {Hoboken}, issn = {0035-9009}, doi = {10.1002/qj.982}, pages = {1388 -- 1399}, year = {2012}, abstract = {We develop a hydrostatic Hamiltonian particle-mesh (HPM) method for efficient long-term numerical integration of the atmosphere. In the HPM method, the hydrostatic approximation is interpreted as a holonomic constraint for the vertical position of particles. This can be viewed as defining a set of vertically buoyant horizontal meshes, with the altitude of each mesh point determined so as to satisfy the hydrostatic balance condition and with particles modelling horizontal advection between the moving meshes. We implement the method in a vertical-slice model and evaluate its performance for the simulation of idealized linear and nonlinear orographic flow in both dry and moist environments. The HPM method is able to capture the basic features of the gravity wave to a degree of accuracy comparable with that reported in the literature. The numerical solution in the moist experiment indicates that the influence of moisture on wave characteristics is represented reasonably well and the reduction of momentum flux is in good agreement with theoretical analysis.}, language = {en} } @article{SeeligRabeMalemShinitskietal.2020, author = {Seelig, Stefan A. and Rabe, Maximilian Michael and Malem-Shinitski, Noa and Risse, Sarah and Reich, Sebastian and Engbert, Ralf}, title = {Bayesian parameter estimation for the SWIFT model of eye-movement control during reading}, series = {Journal of mathematical psychology}, volume = {95}, journal = {Journal of mathematical psychology}, publisher = {Elsevier}, address = {San Diego}, issn = {0022-2496}, doi = {10.1016/j.jmp.2019.102313}, pages = {32}, year = {2020}, abstract = {Process-oriented theories of cognition must be evaluated against time-ordered observations. Here we present a representative example for data assimilation of the SWIFT model, a dynamical model of the control of fixation positions and fixation durations during natural reading of single sentences. First, we develop and test an approximate likelihood function of the model, which is a combination of a spatial, pseudo-marginal likelihood and a temporal likelihood obtained by probability density approximation Second, we implement a Bayesian approach to parameter inference using an adaptive Markov chain Monte Carlo procedure. Our results indicate that model parameters can be estimated reliably for individual subjects. We conclude that approximative Bayesian inference represents a considerable step forward for computational models of eye-movement control, where modeling of individual data on the basis of process-based dynamic models has not been possible so far.}, language = {en} } @article{SchuettRothkegelTrukenbrodetal.2017, author = {Sch{\"u}tt, Heiko Herbert and Rothkegel, Lars Oliver Martin and Trukenbrod, Hans Arne and Reich, Sebastian and Wichmann, Felix A. and Engbert, Ralf}, title = {Likelihood-based parameter estimation and comparison of dynamical cognitive models}, series = {Psychological Review}, volume = {124}, journal = {Psychological Review}, number = {4}, publisher = {American Psychological Association}, address = {Washington}, issn = {0033-295X}, doi = {10.1037/rev0000068}, pages = {505 -- 524}, year = {2017}, abstract = {Dynamical models of cognition play an increasingly important role in driving theoretical and experimental research in psychology. Therefore, parameter estimation, model analysis and comparison of dynamical models are of essential importance. In this article, we propose a maximum likelihood approach for model analysis in a fully dynamical framework that includes time-ordered experimental data. Our methods can be applied to dynamical models for the prediction of discrete behavior (e.g., movement onsets); in particular, we use a dynamical model of saccade generation in scene viewing as a case study for our approach. For this model, the likelihood function can be computed directly by numerical simulation, which enables more efficient parameter estimation including Bayesian inference to obtain reliable estimates and corresponding credible intervals. Using hierarchical models inference is even possible for individual observers. Furthermore, our likelihood approach can be used to compare different models. In our example, the dynamical framework is shown to outperform nondynamical statistical models. Additionally, the likelihood based evaluation differentiates model variants, which produced indistinguishable predictions on hitherto used statistics. Our results indicate that the likelihood approach is a promising framework for dynamical cognitive models.}, language = {en} } @article{ReichWeissmann2021, author = {Reich, Sebastian and Weissmann, Simon}, title = {Fokker-Planck particle systems for Bayesian inference: computational approaches}, series = {SIAM ASA journal on uncertainty quantification}, volume = {9}, journal = {SIAM ASA journal on uncertainty quantification}, number = {2}, publisher = {Society for Industrial and Applied Mathematics}, address = {Philadelphia}, issn = {2166-2525}, doi = {10.1137/19M1303162}, pages = {446 -- 482}, year = {2021}, abstract = {Bayesian inference can be embedded into an appropriately defined dynamics in the space of probability measures. In this paper, we take Brownian motion and its associated Fokker-Planck equation as a starting point for such embeddings and explore several interacting particle approximations. More specifically, we consider both deterministic and stochastic interacting particle systems and combine them with the idea of preconditioning by the empirical covariance matrix. In addition to leading to affine invariant formulations which asymptotically speed up convergence, preconditioning allows for gradient-free implementations in the spirit of the ensemble Kalman filter. While such gradient-free implementations have been demonstrated to work well for posterior measures that are nearly Gaussian, we extend their scope of applicability to multimodal measures by introducing localized gradient-free approximations. Numerical results demonstrate the effectiveness of the considered methodologies.}, language = {en} } @article{Reich2006, author = {Reich, Sebastian}, title = {Linearly implicit time stepping methods for numerical weather prediction}, series = {BIT : numerical mathematics ; the leading applied mathematics journal for all computational mathematicians}, volume = {46}, journal = {BIT : numerical mathematics ; the leading applied mathematics journal for all computational mathematicians}, publisher = {Springer}, address = {Dordrecht}, issn = {0006-3835}, doi = {10.1007/s10543-006-0065-0}, pages = {607 -- 616}, year = {2006}, abstract = {The efficient time integration of the dynamic core equations for numerical weather prediction (NWP) remains a key challenge. One of the most popular methods is currently provided by implementations of the semi-implicit semi-Lagrangian (SISL) method, originally proposed by Robert (J. Meteorol. Soc. Jpn., 1982). Practical implementations of the SISL method are, however, not without certain shortcomings with regard to accuracy, conservation properties and stability. Based on recent work by Gottwald, Frank and Reich (LNCSE, Springer, 2002), Frank, Reich, Staniforth, White and Wood (Atm. Sci. Lett., 2005) and Wood, Staniforth and Reich (Atm. Sci. Lett., 2006) we propose an alternative semi-Lagrangian implementation based on a set of regularized equations and the popular Stormer-Verlet time stepping method in the context of the shallow-water equations (SWEs). Ultimately, the goal is to develop practical implementations for the 3D Euler equations that overcome some or all shortcomings of current SISL implementations.}, language = {en} } @misc{Reich1995, author = {Reich, Sebastian}, title = {Smoothed dynamics of highly oscillatory Hamiltonian systems}, url = {http://nbn-resolving.de/urn:nbn:de:kobv:517-opus-15639}, year = {1995}, abstract = {We consider the numerical treatment of Hamiltonian systems that contain a potential which grows large when the system deviates from the equilibrium value of the potential. Such systems arise, e.g., in molecular dynamics simulations and the spatial discretization of Hamiltonian partial differential equations. Since the presence of highly oscillatory terms in the solutions forces any explicit integrator to use very small step size, the numerical integration of such systems provides a challenging task. It has been suggested before to replace the strong potential by a holonomic constraint that forces the solutions to stay at the equilibrium value of the potential. This approach has, e.g., been successfully applied to the bond stretching in molecular dynamics simulations. In other cases, such as the bond-angle bending, this methods fails due to the introduced rigidity. Here we give a careful analysis of the analytical problem by means of a smoothing operator. This will lead us to the notion of the smoothed dynamics of a highly oscillatory Hamiltonian system. Based on our analysis, we suggest a new constrained formulation that maintains the flexibility of the system while at the same time suppressing the high-frequency components in the solutions and thus allowing for larger time steps. The new constrained formulation is Hamiltonian and can be discretized by the well-known SHAKE method.}, language = {en} } @misc{Reich1992, author = {Reich, Sebastian}, title = {Differential-algebraic equations and applications in circuit theory}, url = {http://nbn-resolving.de/urn:nbn:de:kobv:517-opus-46646}, year = {1992}, abstract = {Technical and physical systems, especially electronic circuits, are frequently modeled as a system of differential and nonlinear implicit equations. In the literature such systems of equations are called differentialalgebraic equations (DAEs). It turns out that the numerical and analytical properties of a DAE depend on an integer called the index of the problem. For example, the well-known BDF method of Gear can be applied, in general, to a DAE only if the index does not exceed one. In this paper we give a geometric interpretation of higherindex DAEs and indicate problems arising in connection with such DAEs by means of several examples.}, language = {en} } @article{Reich2019, author = {Reich, Sebastian}, title = {Data assimilation}, series = {Acta numerica}, volume = {28}, journal = {Acta numerica}, publisher = {Cambridge Univ. Press}, address = {New York}, issn = {0962-4929}, doi = {10.1017/S0962492919000011}, pages = {635 -- 711}, year = {2019}, abstract = {Data assimilation addresses the general problem of how to combine model-based predictions with partial and noisy observations of the process in an optimal manner. This survey focuses on sequential data assimilation techniques using probabilistic particle-based algorithms. In addition to surveying recent developments for discrete- and continuous-time data assimilation, both in terms of mathematical foundations and algorithmic implementations, we also provide a unifying framework from the perspective of coupling of measures, and Schr{\"o}dinger's boundary value problem for stochastic processes in particular.}, language = {en} } @misc{Reich1990, author = {Reich, Sebastian}, title = {On a geometrical interpretation of differential-algebraic equations}, url = {http://nbn-resolving.de/urn:nbn:de:kobv:517-opus-46683}, year = {1990}, abstract = {The subject of this paper is the relation of differential-algebraic equations (DAEs) to vector fields on manifolds. For that reason, we introduce the notion of a regular DAE as a DAE to which a vector field uniquely corresponds. Furthermore, a technique is described which yields a family of manifolds for a given DAE. This socalled family of constraint manifolds allows in turn the formulation of sufficient conditions for the regularity of a DAE. and the definition of the index of a regular DAE. We also state a method for the reduction of higher-index DAEs to lowsr-index ones that can be solved without introducing additional constants of integration. Finally, the notion of realizability of a given vector field by a regular DAE is introduced, and it is shown that any vector field can be realized by a regular DAE. Throughout this paper the problem of path-tracing is discussed as an illustration of the mathematical phenomena.}, language = {en} } @misc{Reich1994, author = {Reich, Sebastian}, title = {Momentum conserving symplectic integrators}, url = {http://nbn-resolving.de/urn:nbn:de:kobv:517-opus-16824}, year = {1994}, abstract = {In this paper, we show that symplectic partitioned Runge-Kutta methods conserve momentum maps corresponding to linear symmetry groups acting on the phase space of Hamiltonian differential equations by extended point transformation. We also generalize this result to constrained systems and show how this conservation property relates to the symplectic integration of Lie-Poisson systems on certain submanifolds of the general matrix group GL(n).}, language = {en} } @misc{Reich1995, author = {Reich, Sebastian}, title = {On the local qualitative behavior of differential-algebraic equations}, url = {http://nbn-resolving.de/urn:nbn:de:kobv:517-opus-46739}, year = {1995}, abstract = {A theoretical famework for the investigation of the qualitative behavior of differential-algebraic equations (DAEs) near an equilibrium point is established. The key notion of our approach is the notion of regularity. A DAE is called regular locally around an equilibrium point if there is a unique vector field such that the solutions of the DAE and the vector field are in one-to-one correspondence in a neighborhood of this equili Drium point. Sufficient conditions for the regularity of an equilibrium point are stated. This in turn allows us to translate several local results, as formulated for vector fields, to DAEs that are regular locally around a g: ven equilibrium point (e.g. Local Stable and Unstable Manifold Theorem, Hopf theorem). It is important that ihese theorems are stated in terms of the given problem and not in terms of the corresponding vector field.}, language = {en} }