New paper published: Reliable genotyping of recombinant genomes using a robust hidden Markov model

Meiotic recombination is an essential mechanism during sexual reproduction and includes the exchange of chromosome segments between homologous chromosomes. New allelic combinations are transmitted to the new generation, introducing novel genetic variation in the offspring genomes. With the improvement of high-throughput whole-genome sequencing technologies, large numbers of recombinant individuals can now be sequenced with low sequencing depth at low costs, necessitating computational methods for reconstructing their haplotypes. The main challenge is the uncertainty in haplotype calling that arises from the low information content of a single genomic position. Straightforward sliding window-based approaches are difficult to tune and fail to place recombination breakpoints precisely. Hidden Markov Model-based approaches, on the other hand, tend to over-segment the genome. Here, we present RTIGER, an HMM-based model that exploits in a mathematically precise way the fact that true chromosome segments typically have a certain minimum length. We further separate the task of identifying the correct haplotype sequence from the accurate placement of haplotype borders, thereby maximizing the accuracy of border positions. By comparing segmentations based on simulated data with known underlying haplotypes, we highlight the reasons for RTIGER outperforming traditional segmentation approaches. We then analyze the meiotic recombination pattern of segregants of two Arabidopsis (Arabidopsis thaliana) accessions and a previously described hyper-recombining mutant. RTIGER is available as an R-package with an efficient Julia implementation of the core algorithm.

Plants Physiology: https://doi.org/10.1093/plphys/kiad191

New paper submitted: On the Theoretical Foundation of Overset Grid Methods for Hyperbolic Problems II: Entropy Bounded Formulations for Nonlinear Conservation Laws

We derive entropy conserving and entropy dissipative overlapping domain formulations for systems of nonlinear hyperbolic equations in conservation form, such as would be approximated by overset mesh methods. The entropy conserving formulation imposes two-way coupling at the artificial interface boundaries through nonlinear penalty functions that vanish when the solutions coincide. The penalty functions are expressed in terms of entropy conserving fluxes originally introduced for finite volume schemes. Entropy dissipation and additional coupling in the overlap region are added through the use of linear penalties.

arXiv: https://arxiv.org/abs/2203.11149

New paper submitted: On the entropy projection and the robustness of high order entropy stable discontinuous Galerkin schemes for under-resolved flows

High order entropy stable schemes provide improved robustness for computational simulations of fluid flows. However, additional stabilization and positivity preserving limiting can still be required for variable-density flows with under-resolved features. We demonstrate numerically that entropy stable DG methods which incorporate an “entropy projection” are less likely to require additional limiting to retain positivity for certain types of flows. We conclude by investigating potential explanations for this observed improvement in robustness.

arXiv: https://arxiv.org/abs/2203.10238

New paper submitted: Entropy-Stable Gauss Collocation Methods for Ideal Magneto-Hydrodynamics

In this paper, we present an entropy-stable Gauss collocation discontinuous Galerkin (DG) method on 3D curvilinear meshes for the GLM-MHD equations: the single-fluid magneto-hydrodynamics (MHD) equations with a generalized Lagrange multiplier (GLM) divergence cleaning mechanism. For the continuous entropy analysis to hold and to ensure Galilean invariance in the divergence cleaning technique, the GLM-MHD system requires the use of non-conservative terms.
Traditionally, entropy-stable DG discretizations have used a collocated nodal variant of the DG method, also known as the discontinuous Galerkin spectral element method (DGSEM) on Legendre-Gauss-Lobatto (LGL) points. Recently, Chan et al. (“Efficient Entropy Stable Gauss Collocation Methods”. SIAM -2019) presented an entropy-stable DGSEM scheme that uses Legendre-Gauss points (instead of LGL points) for conservation laws. Our main contribution is to extend the discretization technique of Chan et al. to the non-conservative GLM-MHD system.
We provide a numerical verification of the entropy behavior and convergence properties of our novel scheme on 3D curvilinear meshes. Moreover, we test the robustness and accuracy of our scheme with a magneto-hydrodynamic Kelvin-Helmholtz instability problem. The numerical experiments suggest that the entropy-stable DGSEM on Gauss points for the GLM-MHD system is more accurate than the LGL counterpart.

arXiv:2203.06062

New paper submitted: Subcell limiting strategies for discontinuous Galerkin spectral element methods

We present a general family of subcell limiting strategies to construct robust high-order accurate nodal discontinuous Galerkin (DG) schemes. The main strategy is to construct compatible low order finite volume (FV) type discretizations that allow for convex blending with the high-order variant with the goal of guaranteeing additional properties, such as bounds on physical quantities and/or guaranteed entropy dissipation. For an implementation of this main strategy, four main ingredients are identified that may be combined in a flexible manner: (i) a nodal high-order DG method on Legendre-Gauss-Lobatto nodes, (ii) a compatible robust subcell FV scheme, (iii) a convex combination strategy for the two schemes, which can be element-wise or subcell-wise, and (iv) a strategy to compute the convex blending factors, which can be either based on heuristic troubled-cell indicators, or using ideas from flux-corrected transport methods.
By carefully designing the metric terms of the subcell FV method, the resulting methods can be used on unstructured curvilinear meshes, are locally conservative, can handle strong shocks efficiently while directly guaranteeing physical bounds on quantities such as density, pressure or entropy. We further show that it is possible to choose the four ingredients to recover existing methods such as a provably entropy dissipative subcell shock-capturing approach or a sparse invariant domain preserving approach.
We test the versatility of the presented strategies and mix and match the four ingredients to solve challenging simulation setups, such as the KPP problem (a hyperbolic conservation law with non-convex flux function), turbulent and hypersonic Euler simulations, and MHD problems featuring shocks and turbulence.

arXiv:2202.00576

New paper submitted: A Discontinuous Galerkin Solver in the FLASH Multi-Physics Framework

In this paper, we present a discontinuous Galerkin solver based on previous work by Markert et al. (2021) for magneto-hydrodynamics in form of a new fluid solver module integrated into the established and well-known multi-physics simulation code FLASH. Our goal is to enable future research on the capabilities and potential advantages of discontinuous Galerkin methods for complex multi-physics simulations in astrophysical settings. We give specific details and adjustments of our implementation within the FLASH framework and present extensive validations and test cases, specifically its interaction with several other physics modules such as (self-)gravity and radiative transfer. We conclude that the new DG solver module in FLASH is ready for use in astrophysics simulations and thus ready for assessments and investigations.

arXiv:2112.11318

New paper submitted: Efficient implementation of modern entropy stable and kinetic energy preserving discontinuous Galerkin methods for conservation laws

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 GPU computing and 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.

arXiv:2112.10517

New paper published: On the theoretical foundation of overset grid methods for hyperbolic problems: Well-posedness and conservation

We use the energy method to study the well-posedness of initial-boundary value problems approximated by overset mesh methods in one and two space dimensions for linear constant-coefficient hyperbolic systems. We show that in one space dimension, for both scalar equations and systems of equations, the problem where one domain partially oversets another is well-posed when characteristic coupling conditions are used. If a system cannot be diagonalized, as is usually the case in multiple space dimensions, then the energy method does not give proper bounds in terms of initial and boundary data. For those problems, we propose a novel penalty approach. We show, by using a global energy that accounts for the energy in the overlap region of the domains, that under well-defined conditions on the coupling matrices the penalized overset domain problems are energy bounded, conservative, well-posed and have solutions equivalent to the original single domain problem.

Published in Journal of Computational Physics (ScienceDirect)

New paper published: An entropy stable nodal discontinuous Galerkin method for the resistive MHD equations. Part II: Subcell finite volume shock capturing

The second paper of this series presents two robust entropy stable shock-capturing methods for discontinuous Galerkin spectral element (DGSEM) discretizations of the compressible magneto-hydrodynamics (MHD) equations. Specifically, we use the resistive GLM-MHD equations, which include a divergence cleaning mechanism that is based on a generalized Lagrange multiplier (GLM). For the continuous entropy analysis to hold, and due to the divergence-free constraint on the magnetic field, the GLM-MHD system requires the use of non-conservative terms, which need special treatment.

Hennemann et al. (2020) [25] recently presented an entropy stable shock-capturing strategy for DGSEM discretizations of the Euler equations that blends the DGSEM scheme with a subcell first-order finite volume (FV) method. Our first contribution is the extension of the method of Hennemann et al. to systems with non-conservative terms, such as the GLM-MHD equations. In our approach, the advective and non-conservative terms of the equations are discretized with a hybrid FV/DGSEM scheme, whereas the visco-resistive terms are discretized only with the high-order DGSEM method. We prove that the extended method is semi-discretely entropy stable on three-dimensional unstructured curvilinear meshes. Our second contribution is the derivation and analysis of a second entropy stable shock-capturing method that provides enhanced resolution by using a subcell reconstruction procedure that is carefully built to ensure entropy stability.

We provide a numerical verification of the properties of the hybrid FV/DGSEM schemes on curvilinear meshes and show their robustness and accuracy with common benchmark cases, such as the Orszag-Tang vortex and the GEM (Geospace Environmental Modeling) reconnection challenge. Finally, we simulate a space physics application: the interaction of Jupiter’s magnetic field with the plasma torus generated by the moon Io.

Published in Journal of Computational Physics (ScienceDirect)

New paper submitted: On the Theoretical Foundation of Overset Grid Methods for Hyperbolic Problems: Well-Posedness and Conservation

We use the energy method to study the well-posedness of initial-boundary value problems approximated by overset mesh methods in one and two space dimensions for linear constant-coefficient hyperbolic systems. We show that in one space dimension, for both scalar equations and systems of equations, the problem where one domain partially oversets another is well-posed when characteristic coupling conditions are used. If a system cannot be diagonalized, as is ususally the case in multiple space dimensions, then the energy method does not give proper bounds in terms of initial and boundary data. For those problems, we propose a novel penalty approach. We show, by using a global energy that accounts for the energy in the overlap region of the domains, that under well-defined conditions on the coupling matrices the penalized overset domain problems are energy bounded, conservative, well-posed and have solutions equivalent to the original single domain problem.

Preprint available at: https://arxiv.org/abs/2105.04664