Florent Bréhard, Nicolas Brisebarre, Mioara Joldes, Warwick Tucker
Abelian integrals play a key role in the infinitesimal version of Hilbert’s 16th problem. Being able to evaluate such integrals – with guaranteed error bounds – is a fundamental step in computer-aided proofs aimed at this problem. Using interpolation by trigonometric polynomials and quasi-Newton-Kantorovitch validation, we develop a validated numerics method for computing Abelian integrals in a quasi-linear number of arithmetic operations. Our approach is both effective, as exemplified on two practical perturbed integrable systems, and amenable to an implementation in a formal proof assistant, which is key to provide fully reliable computer-aided proofs.
{"title":"Efficient and Validated Numerical Evaluation of Abelian Integrals","authors":"Florent Bréhard, Nicolas Brisebarre, Mioara Joldes, Warwick Tucker","doi":"10.1145/3637550","DOIUrl":"https://doi.org/10.1145/3637550","url":null,"abstract":"<p>Abelian integrals play a key role in the infinitesimal version of Hilbert’s 16th problem. Being able to evaluate such integrals – with guaranteed error bounds – is a fundamental step in computer-aided proofs aimed at this problem. Using interpolation by trigonometric polynomials and quasi-Newton-Kantorovitch validation, we develop a validated numerics method for computing Abelian integrals in a quasi-linear number of arithmetic operations. Our approach is both effective, as exemplified on two practical perturbed integrable systems, and amenable to an implementation in a formal proof assistant, which is key to provide fully reliable computer-aided proofs.</p>","PeriodicalId":50935,"journal":{"name":"ACM Transactions on Mathematical Software","volume":"80 1","pages":""},"PeriodicalIF":2.7,"publicationDate":"2023-12-18","publicationTypes":"Journal Article","fieldsOfStudy":null,"isOpenAccess":false,"openAccessPdf":"","citationCount":null,"resultStr":null,"platform":"Semanticscholar","paperid":"138714423","PeriodicalName":null,"FirstCategoryId":null,"ListUrlMain":null,"RegionNum":1,"RegionCategory":"数学","ArticlePicture":[],"TitleCN":null,"AbstractTextCN":null,"PMCID":"","EPubDate":null,"PubModel":null,"JCR":null,"JCRName":null,"Score":null,"Total":0}
Data-flow reversal is at the heart of source-transformation reverse algorithmic differentiation (reverse ST-AD), arguably the most efficient way to obtain gradients of numerical models. However, when the model implementation language uses garbage collection (GC), for instance in Java or Python, the notion of address that is needed for data-flow reversal disappears. Moreover, GC is asynchronous and does not appear explicitly in the source. This paper presents an extension to the model of reverse ST-AD suitable for a language with GC. The approach is validated on a Java implementation of a simple Navier-Stokes solver. Performance is compared with existing AD tools ADOL-C and Tapenade on an equivalent implementation in C and Fortran.
{"title":"Data-flow Reversal and Garbage Collection","authors":"Laurent Hascoët","doi":"10.1145/3627537","DOIUrl":"https://doi.org/10.1145/3627537","url":null,"abstract":"<p>Data-flow reversal is at the heart of source-transformation reverse algorithmic differentiation (reverse ST-AD), arguably the most efficient way to obtain gradients of numerical models. However, when the model implementation language uses garbage collection (GC), for instance in Java or Python, the notion of address that is needed for data-flow reversal disappears. Moreover, GC is asynchronous and does not appear explicitly in the source. This paper presents an extension to the model of reverse ST-AD suitable for a language with GC. The approach is validated on a Java implementation of a simple Navier-Stokes solver. Performance is compared with existing AD tools ADOL-C and Tapenade on an equivalent implementation in C and Fortran.</p>","PeriodicalId":50935,"journal":{"name":"ACM Transactions on Mathematical Software","volume":"203 1","pages":""},"PeriodicalIF":2.7,"publicationDate":"2023-11-18","publicationTypes":"Journal Article","fieldsOfStudy":null,"isOpenAccess":false,"openAccessPdf":"","citationCount":null,"resultStr":null,"platform":"Semanticscholar","paperid":"138537824","PeriodicalName":null,"FirstCategoryId":null,"ListUrlMain":null,"RegionNum":1,"RegionCategory":"数学","ArticlePicture":[],"TitleCN":null,"AbstractTextCN":null,"PMCID":"","EPubDate":null,"PubModel":null,"JCR":null,"JCRName":null,"Score":null,"Total":0}
Timbwoga A. J. Ouermi, Robert M. Kirby, Martin Berzins
Polynomial interpolation is an important component of many computational problems. In several of these computational problems, failure to preserve positivity when using polynomials to approximate or map data values between meshes can lead to negative unphysical quantities. Currently, most polynomial-based methods for enforcing positivity are based on splines and polynomial rescaling. The spline-based approaches build interpolants that are positive over the intervals in which they are defined and may require solving a minimization problem and/or system of equations. The linear polynomial rescaling methods allow for high-degree polynomials but enforce positivity only at limited locations (e.g., quadrature nodes). This work introduces open-source software (HiPPIS) for high-order data-bounded interpolation (DBI) and positivity-preserving interpolation (PPI) that addresses the limitations of both the spline and polynomial rescaling methods. HiPPIS is suitable for approximating and mapping physical quantities such as mass, density, and concentration between meshes while preserving positivity. This work provides Fortran and Matlab implementations of the DBI and PPI methods, presents an analysis of the mapping error in the context of PDEs, and uses several 1D and 2D numerical examples to demonstrate the benefits and limitations of HiPPIS.
{"title":"Algorithm xxxx: HiPPIS A High-Order Positivity-Preserving Mapping Software for Structured Meshes","authors":"Timbwoga A. J. Ouermi, Robert M. Kirby, Martin Berzins","doi":"10.1145/3632291","DOIUrl":"https://doi.org/10.1145/3632291","url":null,"abstract":"Polynomial interpolation is an important component of many computational problems. In several of these computational problems, failure to preserve positivity when using polynomials to approximate or map data values between meshes can lead to negative unphysical quantities. Currently, most polynomial-based methods for enforcing positivity are based on splines and polynomial rescaling. The spline-based approaches build interpolants that are positive over the intervals in which they are defined and may require solving a minimization problem and/or system of equations. The linear polynomial rescaling methods allow for high-degree polynomials but enforce positivity only at limited locations (e.g., quadrature nodes). This work introduces open-source software (HiPPIS) for high-order data-bounded interpolation (DBI) and positivity-preserving interpolation (PPI) that addresses the limitations of both the spline and polynomial rescaling methods. HiPPIS is suitable for approximating and mapping physical quantities such as mass, density, and concentration between meshes while preserving positivity. This work provides Fortran and Matlab implementations of the DBI and PPI methods, presents an analysis of the mapping error in the context of PDEs, and uses several 1D and 2D numerical examples to demonstrate the benefits and limitations of HiPPIS.","PeriodicalId":50935,"journal":{"name":"ACM Transactions on Mathematical Software","volume":"103 17","pages":"0"},"PeriodicalIF":0.0,"publicationDate":"2023-11-10","publicationTypes":"Journal Article","fieldsOfStudy":null,"isOpenAccess":false,"openAccessPdf":"","citationCount":null,"resultStr":null,"platform":"Semanticscholar","paperid":"135136831","PeriodicalName":null,"FirstCategoryId":null,"ListUrlMain":null,"RegionNum":1,"RegionCategory":"数学","ArticlePicture":[],"TitleCN":null,"AbstractTextCN":null,"PMCID":"","EPubDate":null,"PubModel":null,"JCR":null,"JCRName":null,"Score":null,"Total":0}
The Sparse Grids Matlab Kit provides a Matlab implementation of sparse grids, and can be used for approximating high-dimensional functions and, in particular, for surrogate-model-based uncertainty quantification. It is lightweight, high-level and easy to use, good for quick prototyping and teaching; however, it is equipped with some features that allow its use also in realistic applications. The goal of this paper is to provide an overview of the data structure and of the mathematical aspects forming the basis of the software, as well as comparing the current release of our package to similar available software.
{"title":"Algorithm X: The Sparse Grids Matlab Kit - a Matlab implementation of sparse grids for high-dimensional function approximation and uncertainty quantification","authors":"Chiara Piazzola, Lorenzo Tamellini","doi":"10.1145/3630023","DOIUrl":"https://doi.org/10.1145/3630023","url":null,"abstract":"The Sparse Grids Matlab Kit provides a Matlab implementation of sparse grids, and can be used for approximating high-dimensional functions and, in particular, for surrogate-model-based uncertainty quantification. It is lightweight, high-level and easy to use, good for quick prototyping and teaching; however, it is equipped with some features that allow its use also in realistic applications. The goal of this paper is to provide an overview of the data structure and of the mathematical aspects forming the basis of the software, as well as comparing the current release of our package to similar available software.","PeriodicalId":50935,"journal":{"name":"ACM Transactions on Mathematical Software","volume":"40 9","pages":"0"},"PeriodicalIF":0.0,"publicationDate":"2023-11-03","publicationTypes":"Journal Article","fieldsOfStudy":null,"isOpenAccess":false,"openAccessPdf":"","citationCount":null,"resultStr":null,"platform":"Semanticscholar","paperid":"135818721","PeriodicalName":null,"FirstCategoryId":null,"ListUrlMain":null,"RegionNum":1,"RegionCategory":"数学","ArticlePicture":[],"TitleCN":null,"AbstractTextCN":null,"PMCID":"","EPubDate":null,"PubModel":null,"JCR":null,"JCRName":null,"Score":null,"Total":0}
Jonas Kusch, Steffen Schotthöfer, Pia Stammer, Jannick Wolters, Tianbai Xiao
In this paper, we present KiT-RT (Kinetic Transport Solver for Radiation Therapy), an open-source C++ based framework for solving kinetic equations in therapy applications available at https://github.com/CSMMLab/KiT-RT. This software framework aims to provide a collection of classical deterministic solvers for unstructured meshes that allow for easy extendability. Therefore, KiT-RT is a convenient base to test new numerical methods in various applications and compare them against conventional solvers. The implementation includes spherical harmonics, minimal entropy, neural minimal entropy, and discrete ordinates methods. Solution characteristics and efficiency are presented through several test cases ranging from radiation transport to electron radiation therapy. Due to the variety of included numerical methods and easy extendability, the presented open-source code is attractive for both developers, who want a basis to build their numerical solvers, and users or application engineers, who want to gain experimental insights without directly interfering with the codebase.
{"title":"KiT-RT: An extendable framework for radiative transfer and therapy","authors":"Jonas Kusch, Steffen Schotthöfer, Pia Stammer, Jannick Wolters, Tianbai Xiao","doi":"10.1145/3630001","DOIUrl":"https://doi.org/10.1145/3630001","url":null,"abstract":"In this paper, we present KiT-RT (Kinetic Transport Solver for Radiation Therapy), an open-source C++ based framework for solving kinetic equations in therapy applications available at https://github.com/CSMMLab/KiT-RT. This software framework aims to provide a collection of classical deterministic solvers for unstructured meshes that allow for easy extendability. Therefore, KiT-RT is a convenient base to test new numerical methods in various applications and compare them against conventional solvers. The implementation includes spherical harmonics, minimal entropy, neural minimal entropy, and discrete ordinates methods. Solution characteristics and efficiency are presented through several test cases ranging from radiation transport to electron radiation therapy. Due to the variety of included numerical methods and easy extendability, the presented open-source code is attractive for both developers, who want a basis to build their numerical solvers, and users or application engineers, who want to gain experimental insights without directly interfering with the codebase.","PeriodicalId":50935,"journal":{"name":"ACM Transactions on Mathematical Software","volume":"8 14","pages":"0"},"PeriodicalIF":0.0,"publicationDate":"2023-10-27","publicationTypes":"Journal Article","fieldsOfStudy":null,"isOpenAccess":false,"openAccessPdf":"","citationCount":null,"resultStr":null,"platform":"Semanticscholar","paperid":"136234105","PeriodicalName":null,"FirstCategoryId":null,"ListUrlMain":null,"RegionNum":1,"RegionCategory":"数学","ArticlePicture":[],"TitleCN":null,"AbstractTextCN":null,"PMCID":"","EPubDate":null,"PubModel":null,"JCR":null,"JCRName":null,"Score":null,"Total":0}
An interface preserving moving mesh algorithm in two or higher dimensions is presented. It resolves a moving ( d − 1)-dimensional manifold directly within the d -dimensional mesh, which means that the interface is represented by a subset of moving mesh cell-surfaces. The underlying mesh is a conforming simplicial partition that fulfills the Delaunay property. The local remeshing algorithms allow for strong interface deformations. We give a proof that the given algorithms preserve the interface after interface deformation and remeshing steps. Originating from various numerical methods, data is attached cell-wise to the mesh. After each remeshing operation the interface preserving moving mesh retains valid data by projecting the data to the new mesh cells. An open source implementation of the moving mesh algorithm is available at [1].
{"title":"An Interface Preserving Moving Mesh in Multiple Space Dimensions","authors":"Maria Alkämper, Jim Magiera, Christian Rohde","doi":"10.1145/3630000","DOIUrl":"https://doi.org/10.1145/3630000","url":null,"abstract":"An interface preserving moving mesh algorithm in two or higher dimensions is presented. It resolves a moving ( d − 1)-dimensional manifold directly within the d -dimensional mesh, which means that the interface is represented by a subset of moving mesh cell-surfaces. The underlying mesh is a conforming simplicial partition that fulfills the Delaunay property. The local remeshing algorithms allow for strong interface deformations. We give a proof that the given algorithms preserve the interface after interface deformation and remeshing steps. Originating from various numerical methods, data is attached cell-wise to the mesh. After each remeshing operation the interface preserving moving mesh retains valid data by projecting the data to the new mesh cells. An open source implementation of the moving mesh algorithm is available at [1].","PeriodicalId":50935,"journal":{"name":"ACM Transactions on Mathematical Software","volume":"42 6","pages":"0"},"PeriodicalIF":0.0,"publicationDate":"2023-10-26","publicationTypes":"Journal Article","fieldsOfStudy":null,"isOpenAccess":false,"openAccessPdf":"","citationCount":null,"resultStr":null,"platform":"Semanticscholar","paperid":"134908060","PeriodicalName":null,"FirstCategoryId":null,"ListUrlMain":null,"RegionNum":1,"RegionCategory":"数学","ArticlePicture":[],"TitleCN":null,"AbstractTextCN":null,"PMCID":"","EPubDate":null,"PubModel":null,"JCR":null,"JCRName":null,"Score":null,"Total":0}
Ana Budiša, Xiaozhe Hu, Miroslav Kuchta, Kent–André Mardal, Ludmil T. Zikatanov
We introduce the software toolbox HAZniCS for solving interface-coupled multiphysics problems. HAZniCS is a suite of modules that combines the well-known FEniCS framework for finite element discretization with solver and graph library HAZmath. The focus of the paper is on the design and implementation of robust and efficient solver algorithms which tackle issues related to the complex interfacial coupling of the physical problems often encountered in applications in brain biomechanics. The robustness and efficiency of the numerical algorithms and methods is shown in several numerical examples, namely the Darcy-Stokes equations that model flow of cerebrospinal fluid in the human brain and the mixed-dimensional model of electrodiffusion in the brain tissue.
{"title":"HAZniCS – Software Components for Multiphysics Problems","authors":"Ana Budiša, Xiaozhe Hu, Miroslav Kuchta, Kent–André Mardal, Ludmil T. Zikatanov","doi":"10.1145/3625561","DOIUrl":"https://doi.org/10.1145/3625561","url":null,"abstract":"We introduce the software toolbox HAZniCS for solving interface-coupled multiphysics problems. HAZniCS is a suite of modules that combines the well-known FEniCS framework for finite element discretization with solver and graph library HAZmath. The focus of the paper is on the design and implementation of robust and efficient solver algorithms which tackle issues related to the complex interfacial coupling of the physical problems often encountered in applications in brain biomechanics. The robustness and efficiency of the numerical algorithms and methods is shown in several numerical examples, namely the Darcy-Stokes equations that model flow of cerebrospinal fluid in the human brain and the mixed-dimensional model of electrodiffusion in the brain tissue.","PeriodicalId":50935,"journal":{"name":"ACM Transactions on Mathematical Software","volume":"210 1","pages":"0"},"PeriodicalIF":0.0,"publicationDate":"2023-10-17","publicationTypes":"Journal Article","fieldsOfStudy":null,"isOpenAccess":false,"openAccessPdf":"","citationCount":null,"resultStr":null,"platform":"Semanticscholar","paperid":"136034148","PeriodicalName":null,"FirstCategoryId":null,"ListUrlMain":null,"RegionNum":1,"RegionCategory":"数学","ArticlePicture":[],"TitleCN":null,"AbstractTextCN":null,"PMCID":"","EPubDate":null,"PubModel":null,"JCR":null,"JCRName":null,"Score":null,"Total":0}
Alice Le Brigant, Jules Deschamps, Antoine Collas, Nina Miolane
We introduce the information geometry module of the Python package Geomstats. The module first implements Fisher-Rao Riemannian manifolds of widely used parametric families of probability distributions, such as normal, gamma, beta, Dirichlet distributions, and more. The module further gives the Fisher-Rao Riemannian geometry of any parametric family of distributions of interest, given a parameterized probability density function as input. The implemented Riemannian geometry tools allow users to compare, average, interpolate between distributions inside a given family. Importantly, such capabilities open the door to statistics and machine learning on probability distributions. We present the object-oriented implementation of the module along with illustrative examples and show how it can be used to perform learning on manifolds of parametric probability distributions.
{"title":"Parametric information geometry with the package Geomstats","authors":"Alice Le Brigant, Jules Deschamps, Antoine Collas, Nina Miolane","doi":"10.1145/3627538","DOIUrl":"https://doi.org/10.1145/3627538","url":null,"abstract":"We introduce the information geometry module of the Python package Geomstats. The module first implements Fisher-Rao Riemannian manifolds of widely used parametric families of probability distributions, such as normal, gamma, beta, Dirichlet distributions, and more. The module further gives the Fisher-Rao Riemannian geometry of any parametric family of distributions of interest, given a parameterized probability density function as input. The implemented Riemannian geometry tools allow users to compare, average, interpolate between distributions inside a given family. Importantly, such capabilities open the door to statistics and machine learning on probability distributions. We present the object-oriented implementation of the module along with illustrative examples and show how it can be used to perform learning on manifolds of parametric probability distributions.","PeriodicalId":50935,"journal":{"name":"ACM Transactions on Mathematical Software","volume":"74 1","pages":"0"},"PeriodicalIF":0.0,"publicationDate":"2023-10-13","publicationTypes":"Journal Article","fieldsOfStudy":null,"isOpenAccess":false,"openAccessPdf":"","citationCount":null,"resultStr":null,"platform":"Semanticscholar","paperid":"135805180","PeriodicalName":null,"FirstCategoryId":null,"ListUrlMain":null,"RegionNum":1,"RegionCategory":"数学","ArticlePicture":[],"TitleCN":null,"AbstractTextCN":null,"PMCID":"","EPubDate":null,"PubModel":null,"JCR":null,"JCRName":null,"Score":null,"Total":0}
General expressions are derived for the amplitude equation valid at a Turing bifurcation of a system of reaction-diffusion equations in one spatial dimension, with an arbitrary number of components. The normal form is computed up to fifth order, which enables the detection and analysis of codimension-two points where the criticality of the bifurcation changes. The expressions are implemented within a Python package, in which the user needs to specify only expressions for the reaction kinetics and the values of diffusion constants. The code is augmented with a Mathematica routine to compute curves of Turing bifurcations in a parameter plane and automatically detect codimension-two points. The software is illustrated with examples that show the versatility of the method including a case with cross-diffusion, a higher-order scalar equation and a four-component system.
{"title":"Computation of Turing bifurcation normal form for <i>n</i> -component reaction-diffusion systems.","authors":"Edgardo Villar-Sepúlveda, Alan Champneys","doi":"10.1145/3625560","DOIUrl":"https://doi.org/10.1145/3625560","url":null,"abstract":"General expressions are derived for the amplitude equation valid at a Turing bifurcation of a system of reaction-diffusion equations in one spatial dimension, with an arbitrary number of components. The normal form is computed up to fifth order, which enables the detection and analysis of codimension-two points where the criticality of the bifurcation changes. The expressions are implemented within a Python package, in which the user needs to specify only expressions for the reaction kinetics and the values of diffusion constants. The code is augmented with a Mathematica routine to compute curves of Turing bifurcations in a parameter plane and automatically detect codimension-two points. The software is illustrated with examples that show the versatility of the method including a case with cross-diffusion, a higher-order scalar equation and a four-component system.","PeriodicalId":50935,"journal":{"name":"ACM Transactions on Mathematical Software","volume":"57 1","pages":"0"},"PeriodicalIF":0.0,"publicationDate":"2023-09-29","publicationTypes":"Journal Article","fieldsOfStudy":null,"isOpenAccess":false,"openAccessPdf":"","citationCount":null,"resultStr":null,"platform":"Semanticscholar","paperid":"135246930","PeriodicalName":null,"FirstCategoryId":null,"ListUrlMain":null,"RegionNum":1,"RegionCategory":"数学","ArticlePicture":[],"TitleCN":null,"AbstractTextCN":null,"PMCID":"","EPubDate":null,"PubModel":null,"JCR":null,"JCRName":null,"Score":null,"Total":0}
Hendrik Ranocha, Michael Schlottke-Lakemper, Jesse Chan, Andrés M. Rueda-Ramírez, Andrew R. Winters, Florian Hindenlang, Gregor J. Gassner
Many modern discontinuous Galerkin (DG) methods for conservation laws make use of summation by parts operators and flux differencing to achieve kinetic energy preservation or entropy stability. While these techniques increase the robustness of DG methods significantly, they are also computationally more demanding than standard weak form nodal DG methods. We present several implementation techniques to improve the efficiency of flux differencing DG methods that use tensor product quadrilateral or hexahedral elements, in 2D or 3D respectively. Focus is mostly given to CPUs and DG methods for the compressible Euler equations, although these techniques are generally also useful for other physical systems including the compressible Navier-Stokes and magnetohydrodynamics equations. We present results using two open source codes, Trixi.jl written in Julia and FLUXO written in Fortran, to demonstrate that our proposed implementation techniques are applicable to different code bases and programming languages.
{"title":"Efficient implementation of modern entropy stable and kinetic energy preserving discontinuous Galerkin methods for conservation laws","authors":"Hendrik Ranocha, Michael Schlottke-Lakemper, Jesse Chan, Andrés M. Rueda-Ramírez, Andrew R. Winters, Florian Hindenlang, Gregor J. Gassner","doi":"10.1145/3625559","DOIUrl":"https://doi.org/10.1145/3625559","url":null,"abstract":"Many modern discontinuous Galerkin (DG) methods for conservation laws make use of summation by parts operators and flux differencing to achieve kinetic energy preservation or entropy stability. While these techniques increase the robustness of DG methods significantly, they are also computationally more demanding than standard weak form nodal DG methods. We present several implementation techniques to improve the efficiency of flux differencing DG methods that use tensor product quadrilateral or hexahedral elements, in 2D or 3D respectively. Focus is mostly given to CPUs and DG methods for the compressible Euler equations, although these techniques are generally also useful for other physical systems including the compressible Navier-Stokes and magnetohydrodynamics equations. We present results using two open source codes, Trixi.jl written in Julia and FLUXO written in Fortran, to demonstrate that our proposed implementation techniques are applicable to different code bases and programming languages.","PeriodicalId":50935,"journal":{"name":"ACM Transactions on Mathematical Software","volume":"23 1","pages":"0"},"PeriodicalIF":0.0,"publicationDate":"2023-09-27","publicationTypes":"Journal Article","fieldsOfStudy":null,"isOpenAccess":false,"openAccessPdf":"","citationCount":null,"resultStr":null,"platform":"Semanticscholar","paperid":"135536808","PeriodicalName":null,"FirstCategoryId":null,"ListUrlMain":null,"RegionNum":1,"RegionCategory":"数学","ArticlePicture":[],"TitleCN":null,"AbstractTextCN":null,"PMCID":"","EPubDate":null,"PubModel":null,"JCR":null,"JCRName":null,"Score":null,"Total":0}