  In this work, we develop a novel preconditioned method for solving
space-fractional diffusion equations, which both accounts for and improves upon
an ideal preconditioner pioneered in [J. Pestana. Preconditioners for
symmetrized Toeplitz and multilevel Toeplitz matrices. SIAM Journal on Matrix
Analysis and Applications, 40(3):870-887, 2019]. Following standard
discretization on the equation, the resultant linear system is a non-symmetric,
multilevel Toeplitz system. Through a simple symmetrization strategy, we
transform the original linear system into a symmetric multilevel Hankel system.
Subsequently, we propose a symmetric positive definite multilevel Tau
preconditioner for the symmetrized system, which can be efficiently implemented
using discrete sine transforms. Theoretically, we demonstrate that
mesh-independent convergence can be achieved when employing the minimal
residual method. In particular, we prove that the eigenvalues of the
preconditioned matrix are bounded within disjoint intervals containing $\pm 1$,
without any outliers. Numerical examples are provided to critically discuss the
results, showcase the spectral distribution, and support the efficacy of our
preconditioning strategy.
