Fast computation of sub-diffusion and wave equations with frequency dependent attenuation

We describe a new fast and oblivious quadrature for fractional time derivatives. This is based on a new approach to the computation of the inverse Laplace transform, that allows for a single contour and single quadrature to be used for all times. This significantly simplifies the implementation while retaining (and improving) the excellent efficiency of the original oblivious algorithm. As a first application, we describe two ways of solving the time-fractional sub-diffusion equation using these ideas. We end the talk with an application to acoustic wave equation with frequency dependent attenuation .

The talk is based on joint work with Katie Baker and Maria Lopez Fernandez.

A unified theory of fractional and nonlocal calculus and its consequences on nonlocal model discovery

Nonlocal and fractional models capture effects that classical partial differential equations cannot describe; for this reason, they are suitable for a broad class of engineering and scientific applications that feature multiscale or anomalous behavior. This has driven a desire for a vector calculus based on nonlocal and fractional derivatives to derive models of, e.g., subsurface transport, turbulence, and conservation laws. In the literature, several independent definitions and theories of nonlocal and fractional vector calculus have been put forward. This fragmented literature suffers from a lack of rigorous comparison and unified notation, hindering the development of nonlocal modeling. In this work we provide a unified theory and connect all the dots by defining a universal form of nonlocal operators under a theory that includes, as a special case, several well-known proposals for fractional vector calculus in the limit of infinite interactions.

Discovering the subdiffusion model in an unknown medium

The subdiffusion phenomenon is now widely recognized in many engineering and physical applications. The mathematical models for subdiffusion involve several parameters, e.g., diffusion coefficient, potential, initial and boundary conditions along with the order of derivation. Sometimes some of these parameters are not readily available, but one can measure additional information about the solution. Then one natural question is how much we can say about the mathematical model. In this talk, I will discuss several theoretical and computational results on determining the order of derivation and other parameters when the other problem data are not fully specified.

Preconditioning Methods for Fractional Diffusion Equations

In this talk, we discuss some preconditioning methods for fractional diffusion equations. Also preconditioning for time fractional diffusion inverse source problems is studied. Numerical examples are reported to demonstrate these preconditioning techniques.

Numerical methods for spectral fractional diffusion

We present three schemes for the numerical approximation of the spectral fractional Laplacian. The first method uses a rational approximation to approximate negative fractional powers of an operator. The second is a discretization of the so-called Balakrishnan formula, which expresses fractional powers of a positive operator. The third, is a PDE approach that exploits the extension to one higher dimension. We discuss pros and cons of each method, and error estimates. Time permitting, we extend some of these to more general problems, both steady, time dependent; linear, and nonlinear. We show illustrative simulations, applications, and mention challenging open questions.

Efficient and accurate space-time numerical methods for solving time-fractional PDEs

It is well known that solutions of time-fractional ODEs/PDEs exhibit weak singularities at $t=0$ so that usual approximation techniques do not yield accurate results.

In this talk, we present two classes of numerical methods to address this problem: The first one is based on the generalized Jacobi functions with enrichment, while the second one is based on a new set of log orthogonal functions (LOFs). Space-time Galerkin and Petrov-Galerkin approaches will be used. We shall show rigorous error analysis with high-order convergence rate for both methods despite the weak singularities at $t=0$, and present ample numerical results to demonstrate the effectiveness of our methods.

Second-order error analysis of the averaged L1 scheme for a time-fractional subdiffusion problem.

Subdiffusion initial-boundary value problems with a Caputo temporal derivative of order $\alpha\in (0,1)$ are considered. An averaged variant of the well-known L1 scheme is proved to be $O(N^{-2})$ convergent on suitably graded meshes with $N$ points, improving the $O(N^{-(2-\alpha)})$ convergence rate of the basic L1 scheme. The analysis relies on a delicate decomposition of the temporal truncation error that yields a sharp dependence of the order of convergence on the degree of mesh grading used. This averaged L1 scheme can be combined with finite difference or finite element discretisations in space, for which error bounds are derived that have sharp orders of convergence in both time and space. Numerical experiments support our results.

Spline collocation methods for fractional differential equations: theoretical and computational aspects

Spline collocation methods are a powerful tool to discretize fractional differential equations, since they have a high order of convergence and strong stability properties. In this talk, we illustrate the convergence and stability analysis of one and two step spline collocation methods. We pay attention also to the efficient implementation of these methods, which requires the evaluation of fractional integrals. Some numerical experiments are provided to confirm theoretical results and to compare one and two step collocation methods. This is a joint work with B. Paternoster and D. Conte.

Mittag-Leffler functions and their applications in network science

In network analysis measures of centrality identify the most important vertices within a graph. In this talk we will describe a complete theory for walk-based centrality indices defined in terms of Mittag-Leffler functions. This overarching theory includes as special cases well-known centrality measures like subgraph centrality and Katz centrality. We will introduce a family of indices parametrized by two numbers; by letting these vary, one can show that Mittag-Leffler centralities interpolate between degree and eigenvector centrality, as well as between resolvent-based and exponential-based indices. To provide guidelines on parameter selection we will discuss some modeling and computational issues. The theory is then extended to the case of networks that evolve over time.

Hadamard-type fractional heat equations and ultra-slow diffusions

Ultra-slow diffusion processes include a wide class of different stochastic processes characterized by a logarithmic growth of the mean squared displacement. In this talk we will consider diffusion equations involving Hadamard-type time-fractional derivatives related to ultra-slow random models. We first analyze the abstract Cauchy problem involving the Hadamard-type fractional derivative and discuss the stochastic interpretation of its solution, by using the theory of time-changed processes. Then, we consider in detail the particular case of heat-type equations based on the Hadamard-type fractional derivative that can be directly related to ultra-slow diffusions. The ultra-slow behavior emerges from the explicit form of the variance of the random process arising from our analysis. We finally underline the utility of the applications of the more general fractional derivatives w.r.t. another function in diffusive models and the possible future applications.

Variable-order fractional calculus: a change of perspective

Derivatives and integrals of variable fractional order find application in several fields and a number of attempts to generalize constant-order operators have been made. In this talk we present a completely different approach based on the seminal ideas presented in 1972 by the Italian engineer G.B.Scarpi. Instead of generalizing the usual definitions of constant-order operators in the time domain (as it is usually proposed), it is the Laplace transform representation of the standard fractional-order derivative and integral which is extended to cover the variable-order case. This simple but powerful idea has been recently framed in the robust mathematical theory of General Fractional Derivatives and Integrals based on Sonine equations. Since the analytical representation of these variable-order fractional derivatives and integrals is possible only in the Laplace transform domain, the proposed approach requires an extensive use of numerical techniques for the accurate inversion of the Laplace transform which are discussed in this talk.

The semigroup approximation to Fractional Elliptic equations

We present a numerical discretization that combines quadratures rules with finite element methods for the spectral Laplacian associated to different type of boundary conditions. Relying on the integral formulation of the operator via the heat-semigroup formalism, such formulation can handle at once cartesian domains and multi-dimensional (possibly irregular) geometries. We discuss the numerical approximation of the corresponding fractional Poisson problem, and we show super quadratic convergence rates up to $O(h^4)$ for sufficiently regular data, or simply $O(h^{2s})$ for data in $L^2(\Omega)$. We illustrate our findings with some numerical tests.

Fractional Newtonian Gravity and Galactic Dynamics

In this presentation I provide a derivation of some characteristic effects of Milgrom’s modified Newtonian dynamics (MOND) from a fractional version of Newton’s theory based on the fractional Poisson equation. Specifically, I employ the properties of the fractional Laplacian to investigate the features of the fundamental solution of the proposed model. Then, taking advantage of the Tully-Fisher relation, as the fractional order $s$ approaches 3/2, I relate the typical length scale emerging from this modification of Newton’s gravity with the critical acceleration $a_0$ of MOND. Finally, I explain the need for a variable-order version of the proposed model in order to properly capture the phenomenology of galactic dynamics.

Exact and Numerical Solutions of Fractional Reaction–Diffusion Equations

We present a procedure that combines the Lie symmetry analysis with the numerical methods to get exact and numerical solutions of fractional reaction–diffusion equations. By the Lie symmetries, the original model, defined in terms of the Riemann-Liouville fractional derivatives, is reduced to fractional differential equations with smaller number of independent variables or, also, to fractional ordinary differential equations. By introducing the Caputo derivative, the reduced equations are numerically solved by the classical numerical methods. Numerical results show the efficiency and the reliability of the proposed approach for solving a wide class of fractional models.

Rational Krylov for evaluating inverse fractional powers: convergence and pole selection

Evaluating the action of a matrix function on a vector, that is $x=f(\mathcal M)v$, is a ubiquitous task in applications. When the matrix $\mathcal M$ is large, subspace projection methods, such as the rational Krylov method, are usually employed. In this work, we provide a quasi-optimal pole choice for rational Krylov methods when $f$ is either Cauchy-Stieltjes or Laplace-Stieltjes for a positive definite matrix $\mathcal M$. This includes the case $f(z)=z^{-\alpha}$ with $\alpha>0$, which arises when solving fractional diffusion problems.

Then, we consider the case when the argument $\mathcal M$ has the Kronecker structure $\mathcal M=I \otimes A + B \otimes I$, and $v$ is obtained by vectorizing a low-rank matrix. This finds application, for instance, in solving two-dimensional FDEs on rectangular domains. We introduce an error analysis for the numerical approximation of $x$. Pole choices and explicit convergence bounds are given also in this case.

On B-spline collocation matrices for (Riemann-Liouville or Caputo) Riesz fractional operator and their spectral properties

In this work, we focus on a fractional differential equation in Riesz form with Riemann-Liouville fractional derivatives discretized by a polynomial B-spline collocation method. For an arbitrary polynomial degree $p$, we show that the resulting coefficient matrices possess a Toeplitz-like structure. We investigate their spectral properties via their symbol and we prove that, like for second order differential problems, also in this case the given matrices are ill-conditioned both in the low and high frequencies for large $p$. More precisely, in the fractional scenario the symbol has a single zero at $0$ of order $\alpha$, with $\alpha$ the fractional derivative order that ranges from $1$ to $2$, and it presents an exponential decay to zero at $\pi$ for increasing $p$ that becomes faster as $\alpha$ approaches $1$. This translates in a mitigated conditioning in the low frequencies and in a deterioration in the high frequencies when compared to second order problems. Moreover, we perform a numerical study of the approximation behavior of polynomial B-spline collocation. This study suggests that, in line with non-fractional diffusion problems, the approximation order for smooth solutions in the fractional case is $p+2-\alpha$ for even $p$, and $p+1-\alpha$ for odd $p$. Finally, we discuss what changes in terms of matrices and related spectral distribution when replacing the Riemann-Liouville fractional derivatives with the Caputo ones.

A truncated Laguerre approach for the fractional power of operators

In this work we consider the numerical approximation of $\mathcal{L}^{-\alpha },$ $0<\alpha <1$, where $\mathcal{L}$ is a self-adjoint positive operator acting on a separable Hilbert space $\mathcal{H}$, with spectrum $\sigma (\mathcal{L})\subseteq \lbrack c,+\infty ), c>0$. This problem finds immediate application when solving equations involving a fractional diffusion term like $(-\Delta )^{\alpha }$ where $\Delta $ denotes the standard Laplacian.

By exploiting the existing representations of the function $\lambda^{-\alpha}$ in terms of contour integrals, after suitable changes of variable and quadrature rules one typically finds rational approximations of the type

$\mathcal{L}^{-\alpha }\approx \mathcal{R}_{n-1,n}(\mathcal{L}),$ $\mathcal{R}_{n-1,n}(\lambda )=\frac{p_{n-1}(\lambda )}{q_{n}(\lambda )}$, $p_{n-1}\in \Pi _{n-1},\,q_{n}\in \Pi _{n}$

where $n$ is equal or closely related to the number of points of the quadrature formula.

In this work we consider a particular integral representation that leads to the the use of a truncated Gauss-Laguerre rule. We show that the error decays exponentially and that is well estimated by

$ 8\sin (\alpha \pi )\exp \left( -5\alpha ^{1/2}k_{n}^{1/2}\right) $

where $k_{n}$ is number of operator inversions. The method is then very fast and the quality of the above estimate allows to define a-priori $n$ and $k_{n}$ in order to fulfil a given tolerance. The poles can also be fruitfully used to construct a rational Krylov method for applications where only the action of $\mathcal{L}^{-\alpha }$ is needed.

Fast Solvers for Optimization Problems with FDE Constraints

In this talk we consider fast iterative solvers and preconditioners for optimization problems where fractional differential equations act as constraints. A focus is on the multilevel Toeplitz structure of a number of such problems, for which it is found to be valuable to apply multilevel circulant preconditioners to accelerate iterative methods. We also examine problems with additional bound constraints imposed on the solution, and we devise an Alternating Direction Method of Multipliers (ADMM) framework coupled with iterative solvers and multilevel circulant preconditioners for the resulting linear systems. Although we work with dense linear systems, we are able to keep the storage requirements linear with respect to grid size N, while ensuring order N log N computational complexity. We also briefly discuss some other approaches for preconditioning FDE-constrained optimization problems, including a low-rank tensor-product solver which will be outlined in detail by another speaker.

Spectral methods in Peridynamics models

Peridynamics is a non local theory for dynamic fracture analysis consisting in a second order in time partial integro-differential equation. In this talk, we describe numerical methods for both the linear and the non linear models. In particular, we implement a spectral method with volume penalization for the space discretization based on the Fourier expansion of the solution, while we consider the Newmark-β and the Stormer-Verlet method for the time marching. This computational approach takes advantages from the convolutional form of the peridynamic operator and from the use of the discrete Fourier transform. We show a convergence result for the fully discrete approximation and study the stability of the method applied to the linear peridynamic model. We additionally study the eigenvalues of the peridynamic operator and perform several numerical tests and comparisons to validate our results.

On the numerical solution of fractional differential equations by spline quasi-interpolant operators

The nonlocal behavior of the fractional derivative poses a challenge when constructing numerical methods to solve fractional differential problems. Collocation methods seem to provide a good compromise between accuracy and efficiency. In fact, collocation methods are global methods that approximate the solution in the whole discretization interval, thus taking into account the nonlocal behavior of the fractional derivative. On the other side, the use of local operators to approximate the solution allows us to keep a low computational cost. In this talk, I will present some recent results on the numerical solution of some fractional differential problems by a collocation method based on spline quasi-interpolant operators. Quasi-interpolant operators are approximating operators that reproduce polynomials up to a given degree. For this reason, they have a great flexibility that can be used to preserve special properties of the function to be approximated. I will show the performance of the proposed method through numerical experiments.

On the efficient numerical computation of the matrix Mittag-Leffler function

Important applications in fractional calculus require the numerical computation of the Mittag-Leffler function with matrix arguments. This topic presents sensitive issues that need to be properly addressed. Furthermore, ad hoc numerical strategies are required when large matrices are involved. Numerical methods will be presented, along with numerical tests to demonstrate their effectiveness

Rational Krylov methods for fractional diffusion problems on graphs

The solution to the fractional diffusion equation on a directed network can be expressed as $f(L^T) b$, where $L$ is the graph Laplacian and $f$ is a function involving fractional powers, with a branch cut on the negative real line. Since $L$ is a singular matrix, Krylov methods for $f(L^T) b$ may converge slowly. To achieve faster convergence, we propose to use rational Krylov methods applied to a desingularized version of the graph Laplacian. We present three different desingularization strategies and we compare their effectiveness with numerical experiments on real-world networks.

Low-rank solvers for fractional equations in simulation and optimization

When discretized with finite difference schemes FDEs lead to dense but highly structured linear systems. In particular, we see a Kronecker product structure involving Toeplitz matrices. As the curse of dimensionality hits we want to use this structure to approximate the solution to the linear system in low rank form where we approximate the solution with high accuracy but a fraction of the storage cost. We will illustrate this both in the simulation as well as the optimization setup.

Anisotropic multigrid preconditioners for space-fractional diffusion equations

In this talk, we focus on a two-dimensional space-fractional diffusion problem discretized by means of a second order finite difference scheme obtained as combination of the Crank-Nicolson scheme and the so-called weighted and shifted Grünwald formula. Efficient multigrid strategies for the resulting Toeplitz-like linear systems have shown to be effective when the fractional orders are both close to 2. Here we seek to investigate how multigrid approaches can be efficiently extended to the case where only one of the two fractional order is close to 2, while the other is close to 1. In other words, we consider space-fractional diffusion problems that involve an intrinsic anisotropy in the direction corresponding to the minimum fractional order. We design a multigrid preconditioner where the grid transfer operator is obtained with a semicoarsening technique and the smoothing is performed with relaxed Jacobi whose damping parameter is accurately estimated by using a symbol-based approach. Moreover, for large-sized problems a further improvement in the robustness of the multigrid method can be reached using a V-cycle with semicoarsening as smoother. The scaled-Laplacian matrix is used in the direction where the fractional derivative order is close to 2, while in the other direction a rediscretization is adopted. Several numerical results confirm that the resulting multigrid preconditioner is computationally effective and outperforms current state of the art techniques.

Fractional exponential sums

The computation of the inverse of the $\alpha$-th power of a Kroneker sum is an important task for many application; for instance in the numerical resolution of fractional PDEs is often needed the computation of $(I\otimes A+B\otimes I)^{-\alpha},$ for some matrices $A,B$. We propose a way to approximate the function $f(\xi)=\xi^{-\alpha}$, for $\alpha \in (0,1)$ with a sum of exponentials, that is \(f(\xi) \approx \sum_{j = 1}^k w_j e^{-\beta_j \xi}.\) This allows an easy evaluation in a Kroneker sum, by exploiting the Kroneker form of the exponential of a Kroneker sum (i.e., $e^{I\otimes A + I\otimes B}=e^A\otimes e^B$). Moreover, we provide upper bounds on the errors computed using such approximations.

Sum of Independent Generalized Mittag-Leffler Random Variables and Related Fractional Processes

We study the distribution of the sum of independent and non-identically distributed generalized Mittag-Leffler random variables by means of the Laplace transform.

The probability density can be expressed in terms of a multivariate analogue of the MittagLeffler function. We then apply the result on the convolution of generalized Mittag-Leffler to study some state-dependent fractional point processes. We present their explicit probability mass functions as well as their connections with the fractional integral/differential equations. In particular, we focus on the point process whose waiting times are given by independent Mittag-Leffler random variables which alternates $N$ couples of parameters $(\nu_i,\lambda_i) \in (0,1] \times (0,+\infty)$, $i=1,\dots,N$.

In the case of two alternating kinds of Mittag-Leffler distributed waiting times we present the conditional arrival times and we show an application to the telegraph process.

Hadamard-Type Fractional Heat Equations on bounded domains

Since the first time a fractional Cauchy problem was considered [Nigmatullin, 1986], the concern about this topic has increased. Nowadays many fractional operators have been defined: Caputo, Riemann-Liouville, Prabhakar are the most famous cases which can describe different kinds of physical phenomena. In this work, we focus on the diffusion problem led by Hadamard type fractional derivative operator on a bounded domain in the d-dimensional euclidean space. The results will be given by a probabilistic and analytic point of view.

Fractional-Order Perona–Malik Diffusion for Image Denoising

In 1987, Perona and Malik presented a denoising technique expressed in terms of a nonlinear partial differential equation that includes a diffusion coefficient given as a function of the gradient of the image. Their approach has recently been improved by means of the combination of integer-order operators with a diffusion coefficient that is a function of the fractional gradient of the image. In this talk, we consider the Grünwald–Letnikov definition of fractional derivatives and we give an idea on how to discretize the fractional-order anisotropic diffusion equation on a square lattice, using central finite differences with 8 points.

Continuous extensions of numerical methods for Stochastic Fractional Differential Equations

We present a technique to provide continuous-time extension of numerical methods solving Stochastic Fractional Differential Equations (SFDEs). The basic idea we follow is closely related to the classic scenario of deterministic collocation methods for ordinary differential equations, useful to provide accurate error estimates and to perform a variable step-size implementation. The building block of this analysis is the continuous extension of Euler Maruyama method, whose effectiveness is also confirmed by selected numerical experiments. This is a joint work with B. Paternoster, R. D’Ambrosio and D. Conte.

Equation dependent numerical methods for FDEs

In this talk we describe techniques that allow to enlarge the absolute stability regions of classical explicit numerical methods for Ordinary Differential Equations (ODEs). The basic idea of these techniques consists in the modification of method coefficients, which result in depending on the Jacobian of the ODE to be solved. We then analyze the possibility to apply these methodologies to numerical methods for Fractional Differential Equations (FDEs), in order to obtain an improvement in terms of accuracy and stability properties. This is a joint work with Prof. Beatrice Paternoster and Dajana Conte.