New preprints on efficient solution of multiple-network porous media models

Our new preprints on Parameter robust preconditioning for multi-compartmental Darcy equations and Parameter robust preconditioning by congruence for multiple-network poroelasticity are out on arXiv!

In these studies, we consider systems of equations describing flow in multiple, interacting fluid networks with different permeabilities. The resulting discrete systems of equations can easily become very large, and efficient solution techniques are key. In cases where the inter-network interactions dominate the intra-network permeabilities, we observed that the systems of equations became increasingly hard to solve. In particular, the number of solver iterations would increase with increasing discrete resolution and interaction parameters. However, by creating an algorithm for defining tailored changes of variables, we were able to define robust preconditioners with uniformly bounded iteration counts.

Personal favorite: The lovely trick for how to define parameter-tailored changes of variables via classical matrix analysis including simultaneous diagonalization by congruence.

Fast uncertainty quantification of tracer distribution in the brain interstitial fluid

Our new preprint on Fast uncertainty quantification of tracer distribution in the brain interstitial fluid with multilevel and quasi Monte Carlo is now out on arXiv.

Meshes of brains with regional markers for defining output quantities of interest and bounding boxes for creating field representations
Figure: Meshes of brains with regional markers for defining output quantities of interest and bounding boxes for creating field representations

Mathematical models in biology involve many parameters that are uncertain or in some cases unknown. Over the last years, increased computing power has expanded the complexity and increased the number of degrees of freedom of many such models. For this reason, efficient uncertainty quantification algorithms are now needed to explore the often large parameter space of a given model. Advanced Monte Carlo methods such as quasi Monte Carlo (QMC) and multilevel Monte Carlo (MLMC) have become very popular in the mathematical, engineering, and financial literature for the quantification of uncertainty in model predictions. However, applying these methods to physiologically relevant simulations is a difficult task given the typical complexity of the models and geometries involved. In this paper, we design and apply QMC and MLMC methods to quantify uncertainty in a convection-diffusion model for tracer transport within the brain. We show that QMC outperforms standard Monte Carlo simulations when the number of random inputs is small. MLMC considerably outperforms both QMC and standard Monte Carlo methods and should therefore be preferred for brain transport models.

Finite element simulation of ionic electrodiffusion in cellular geometries

Our preprint on Finite element simulation of ionic electrodiffusion in cellular geometries is now out on arxiv!

Our mathematical model accounting for ionic electrodiffusion predicts weaker ephaptic coupling effects compared to previous models.

Mathematical models for excitable cells are commonly based on cable theory, which considers a homogenized domain and spatially constant ionic concentrations. Although such models provide valuable insight, the effect of altered ion concentrations or detailed cell morphology on the electrical potentials cannot be captured. In this paper, we discuss an alternative approach to detailed modelling of electrodiffusion in neural tissue. The mathematical model describes the distribution and evolution of ion concentrations in a geometrically-explicit representation of the intra- and extracellular domains. As a combination of the electroneutral Kirchhoff-Nernst-Planck (KNP) model and the Extracellular-Membrane-Intracellular (EMI) framework, we refer to this model as the KNP-EMI model. Here, we introduce and numerically evaluate a new, finite element-based numerical scheme for the KNP-EMI model, capable of efficiently and flexibly handling geometries of arbitrary dimension and arbitrary polynomial degree. Moreover, we compare the electrical potentials predicted by the KNP-EMI and EMI models. Finally, we study ephaptic coupling induced in an unmyelinated axon bundle and demonstrate how the KNP-EMI framework can give new insights in this setting.

Abstractions and algorithms for mixed domain and mixed dimensional finite elements

Our preprint on Abstractions and automated algorithms for mixed domain finite element methods is now out on arXiv!

Mixed dimensional partial differential equations (PDEs) are equations coupling unknown fields defined over domains of differing topological dimension. Such equations naturally arise in a wide range of scientific fields including geology, physiology, biology and fracture mechanics. Mixed dimensional PDEs are also commonly encountered when imposing non-standard conditions over a subspace of lower dimension e.g. through a Lagrange multiplier. In this paper, we present general abstractions and algorithms for finite element discretizations of mixed domain and mixed dimensional PDEs of co-dimension up to one (i.e. nD-mD with |n-m| <= 1). We introduce high level mathematical software abstractions together with lower level algorithms for expressing and efficiently solving such coupled systems. The concepts introduced here have also been implemented in the context of the FEniCS finite element software. We illustrate the new features through a range of examples, including a constrained Poisson problem, a set of Stokes-type flow models and a model for ionic electrodiffusion.

Paper quantifying uncertainties of tracer distribution in the brain published in FBCNS

Our paper on “Uncertainty quantification of parenchymal tracer distribution using random diffusion and convective velocity fields” is now published in Fluids and Barrier of the Central Nervous System. Our main findings are:

  • Uncertainty in diffusion parameters substantially impact the amount of tracer in gray and white matter, and the average tracer concentration in gray and white subregions.
  • Even with an uncertainty in the diffusion coefficient of a factor three, and a resulting fourfold variation in white matter tracer enhancement, discrepancies between simulations of diffusion and experimental data are too large to be attributed to uncertainties in the diffusion coefficient alone.
  • A convective velocity field modelling the glymphatic theory did not increase tracer enhancement in the brain parenchyma compared to pure diffusion. However, when a large-scale directional structure was added to this glymphatic velocity field, tracer inflow increased.
Illustration of velocity field modelling large-scale directionality associated with a glymphatic circulation

Intracranial pressure elevation alters CSF clearance pathways

Our preprint on Intracranial pressure elevation alters CSF clearance pathways is now out on bioRxiv! 

Infusion testing is a common procedure to determine whether shunting will be beneficial in patients with normal pressure hydrocephalus. The method has a well-developed theoretical foundation and corresponding mathematical models that describe the CSF circulation from the choroid plexus to the arachnoid granulations. Here, we investigate to what extent the proposed glymphatic or paravascular pathway (or similar pathways) modifies the results of the traditional mathematical models.

We used a two-compartment model consisting of the subarachnoid space and the paravascular spaces. For the arachnoid granulations, the cribriform plate, capillaries and paravascular spaces, resistances were calculated and used to estimate flow before and during an infusion test. Next, pressure in the subarachnoid space and paravascular spaces were computed. Finally, different variations to the model were tested to evaluate the sensitivity of selected parameters.

At baseline, we found a very small paravascular flow directed into the subarachnoid space, while 60% of the fluid left through the arachnoid granulations and 40% left through the cribriform plate. However, during the infusion, paravascular flow reversed and 25% of the fluid left through these spaces, while 60% went through the arachnoid granulations and only 15% through the cribriform plate.

The relative distribution of CSF flow to different clearance pathways depends on intracranial pressure (ICP), with the arachnoid granulations as the main contributor to outflow. As such, ICP increase is an important factor that should be addressed when determining the pathways of injected substances in the subarachnoid space.

New (old) paper published on automated adjoints of coupled PDE-ODE systems!

Our paper on Automated adjoints of coupled PDE-ODE systems is now published online in the SIAM Journal on Scientific Computing (SISC). Mathematical models that couple partial differential equations (PDEs) and spatially distributed ordinary differential equations (ODEs) arise in biology, medicine, chemistry, and many other fields. In this paper we discuss an extension to the FEniCS finite element software for expressing and efficiently solving such coupled systems. Given an ODE described using an augmentation of the Unified Form Language (UFL) and a discretization described by an arbitrary Butcher tableau, efficient code is automatically generated for the parallel solution of the ODE. The high-level description of the solution algorithm also facilitates the automatic derivation of the adjoint and tangent linearization of coupled PDE-ODE solvers. We demonstrate the capabilities of the approach on examples from cardiac electrophysiology and mitochondrial swelling.