Stefano Serra-Capizzano

Carlo Garoni, Stefano Serra-Capizzano

The theory of generalized locally Toeplitz (GLT) sequences is an apparatus for computing the spectral and singular value distribution of sequences of matrices that possess a (possibly hidden) Toeplitz-like structure. These sequences, which are known as GLT sequences, arise in several applications, including the discretization of differential equations. Associated with any GLT sequence is a special function called symbol. In this paper, we prove that, if $\{A_n\}_n$ is a GLT sequence with symbol $κ$ and $P_n$ is any block Jacobi or block Gauss-Seidel preconditioner for $A_n$ with a fixed number of blocks independent of $n$, then $\{P_n\}_n$ is a GLT sequence with symbol $κ$, just like $\{A_n\}_n$. This result allows us to predict a remarkable efficiency of block Jacobi/Gauss-Seidel preconditioning for GLT sequences, which is in fact illustrated through numerical experiments. It also allows us to extend the Fasino-Tilli theorem on the zero distribution of Hankel matrix sequences generated by $L^1$ functions to a larger class of matrix sequences called generalized locally Hankel (GLH) sequences.

Mariarosa Mazza, Carla Manni, Ahmed Ratnani et al.

Alfvén-like operators are of interest in magnetohydrodynamics, which is used in plasma physics to study the macroscopic behavior of plasma. Motivated by this important and complex application, we focus on a parameter-dependent curl-div problem that can be seen as a prototype of an Alfvén-like operator, and we discretize it using isogeometric analysis based on tensor-product B-splines. The involved coefficient matrices can be very ill-conditioned, so that standard numerical solution methods perform quite poorly here. In order to overcome the difficulties caused by such ill-conditioning, a two-step strategy is proposed. First, we conduct a detailed spectral study of the coefficient matrices, highlighting the critical dependence on the different physical and approximation parameters. Second, we exploit such spectral information to design fast iterative solvers for the corresponding linear systems. For the first goal we apply the theory of (multilevel block) Toeplitz and generalized locally Toeplitz sequences, while for the second we use a combination of multigrid techniques and preconditioned Krylov solvers. Several numerical tests are provided both for the study of the spectral problem and for the solution of the corresponding linear systems.

Marco Donatelli, Stefano Serra-Capizzano, Debora Sesana

Starting from the spectral analysis of g-circulant matrices, we consider a new multigrid method for circulant and Toeplitz matrices with given generating function. We assume that the size n of the coefficient matrix is divisible by g \geq 2 such that at the lower level the system is reduced to one of size n/g by employing g-circulant based projectors. We perform a rigorous two-grid convergence analysis in the circulant case and we extend experimentally the results to the Toeplitz setting, by employing structure preserving projectors. The optimality of the proposed two-grid method and of the multigrid method is proved, when the number theta \in N of recursive calls is such that 1 < theta < g. The previous analysis is used also to overcome some pathological cases, in which the generating function has zeros located at "mirror points" and the standard two-grid method with g = 2 is not optimal. The numerical experiments show the correctness and applicability of the proposed ideas both for circulant and Toeplitz matrices.

Minghua Chen, Sven-Erik Ekström, Stefano Serra-Capizzano

The nonlocal problems have been used to model very different applied scientific phenomena, which involve the fractional Laplacian when one looks at the Lévy processes and stochastic interfaces. This paper deals with the nonlocal problems on a bounded domain, where the stiffness matrices of the resulting systems are Toeplitz-plus-tridiagonal and far from being diagonally dominant, as it occurs when dealing with linear finite element approximations. By exploiting a weakly diagonally dominant Toeplitz property of the stiffness matrices, the optimal convergence of the two-grid method is well established [Fiorentino and Serra-Capizzano, {\em SIAM J. Sci. Comput.}, {17} (1996), pp. 1068--1081; Chen and Deng, {\em SIAM J. Matrix Anal. Appl.}, {38} (2017), pp. 869--890]; and there are still questions about best ways to define coarsening and interpolation operator when the stiffness matrix is far from being weakly diagonally dominant [Stüben, {\em J. Comput. Appl. Math.}, {128} (2001), pp. 281--309]. In this work, using spectral indications from our analysis of the involved matrices, the simple (traditional) restriction operator and prolongation operator are employed in order to handle general algebraic systems which are {\em neither Toeplitz nor weakly diagonally dominant} corresponding to the fractional Laplacian kernel and the constant kernel, respectively. We focus our efforts on providing the detailed proof of the convergence of the two-grid method for such situations. Moreover, the convergence of the full multigrid is also discussed with the constant kernel. The numerical experiments are performed to verify the convergence with only $\mathcal{O}(N \mbox{log} N)$ complexity by the fast Fourier transform, where $N$ is the number of the grid points.

Zheng-Jian Bai, Marco Donatelli, Stefano Serra-Capizzano

In recent works several authors have proposed the use of precise boundary conditions (BCs) for blurring models and they proved that the resulting choice (Neumann or reflective, anti-reflective) leads to fast algorithms both for deblurring and for detecting the regularization parameters in presence of noise. When considering a symmetric point spread function, the crucial fact is that such BCs are related to fast trigonometric transforms. In this paper we combine the use of precise BCs with the Total Variation (TV) approach in order to preserve the jumps of the given signal (edges of the given image) as much as possible. We consider a classic fixed point method with a preconditioned Krylov method (usually the conjugate gradient method) for the inner iteration. Based on fast trigonometric transforms, we propose some preconditioning strategies which are suitable for reflective and anti-reflective BCs. A theoretical analysis motivates the choice of our preconditioners and an extensive numerical experimentation is reported and critically discussed. The latter shows that the TV regularization with anti-reflective BCs implies not only a reduced analytical error, but also a lower computational cost of the whole restoration procedure over the other BCs.

Danyal Ahmad, Stefano Serra-Capizzano, Muhammad Sohaib et al.

We consider an initial boundary value problem of the space fractional Allen-Cahn equation with logarithmic Flory-Huggins potential. As an approximation technique, first-order weighted and shifted Grunwald difference formulae of the left and right fractional derivatives are used. The main focus of the present work is to study the spectral features of the underlying matrices and matrix sequences and to design proper preconditioners based on the spectral information. Then a computational and spectral analysis of the resulting preconditioned matrix sequences is performed. Numerical evidence and a short list of open problems complete the current study.

Pietro Dell'Acqua, Stefano Serra-Capizzano, Cristina Tablino Possio

Inspired by the theoretical results on optimal preconditioning stated by Ng, R.Chan, and Tang in the framework of Reflective boundary conditions (BCs), in this paper we present analogous results for Anti-Reflective BCs, where an additional technical difficulty is represented by the non orthogonal character of the Anti-Reflective transform and indeed the technique of Ng, R.Chan, and Tang can not be used. Nevertheless, in both cases, the optimal preconditioner is the blurring matrix associated to the symmetrized Point Spread Function (PSF). The geometrical idea on which our proof is based is very simple and general, so it may be useful in the future to prove theoretical results for new proposed boundary conditions. Computational results show that the preconditioning strategy is effective and it is able to give rise to a meaningful acceleration both for slightly and highly non-symmetric PSFs.

Allal Guessab, Federico Nudo, Stefano Serra-Capizzano

In this paper we study a weighted histopolation problem on $[-1,1]$ associated with Jacobi weights. In the first part of the present work we prove results in approximation theory, while in the second we analyze the resulting matrices from an asymptotic linear algebra perspective. More in detail, in the first part, given weighted cell averages, we construct a reconstruction operator based on weighted primitives of Jacobi polynomials and investigate the resulting discretization matrices. At any fixed discretization level, we derive an exact factorization of the histopolation matrix through a backward-difference operator and a sampling operator of Jacobi weighted primitives. Combining a sharp integration by parts identity with the three-term recurrence of Jacobi polynomials, we further show that the primitive sampling operator admits an explicit decomposition involving a tridiagonal coupling matrix in the Jacobi spectral index. This yields a tridiagonal factor representation of the histopolation matrix. In the second part, under standard mesh-regularity assumptions, we show that all the various induced matrix sequences belong to the Generalized Locally Toeplitz (GLT) class, by describing in detail the related GLT symbols. As a consequence, we provide the corresponding spectral distributions and discuss their implications for numerical stability when solving the associated linear systems.

Stefano Serra-Capizzano, Cristina Tablino Possio

Given a multigrid procedure for linear systems with coefficient matrices $A_n$, we discuss the optimality of a related multigrid procedure with the same smoother and the same projector, when applied to properly related algebraic problems with coefficient matrices $B_n$: we assume that both $A_n$ and $B_n$ are positive definite with $A_n\le \vartheta B_n$, for some positive $\vartheta$ independent of $n$. In this context we prove the Two-Grid method optimality. We apply this elementary strategy for designing a multigrid solution for modifications of multilevel structured (Toeplitz, circulants, Hartley, sine ($τ$ class) and cosine algebras) linear systems, in which the coefficient matrix is banded in a multilevel sense and Hermitian positive definite. In such a way, several linear systems arising from the approximation of integro-differential equations with various boundary conditions can be efficiently solved in linear time (with respect to the size of the algebraic problem). Some numerical experiments are presented and discussed, both with respect to Two-Grid and multigrid procedures.

Samuele Ferri, Chiara Giraudo, Valerio Loi et al.

In the present paper, we analyze in detail the spectral features of the matrix sequences arising from the Taylor-Hood $\mathbb{P}_2$-$\mathbb{P}_1$ approximation of variable viscosity for $2d$ Stokes problem under weak assumptions on the regularity of the diffusion. Localization and distributional spectral results are provided, accompanied by numerical tests and visualizations. A preliminary study of the impact of our findings on the preconditioning problem is also presented. A final section with concluding remarks and open problems ends the current work.

Ayoub Lailoune, Valerio Loi, Stefano Serra-Capizzano

The present work contains a comprehensive treatment of Weyl eigenvalue and singular value distributions for single-axis quaternion block multilevel Toeplitz matrix sequences generated by $s\times t$ quaternion matrix-valued, $d$-variate, Lebesgue integrable generating functions. Furthermore, in view of concrete applications, we are interested in preconditioning and matrix approximation results. To this end, a crucial step is the extension of the notion of an approximating class of sequences (a.c.s.) to the case of matrix sequences with quaternion entries, since it allows us to decompose the difference between a matrix and its preconditioner into low-norm plus (relatively) low-rank terms. As a specific example, we consider classes of quaternion block multilevel circulant matrix sequences as an a.c.s. for quaternion block multilevel Toeplitz matrix sequences. These approximation results lay the foundations for fast preconditioning methods when dealing with large quaternion linear systems stemming from modern applications. We conclude our study with numerical experiments and directions for future research.

Muhammad Faisal Khan, Asim Ilyas, Rolf Krause et al.

This paper studies the spectral properties of large matrices and the preconditioning of linear systems, arising from the finite difference discretization of a time-dependent space-fractional diffusion equation with a variable coefficient $a(x)$ defined on $Ω\subset \mathbb{R}^d$, $d=1,2$. The model involves a one-sided Riemann-Liouville fractional derivative multiplied by the function $a(x)$, discretized by the shifted Gr"unwald formula in space and the $θ$-method in time. The resulting all-at-once linear systems exhibit a $(d+1)$-level Toeplitz-like matrix structure, with $d=1,2$ denoting the space dimension, while the additional level is due to the time variable. A preconditioning strategy is developed based on the structural properties of the discretized operator. Using the generalized locally Toeplitz (GLT) theory, we analyze the spectral distribution of the unpreconditioned and preconditioned matrix sequences. The main novelty is that the analysis fully covers the case where the variable coefficient $a$ is nonconstant. Numerical results are provided to support the GLT based theoretical findings, and some possible extensions are briefly discussed.

Ludovico Bruni Bruno, Stefano Serra-Capizzano

Histopolation is the approximation procedure that associates a degree $ d-1 $ polynomial $ p_{d-1} \in \mathscr{P}_{d-1} (I) $ with a locally integrable function $ f $ imposing that the integral (or, equivalently, the average) of $p$ coincides with that of $f$ on a collection of $ d $ distinct segments $s_i$. In this work we discuss unisolvence and conditioning of the associated matrices, in an asymptotic linear algebra perspective, i.e., when the matrix-size $d$ tends to infinity. While the unisolvence is a rather sparse topic, the conditioning in the unisolvent setting has a uniform behavior: as for the case of standard Vandermonde matrix-sequences with real nodes, the conditioning is inherently exponential as a function of $d$ when the monomial basis is chosen. In contrast, for an appropriate selection of supports, the Chebyshev polynomials of second kind exhibit a bounded conditioning. A linear behavior is also observed in the Frobenius norm.

Paola Ferrari, Isabella Furci, Sean Hon et al.

Given a Lebesgue integrable function $f$ over $[0,2π]$, we consider the sequence of matrices $\{Y_nT_n[f]\}_n$, where $T_n[f]$ is the $n$-by-$n$ Toeplitz matrix generated by $f$ and $Y_n$ is the flip permutation matrix, also called the anti-identity matrix. Because of the unitary character of $Y_n$, the singular values of $T_n[f]$ and $Y_n T_n[f]$ coincide. However, the eigenvalues are affected substantially by the action of the matrix $Y_n$. Under the assumption that the Fourier coefficients are real, we prove that $\{Y_nT_n[f]\}_n$ is distributed in the eigenvalue sense as \[ ϕ_g(θ)=\left\{ \begin{array}{cc} g(θ), & θ\in [0,2π], -g(-θ), & θ\in [-2π,0), \end{array} \right.\, \] with $g(θ)=|f(θ)|$. We also consider the preconditioning introduced by Pestana and Wathen and, by using the same arguments, we prove that the preconditioned sequence is distributed in the eigenvalue sense as $ϕ_1$, under the mild assumption that $f$ is sparsely vanishing. We emphasize that the mathematical tools introduced in this setting have a general character and in fact can be potentially used in different contexts. A number of numerical experiments are provided and critically discussed.

Sean Hon, Mohammad Ayman Mursaleen, Stefano Serra-Capizzano

The singular value and spectral distribution of Toeplitz matrix sequences with Lebesgue integrable generating functions is well studied. Early results were provided in the classical Szeg{ő} theorem and the Avram-Parter theorem, in which the singular value symbol coincides with the generating function. More general versions of the theorem were later proved by Zamarashkin and Tyrtyshnikov, and Tilli. Considering (real) nonsymmetric Toeplitz matrix sequences, we first symmetrize them via a simple permutation matrix and then we show that the singular value and spectral distribution of the symmetrized matrix sequence can be obtained analytically, by using the notion of approximating class of sequences. In particular, under the assumption that the symbol is sparsely vanishing, we show that roughly half of the eigenvalues of the symmetrized Toeplitz matrix (i.e. a Hankel matrix) are negative/positive for sufficiently large dimension, i.e. the matrix sequence is symmetric (asymptotically) indefinite.

Minghua Chen, Weihua Deng, Stefano Serra-Capizzano

In this paper we derive new uniform convergence estimates for the V-cycle MGM applied to symmetric positive definite Toeplitz block tridiagonal matrices, by also discussing few connections with previous results. More concretely, the contributions of this paper are as follows: (1) It tackles the Toeplitz systems directly for the elliptic PDEs. (2) Simple (traditional) restriction operator and prolongation operator are employed in order to handle general Toeplitz systems at each level of the recursion. Such a technique is then applied to systems of algebraic equations generated by the difference scheme of the two-dimensional fractional Feynman-Kac equation, which describes the joint probability density function of non-Brownian motion. In particular, we consider the two coarsening strategies, i.e., doubling the mesh size (geometric MGM) and Galerkin approach (algebraic MGM), which lead to the distinct coarsening stiffness matrices in the general case: however, several numerical experiments show that the two algorithms produce almost the same error behaviour.

Davide Bianchi, Alessandro Buccini, Marco Donatelli et al.

Fractional Tikhonov regularization methods have been recently proposed to reduce the oversmoothing property of the Tikhonov regularization in standard form, in order to preserve the details of the approximated solution. Their regularization and convergence properties have been previously investigated showing that they are of optimal order. This paper provides saturation and converse results on their convergence rates. Using the same iterative refinement strategy of iterated Tikhonov regularization, new iterated fractional Tikhonov regularization methods are introduced. We show that these iterated methods are of optimal order and overcome the previous saturation results. Furthermore, nonstationary iterated fractional Tikhonov regularization methods are investigated, establishing their convergence rate under general conditions on the iteration parameters. Numerical results confirm the effectiveness of the proposed regularization iterations.

Matteo Semplice, Marco Donatelli, Stefano Serra-Capizzano

The novel contribution of this paper relies in the proposal of a fully implicit numerical method designed for nonlinear degenerate parabolic equations, in its convergence/stability analysis, and in the study of the related computational cost. In fact, due to the nonlinear nature of the underlying mathematical model, the use of a fixed point scheme is required and every step implies the solution of large, locally structured, linear systems. A special effort is devoted to the spectral analysis of the relevant matrices and to the design of appropriate iterative or multi-iterative solvers, with special attention to preconditioned Krylov methods and to multigrid procedures: in particular we investigate the mutual benefit of combining in various ways suitable preconditioners with V-cycle algorithms. Numerical experiments in one and two spatial dimensions for the validation of our multi-facet analysis complement this contribution.

Eric Ngondiep, Stefano Serra-Capizzano, Debora Sesana

For a given nonnegative integer alpha, a matrix A_{n} of size n is called alpha-Toeplitz if its entries obey the rule A_{n}=[a_{r-alpha*s}]_{r,s=0}^{n-1}. Analogously, a matrix A_{n} again of size n is called alpha-circulant if A_{n}= [a_{(r-alpha*s)mod n}]_{r,s=0}^{n-1}. Such kind of matrices arises in wavelet analysis, subdivision algorithms and more generally when dealing with multigrid/multilevel methods for structured matrices and approximations of boundary value problems. In this paper we study the singular values of alpha-circulants and we provide an asymptotic analysis of the distribution results for the singular values of alpha-Toeplitz sequences in the case where {a_{k}} can be interpreted as the sequence of Fourier coeffcients of an integrable function f over the domain (-pi;pi). Some generalizations to the block, multilevel case, amounting to choose f multivariate and matrix valued, are briefly considered.

Leonid Golinskii, Stefano Serra-Capizzano

Under the mild trace-norm assumptions we show that the eigenvalues of a generic (non Hermitian) complex perturbation of a Jacobi matrix sequence (not necessarily real) are still distributed as the real-valued function $2\cos t$ on $[0,π]$, which characterizes the nonperturbed case. In this way the real interval $[-2,2]$ is still a cluster for the asymptotic joint spectrum and, moreover, $[-2,2]$ attracts strongly (with infinite order) the perturbed matrix sequence. The results follow in a straightforward way from more general facts that we prove in an asymptotic linear algebra framework and are plainly generalized to the case of matrix-valued symbol, which arises when dealing with orthogonal polynomials with asymptotically periodic recurrence coefficients.