NySALT: Nystr\"{o}m-type inference-based schemes adaptive to large time-stepping

Large time-stepping is important for efficient long-time simulations of deterministic and stochastic Hamiltonian dynamical systems. Conventional structure-preserving integrators, while being successful for generic systems, have limited tolerance to time step size due to stability and accuracy constraints. We propose to use data to innovate classical integrators so that they can be adaptive to large time-stepping and are tailored to each specific system. In particular, we introduce NySALT, Nystr\"{o}m-type inference-based schemes adaptive to large time-stepping. The NySALT has optimal parameters for each time step learnt from data by minimizing the one-step prediction error. Thus, it is tailored for each time step size and the specific system to achieve optimal performance and tolerate large time-stepping in an adaptive fashion. We prove and numerically verify the convergence of the estimators as data size increases. Furthermore, analysis and numerical tests on the deterministic and stochastic Fermi-Pasta-Ulam (FPU) models show that NySALT enlarges the maximal admissible step size of linear stability, and quadruples the time step size of the St\"{o}rmer--Verlet and the BAOAB when maintaining similar levels of accuracy.Comment: 26 pages, 7 figure

Exploring efficient: numerical methods for differential equations

Numerical analysis is a way to do higher mathematical problems on a computer, a technique widely used by scientists and engineers to solve their problems. A major advantage of numerical analysis is that a numerical answer can be obtained even when a problem has no “analytical” solution. Results from numerical analysis are an approximation, which can be made as accurate as desired. The analysis of errors in numerical methods is a critically important part of the study of numerical analysis. Hence, we will see in this research that computation of the error is a must as it is a way to measure the efficiency of the numerical methods developed. Numerical methods require highly tedious and repetitive computations that can only be done using the computer. Hence in this research, it is shown that computer programs must be written for the implementation of numerical methods. In the early part of related research the computer language used was Fortran. Subsequently more and more computer programs used the C programming language. Additionally, now computations can also be carried out using softwares like MATLAB, MATHEMATICA and MAPLE. Many physical problems that arise from ordinary differential equations (ODEs) have magnitudes of eigenvalues which vary greatly, and such systems are commonly known as stiff systems. Stiff systems usually consist of a transient solution, that is, a solution which varies rapidly at the beginning of the integration. This phase is referred to as the transient phase and during this phase, accuracy rather than stability restricts the stepsize of the numerical methods used. Thus the generally the structure of the solutions suggests application of specific methods for non-stiff equations in the transient phase and specific methods for stiff equations during the steady-state phase in a manner whereby computational costs can be reduced. Consequently, in this research we developed embedded Runge-Kutta methods for solving stiff differential equations so that variable stepsize codes can be used in its implementation. We have also included intervalwise partitioning, whereby the system is considered as non-stiff first, and solved using the method with simple iterations, and once stiffness is detected, the system is solved using the same method, but with Newton iterations. By using variable stepsize code and intervalwise partitioning, we have been able to reduce the computational costs. With the aim of increasing the computational efficiency of the Runge-Kutta methods, we have also developed methods of higher order with less number of stages or function evaluations. The method used is an extension of the classical Runge-Kutta method and the approximation at the current point is based on the information at the current internal stage as well as the previous internal stage. This is the idea underlying the construction of Improved Runge-Kutta methods, so that the resulting method will give better accuracy. Usually higher order ordinary differential equations are solved by converting them into a system of first order ODEs and using numerical methods suitable for first order ODEs. However it is more efficient, in terms of accuracy, number of function evaluations as well as computational time, if the higher order ODEs can be solved directly (without being converted to a system of first order ODEs), using numerical methods. In this research we developed numerical methods, particularly Runge-Kutta type methods, which can directly solve special third order and fourth order ODEs. Special second order ODE is an ODE which does not depend on the first derivative. The solution from this type of ODE often exhibits a pronounced oscillatory character. It is well known that it is difficult to obtain accurate numerical results if the ODEs are oscillatory in nature. In order to address this problem a lot of research has been focused on developing methods which have high algebraic order, reduced phase-lag or dispersion and reduced dissipation. Phaselag is the angle between the true and approximate solution, while dissipation is the difference between the approximate solution and the standard cyclic solution. If a method has high algebraic order, high order of dispersion and dissipation, then the numerical solutions obtained will be very accurate. Hence in this research we have developed numerical methods, specifically hybrid methods which have all the above mentioned properties. If the solutions are oscillatory in nature, it means that the solutions will have components which are trigonometric functions, that is, sine and cosine functions. In order to get accurate numerical solutions we thus phase-fitted the methods using trigonometric functions. In this research, it is proven that trigonometrically-fitting the hybrid methods and applying them to solve oscillatory delay differential equations result in better numerical results. These are the highlights of my research journey, though a lot of work has also been done in developing numerical methods which are multistep in nature, for solving higher order ODEs, as well as implementation of methods developed for solving fuzzy differential equations and partial differential equations, which are not covered here

Connecting the Dots: Towards Continuous Time Hamiltonian Monte Carlo

Continuous time Hamiltonian Monte Carlo is introduced, as a powerful alternative to Markov chain Monte Carlo methods for continuous target distributions. The method is constructed in two steps: First Hamiltonian dynamics are chosen as the deterministic dynamics in a continuous time piecewise deterministic Markov process. Under very mild restrictions, such a process will have the desired target distribution as an invariant distribution. Secondly, the numerical implementation of such processes, based on adaptive numerical integration of second order ordinary differential equations is considered. The numerical implementation yields an approximate, yet highly robust algorithm that, unlike conventional Hamiltonian Monte Carlo, enables the exploitation of the complete Hamiltonian trajectories (hence the title). The proposed algorithm may yield large speedups and improvements in stability relative to relevant benchmarks, while incurring numerical errors that are negligible relative to the overall Monte Carlo errors

High-order, Dispersionless "Fast-Hybrid" Wave Equation Solver. Part I: O(1)\mathcal{O}(1) Sampling Cost via Incident-Field Windowing and Recentering

This paper proposes a frequency/time hybrid integral-equation method for the time dependent wave equation in two and three-dimensional spatial domains. Relying on Fourier Transformation in time, the method utilizes a fixed (time-independent) number of frequency-domain integral-equation solutions to evaluate, with superalgebraically-small errors, time domain solutions for arbitrarily long times. The approach relies on two main elements, namely, 1) A smooth time-windowing methodology that enables accurate band-limited representations for arbitrarily-long time signals, and 2) A novel Fourier transform approach which, in a time-parallel manner and without causing spurious periodicity effects, delivers numerically dispersionless spectrally-accurate solutions. A similar hybrid technique can be obtained on the basis of Laplace transforms instead of Fourier transforms, but we do not consider the Laplace-based method in the present contribution. The algorithm can handle dispersive media, it can tackle complex physical structures, it enables parallelization in time in a straightforward manner, and it allows for time leaping---that is, solution sampling at any given time $T$ at $\mathcal{O}(1)$-bounded sampling cost, for arbitrarily large values of $T$, and without requirement of evaluation of the solution at intermediate times. The proposed frequency-time hybridization strategy, which generalizes to any linear partial differential equation in the time domain for which frequency-domain solutions can be obtained (including e.g. the time-domain Maxwell equations), and which is applicable in a wide range of scientific and engineering contexts, provides significant advantages over other available alternatives such as volumetric discretization, time-domain integral equations, and convolution-quadrature approaches.Comment: 33 pages, 8 figures, revised and extended manuscript (and now including direct comparisons to existing CQ and TDIE solver implementations) (Part I of II

The numerical solution of neural field models posed on realistic cortical domains

The mathematical modelling of neural activity is a hugely complex and prominent area of exploration that has been the focus of many researchers since the mid 1900s. Although many advancements and scientific breakthroughs have been made, there is still a great deal that is not yet understood about the brain. There have been a considerable amount of studies in mathematical neuroscience that consider the brain as a simple one-dimensional or two-dimensional domain; however, this is not biologically realistic and is primarily selected as the domain of choice to aid analytical progress. The primary aim of this thesis is to develop and provide a novel suite of codes to facilitate the computationally efficient numerical solution of large-scale delay differential equations, and utilise this to explore both neural mass and neural field models with space-dependent delays. Through this, we seek to widen the scope of models of neural activity by posing them on realistic cortical domains and incorporating real brain data to describe non-local cortical connections. The suite is validated using a selection of examples that compare numerical and analytical results, along with recreating existing results from the literature. The relationship between structural connectivity and functional connectivity is then analysed as we use an eigenmode fitting approach to inform the desired stability regimes of a selection of neural mass models with delays. Here, we explore a next-generation neural mass model developed by Coombes and Byrne [36], and compare results to the more traditional Wilson-Cowan formulation [180, 181]. Finally, we examine a variety of solutions to three different neural field models that incorporate real structural connectivity, path length, and geometric surface data, using our NFESOLVE library to efficiently compute the numerical solutions. We demonstrate how the field version of the next-generation model can yield intricate and detailed solutions which push us closer to recreating observed brain dynamics

Mathematical Modeling with Differential Equations in Physics, Chemistry, Biology, and Economics

This volume was conceived as a Special Issue of the MDPI journal Mathematics to illustrate and show relevant applications of differential equations in different fields, coherently with the latest trends in applied mathematics research. All the articles that were submitted for publication are valuable, interesting, and original. The readers will certainly appreciate the heterogeneity of the 10 papers included in this book and will discover how helpful all the kinds of differential equations are in a wide range of disciplines. We are confident that this book will be inspirational for young scholars as well

Conservation of phase space properties using exponential integrators on the cubic Schrödinger equation

The cubic nonlinear Schrödinger (NLS) equation with periodic boundary conditions is solvable using Inverse Spectral Theory. The nonlinear spectrum of the associated Lax pair reveals topological properties of the NLS phase space that are difficult to assess by other means. In this paper we use the invariance of the nonlinear spectrum to examine the long time behavior of exponential and multisymplectic integrators as compared with the most commonly used split step approach. The initial condition used is a perturbation of the unstable plane wave solution, which is difficult to numerically resolve. Our findings indicate that the exponential integrators from the viewpoint of efficiency and speed have an edge over split step, while a lower order multisymplectic is not as accurate and too slow to compete. © 2006 Elsevier Inc. All rights reserved

Magnus-based geometric integrators for dynamical systems with time-dependent potentials

[ES] Esta tesis trata sobre la integración numérica de sistemas hamiltonianos con potenciales explícitamente dependientes del tiempo. Los problemas de este tipo son comunes en la física matemática, porque provienen de la mecánica cuántica, clásica y celestial. La meta de la tesis es construir integradores para unos problemas relevantes no autónomos: la ecuación de Schrödinger, que es el fundamento de la mecánica cuántica; las ecuaciones de Hill y de onda, que describen sistemas oscilatorios; el problema de Kepler con la masa variante en el tiempo. El Capítulo 1 describe la motivación y los objetivos de la obra en el contexto histórico de la integración numérica. En el Capítulo 2 se introducen los conceptos esenciales y unas herramientas fundamentales utilizadas a lo largo de la tesis. El diseño de los integradores propuestos se basa en los métodos de composición y escisión y en el desarrollo de Magnus. En el Capítulo 3 se describe el primero. Su idea principal consta de una recombinación de unos integradores sencillos para obtener la solución del problema. El concepto importante de las condiciones de orden se describe en ese capítulo. En el Capítulo 4 se hace un resumen de las álgebras de Lie y del desarrollo de Magnus que son las herramientas algebraicas que permiten expresar la solución de ecuaciones diferenciales dependientes del tiempo. La ecuación lineal de Schrödinger con potencial dependiente del tiempo está examinada en el Capítulo 5. Dado su estructura particular, nuevos métodos casi sin conmutadores, basados en el desarrollo de Magnus, son construidos. Su eficiencia es demostrada en unos experimentos numéricos con el modelo de Walker-Preston de una molécula dentro de un campo electromagnético. En el Capítulo 6, se diseñan los métodos de Magnus-escisión para las ecuaciones de onda y de Hill. Su eficiencia está demostrada en los experimentos numéricos con varios sistemas oscilatorios: con la ecuación de Mathieu, la ec. de Hill matricial, las ecuaciones de onda y de Klein-Gordon-Fock. El Capítulo 7 explica cómo el enfoque algebraico y el desarrollo de Magnus pueden generalizarse a los problemas no lineales. El ejemplo utilizado es el problema de Kepler con masa decreciente. El Capítulo 8 concluye la tesis, reseña los resultados y traza las posibles direcciones de la investigación futura.[CA] Aquesta tesi tracta de la integració numèrica de sistemes hamiltonians amb potencials explícitament dependents del temps. Els problemes d'aquest tipus són comuns en la física matemàtica, perquè provenen de la mecànica quàntica, clàssica i celest. L'objectiu de la tesi és construir integradors per a uns problemes rellevants no autònoms: l'equació de Schrödinger, que és el fonament de la mecànica quàntica; les equacions de Hill i d'ona, que descriuen sistemes oscil·latoris; el problema de Kepler amb la massa variant en el temps. El Capítol 1 descriu la motivació i els objectius de l'obra en el context històric de la integració numèrica. En Capítol 2 s'introdueixen els conceptes essencials i unes ferramentes fonamentals utilitzades al llarg de la tesi. El disseny dels integradors proposats es basa en els mètodes de composició i escissió i en el desenvolupament de Magnus. En el Capítol 3, es descriu el primer. La seua idea principal consta d'una recombinació d'uns integradors senzills per a obtenir la solució del problema. El concepte important de les condicions d'orde es descriu en eixe capítol. El Capítol 4 fa un resum de les àlgebres de Lie i del desenvolupament de Magnus que són les ferramentes algebraiques que permeten expressar la solució d'equacions diferencials dependents del temps. L'equació lineal de Schrödinger amb potencial dependent del temps està examinada en el Capítol 5. Donat la seua estructura particular, nous mètodes quasi sense commutadors, basats en el desenvolupament de Magnus, són construïts. La seua eficiència és demostrada en uns experiments numèrics amb el model de Walker-Preston d'una molècula dins d'un camp electromagnètic. En el Capítol 6 es dissenyen els mètodes de Magnus-escissió per a les equacions d'onda i de Hill. El seu rendiment està demostrat en els experiments numèrics amb diversos sistemes oscil·latoris: amb l'equació de Mathieu, l'ec. de Hill matricial, les equacions d'onda i de Klein-Gordon-Fock. El Capítol 7 explica com l'enfocament algebraic i el desenvolupament de Magnus poden generalitzar-se als problemes no lineals. L'exemple utilitzat és el problema de Kepler amb massa decreixent. El Capítol 8 conclou la tesi, ressenya els resultats i traça les possibles direccions de la investigació futura.[EN] The present thesis addresses the numerical integration of Hamiltonian systems with explicitly time-dependent potentials. These problems are common in mathematical physics because they come from quantum, classical and celestial mechanics. The goal of the thesis is to construct integrators for several import ant non-autonomous problems: the Schrödinger equation, which is the cornerstone of quantum mechanics; the Hill and the wave equations, that describe oscillating systems; the Kepler problem with time-variant mass. Chapter 1 describes the motivation and the aims of the work in the historical context of numerical integration. In Chapter 2 essential concepts and some fundamental tools used throughout the thesis are introduced. The design of the proposed integrators is based on the composition and splitting methods and the Magnus expansion. In Chapter 3, the former is described. Their main idea is to recombine some simpler integrators to obtain the solution. The salient concept of order conditions is described in that chapter. Chapter 4 summarises Lie algebras and the Magnus expansion ¿ algebraic tools that help to express the solution of time-dependent differential equations. The linear Schrödinger equation with time-dependent potential is considered in Chapter 5. Given its particular structure, new, Magnus-based quasi-commutator-free integrators are build. Their efficiency is shown in numerical experiments with the Walker-Preston model of a molecule in an electromagnetic field. In Chapter 6, Magnus-splitting methods for the wave and the Hill equations are designed. Their performance is demonstrated in numerical experiments with various oscillatory systems: the Mathieu equation, the matrix Hill eq., the wave and the Klein-Gordon-Fock eq. Chapter 7 shows how the algebraic approach and the Magnus expansion can be generalised to non-linear problems. The example used is the Kepler problem with decreasing mass. The thesis is concluded by Chapter 8, in which the results are reviewed and possible directions of future work are outlined.Kopylov, N. (2019). Magnus-based geometric integrators for dynamical systems with time-dependent potentials [Tesis doctoral no publicada]. Universitat Politècnica de València. https://doi.org/10.4995/Thesis/10251/118798TESI

The numerical solution of neural field models posed on realistic cortical domains

The mathematical modelling of neural activity is a hugely complex and prominent area of exploration that has been the focus of many researchers since the mid 1900s. Although many advancements and scientific breakthroughs have been made, there is still a great deal that is not yet understood about the brain. There have been a considerable amount of studies in mathematical neuroscience that consider the brain as a simple one-dimensional or two-dimensional domain; however, this is not biologically realistic and is primarily selected as the domain of choice to aid analytical progress. The primary aim of this thesis is to develop and provide a novel suite of codes to facilitate the computationally efficient numerical solution of large-scale delay differential equations, and utilise this to explore both neural mass and neural field models with space-dependent delays. Through this, we seek to widen the scope of models of neural activity by posing them on realistic cortical domains and incorporating real brain data to describe non-local cortical connections. The suite is validated using a selection of examples that compare numerical and analytical results, along with recreating existing results from the literature. The relationship between structural connectivity and functional connectivity is then analysed as we use an eigenmode fitting approach to inform the desired stability regimes of a selection of neural mass models with delays. Here, we explore a next-generation neural mass model developed by Coombes and Byrne [36], and compare results to the more traditional Wilson-Cowan formulation [180, 181]. Finally, we examine a variety of solutions to three different neural field models that incorporate real structural connectivity, path length, and geometric surface data, using our NFESOLVE library to efficiently compute the numerical solutions. We demonstrate how the field version of the next-generation model can yield intricate and detailed solutions which push us closer to recreating observed brain dynamics

Efficient implicit spectral/hp element DG techniques for compressible flows

In the simulation of stiff problems, such as fluid flows at high Reynolds numbers, the efficiency of explicit time integration is significantly limited by the need to use very small time steps. To alleviate this limitation and to accelerate compressible flow simulations based on high-order spectral/$hp$ element methods, an implicit time integration method is developed using singly diagonally implicit Runge-Kutta temporal discretization schemes combined with a Jacobian-free Newton Krylov (JFNK) method. This thesis studies several topics influencing the efficiency, accuracy and robustness of the solver. Firstly, an efficient partially matrix-free block relaxed Jacobi (BRJ) preconditioner is proposed, in which the Jacobian matrix and preconditioning matrices are properly approximated based on studies of their influences on convergence. The preconditioner only forms and stores the diagonal part of the Jacobian matrix while the off-diagonal operators are calculated on the fly. Used together with techniques such as using single precision data, the BRJ can largely reduce the memory consumption when compared with matrix-based ones like incomplete LU factorization preconditioners (ILU). To further accelerate the solver, influences of different parts of the flux Jacobian on the preconditioning effects are studied and terms with minor influences are neglected. This reduces the computational cost of the BRJ preconditioner by about 3 times while maintaining similar preconditioning effects. Secondly, adaptive strategies for a suitable choice of some free parameters are designed to maintain temporal accuracy and relatively high efficiency. The several free parameters in the implicit solver have significant influences on the accuracy, efficiency and stability. Therefore, designing proper strategies in choosing them is essential for developing a robust general purpose solver. Based on the idea of constructing proper relations between the temporal, spatial and iterative errors, adaptive strategies are designed for determining the time step and Newton tolerance. These parameters maintain temporal accuracy of the solver in the sense that further decreasing temporal and iterative errors will not obviously improve the efficiency. Meanwhile, they maintain relatively efficient by avoiding excessively small time step and Newton tolerance. The strategies are tested in different types of cases to illustrate their performance and generality. Finally, the implicit solver is studied in high-fidelity simulations of turbulent flows based on a hierarchical implementation in the open-source spectral/$hp$ element framework Nektar++. The solver is applied to large-eddy simulations of Taylor-Green vortex flow, turbulent channel flow and flow over a circular cylinder cases. The efficiency of the solver and the prediction accuracy of these problems are studied. The results show that the solver yields good predictions in turbulence simulations whilst keeping good stability and high efficiency.Open Acces