# Laplacian Filtered Loop-Star Decompositions and Quasi-Helmholtz Filters: Definitions, Analysis, and Efficient Algorithms

Adrien Merlini, Clément Henry, Davide Consoli, Lyes Rahmouni, Alexandre Dély, Francesco Andriulli

<figure><img src="image_1.png" /></figure>

## ► To cite this version:

Adrien Merlini, Clément Henry, Davide Consoli, Lyes Rahmouni, Alexandre Dély, et al., Laplacian Filtered Loop-Star Decompositions and Quasi-Helmholtz Filters: Definitions, Analysis, and Efficient Algorithms. IEEE Transactions on Antennas and Propagation, 2023, 71 (12), pp.9289 - 9302. 10.1109/TAP.2023.3283043. [hal-04124349](https://HAL Id: hal-04124349)

HAL Id: hal-04124349

https://imt-atlantique.hal.science/hal-04124349v1

Submitted on 4 Dec 2024

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_2.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
# Laplacian Filtered Loop-Star Decompositions and Quasi-Helmholtz Laplacian Filters: Definitions, Analysis, and Efficient Algorithms

Adrien Merlini, Member, IEEE, Clément Henry, Member, IEEE, Davide Consoli, Student Member, IEEE, Lyes Rahmouni, Alexandre Dély, and Francesco P. Andriulli, Senior Member, IEEE

*Abstract*—Quasi-Helmholtz decompositions are fundamental tools in integral equation modeling of electromagnetic problems because of their ability of rescaling solenoidal and non-solenoidal components of solutions, operator matrices, and radiated fields. These tools are however incapable, *per se*, of modifying the refinement-dependent spectral behavior of the different operators and often need to be combined with other preconditioning strategies. This paper introduces the new concept of filtered quasi-Helmholtz decompositions proposing them in two incarnations: the filtered Loop-Star functions and the quasi-Helmholtz Laplacian filters. Because they are capable of manipulating large parts of the operators' spectra, new families of preconditioners and fast solvers can be derived from these new tools. A first application to the case of the frequency and *h*-refinement preconditioning of the electric field integral equation is presented together with numerical results showing the practical effectiveness of the newly proposed decompositions.

*Index Terms*—Integral equations, quasi-Helmholtz decompositions, quasi-Helmholtz projectors, preconditioning, EFIE.

## I. INTRODUCTION

INTEGRAL equation formulations are effective numerical strategies for modeling radiation and scattering by perfectly electrically conducting objects [1]–[3]. Their effectiveness primarily derives from the fact that they only require the scatterers’ surfaces to be discretized, automatically impose radiation conditions and, thanks to the advent of fast algorithms [4], give rise to linear-in-complexity approaches when solved with iterative schemes—provided that the conditioning of the linear system matrices resulting from their discretizations is independent of the number of unknowns [5]. Among the well-established formulations, the electric field integral equation (EFIE) plays a crucial role, both in itself and within combined field formulations [6]. The EFIE, lamentably, becomes ill-conditioned when the frequency is low or the discretization density high [7]. These phenomena—respectively known as the low-frequency and *h*-refinement breakdowns—cause the solution of the EFIE to become increasingly challenging to obtain, as the number of iterations of the solution process grows unbounded, which jeopardizes the possibility of achieving an overall linear complexity.

D. Consoli, L. Rahmouni, A. Dély, and F. P. Andriulli are with the Department of Electronics and Telecommunications, Politecnico di Torino, 10129 Torino, Italy; e-mail: name.surname@polito.it.

A. Merlini and C. Henry are with the Microwave department, IMT Atlantique, 29238 Brest cedex 03, France; e-mail: name.surname@imt-atlantique.fr. Manuscript received April 19, 2005.

Traditional approaches to tackle the low-frequency breakdown rely on standard quasi-Helmholtz decompositions such as Loop-Star/Tree bases [8]–[12] that, despite curing the low-frequency behavior, worsen the *h*-refinement ill-conditioning of the EFIE [11] because of the derivative nature of the change of basis [12]. A way to circumvent the issue is the use of hierarchical strategies both on structured [13], [14] and unstructured meshes [15]–[18]. These schemes, when designed properly, can solve both the low-frequency and the *h*-refinement problems but still rely on the construction on an explicit, basis-based, quasi-Helmholtz decomposition that requires the cumbersome detection of topological loops whenever handles are present in the geometry. A popular alternative strategy leverages Calderón identities to form a second kind integral equation out of the EFIE. Calderón approaches concurrently solve the low-frequency and the *h*-refinement breakdowns without calling for an explicit quasi-Helmholtz decomposition [19]–[26]. In their standard incarnations they do, however, require the use of a dual discretization and global loop handling, because global loops reside in the static null-space of the Calderón operator. The introduction of implicit quasi-Helmholtz decompositions via the so called quasi-Helmholtz projectors [27], when combined with Calderón approaches, led to the design of several well-conditioned formulations, free from static nullspaces (see [7], [27]–[29] and references therein) and, in some incarnations, free from the need of performing a barycentric refinement [30]. Quasi-Helmholtz projectors have shown to be an effective and efficiently computable tool for performing quasi-Helmholtz decompositions, but, by themselves, they can only tackle the low-frequency breakdown and must be combined with Calderón-like strategies that involve multiple operators, to obtain *h*-refinement spectral preconditioning effects. A set of tools as versatile as the projectors that could also manipulate the operator spectra beyond a simple rescaling would thus be desirable.

This paper introduces such a new family of tools. The contribution of this work is in fact threefold: (i) we will introduce the concept of Laplacian-filtered Loop-Star decompositions, a new quasi-Helmholtz decomposition approach that will allow for a finer tuning of the operator spectrum with respect to their standard Loop-Star counterparts. (ii) Just like standard Loop-Star bases give rise to the quasi-Helmholtz projectors, a suitable choice of projections on the filtered Loop-Star spaces will give rise to a new family of mathematical objects,
the quasi-Helmholtz Laplacian filters, that will not require
handling basis functions and explicit decompositions, while
still providing the spectral tuning properties of (i). (iii) We
will obtain new frequency and h-refinement preconditioners
for the EFIE based on (i) and (ii) that represent a natural
first application of the newly proposed techniques. We believe,
however, that the applicability of the new spectral filters will
extend beyond EFIE preconditioning in further investigations.
The contribution will be further enriched by a section devoted
to efficient implementations of the newly defined tools that
will be obtained by leveraging strategies developed in the
context of polynomial preconditioning approaches and in
graph wavelet theory [31]–[34]. Numerical results will then
corroborate and confirm our theoretical considerations.

The paper is organized as follows: the background material
and the notation are presented in Section II, the new Laplacian
filtered Loop-Star decompositions and the quasi-Helmholtz
Laplacian filters are presented in Section III and Section IV,
respectively, along with their main properties. Various strate-
gies for computing the filters, in practical scenarios, are
detailed in Section V. Preconditioners tackling simultaneously
the low-frequency and h-refinement breakdowns of the EFIE
are then derived in Section VI. Finer details relating to the
implementation and computation of the filters and precon-
ditioners are then presented in Section VII. Illustrations of
the effectiveness of the schemes are provided in Section VIII,
before concluding in Section IX.

Preliminary results from this work were presented in the
conference contributions [35], [36].

II. NOTATION AND BACKGROUND

Let $\Gamma$ be a smooth surface modeling the boundary of a perfectly electrically conducting (PEC), closed scatterer enclosed in a homogeneous background medium with permittivity $\epsilon$ and permeability $\mu$. The boundary $\Gamma$ can be multiply connected and contain holes. We denote by $\hat{n}(\mathbf{r})$ the outward pointing normal field at $\mathbf{r}$. When illuminated by a time-harmonic incident electric field $\mathbf{E}^i$, a surface current density $\mathbf{J}$ is induced on $\Gamma$ that satisfies the electric field integral equation (EFIE)

$$
\mathcal{T} \mathbf{J} = \mathcal{T}_s \mathbf{J} + \mathcal{T}_h \mathbf{J} = -\hat{\mathbf{n}} \times \mathbf{E}^i, \quad (1)
$$

where

$$
\mathcal{T}_s \mathbf{J} = \hat{\mathbf{n}}(\mathbf{r}) \times i\kappa \int_{\Gamma} \frac{e^{ik\|\mathbf{r}-\mathbf{r}'\|}}{4\pi\|\mathbf{r}-\mathbf{r}'\|} \mathbf{J}(\mathbf{r}') \mathrm{d}S(\mathbf{r}'), \quad (2)
$$

$$
\mathcal{T}_h \mathbf{J} = -\hat{\mathbf{n}}(\mathbf{r}) \times \frac{1}{i\kappa} \nabla \int_{\Gamma} \frac{e^{ik\|\mathbf{r}-\mathbf{r}'\|}}{4\pi\|\mathbf{r}-\mathbf{r}'\|} \nabla' \cdot \mathbf{J}(\mathbf{r}') \mathrm{d}S(\mathbf{r}'), \quad (3)
$$

and $k$ is the wavenumber of the electromagnetic wave in
the background medium. Equation (1) can be solved nu-
merically by approximating $\Gamma$ with triangular elements of
average edge length $h$ and by approximating the current
density as $\mathbf{J} \approx \sum_{n=1}^{N} [\mathbf{j}]_n \mathbf{f}_n$ with the Rao-Wilton-Glisson
basis functions $\{\mathbf{f}_n\}_n$ [37], in which $N$ is the number of

edges in the mesh, $j$ is the vector of the coefficients of the
expansion, and $\mathbf{f}_n$ is defined as

$$
\mathbf{f}_n(\mathbf{r}) = \begin{cases} \frac{\mathbf{r} - \mathbf{r}_n^+}{2A_n^+} & \text{if } \mathbf{r} \in c_n^+ \\ -\frac{\mathbf{r} - \mathbf{r}_n^-}{2A_n^-} & \text{if } \mathbf{r} \in c_n^- \end{cases} \quad (4)
$$

where the notation of Fig. 1 was employed and where $A_n^{\pm}$ is
the area of the cell $c_n^{\pm}$.

The final step to obtain the discretized EFIE is to test (1)
with the rotated RWG functions {$\hat{\mathbf{n}} \times \mathbf{f}_n$}, which results in
the linear system

$$
\mathbf{T}j = (\mathbf{T}_s + \mathbf{T}_h)j = \mathbf{v}, \quad (5)
$$

in which $[\mathbf{T}_s]_{mn} = \langle \hat{\mathbf{n}} \times \mathbf{f}_m, \mathcal{T}_s(\mathbf{f}_n) \rangle$, $[\mathbf{T}_h]_{mn} = \langle \hat{\mathbf{n}} \times \mathbf{f}_m, \mathcal{T}_h(\mathbf{f}_n) \rangle$, $[\mathbf{v}]_m = \langle \hat{\mathbf{n}} \times \mathbf{f}_m, -\hat{\mathbf{n}} \times \mathbf{E}^i \rangle$, and $\langle \mathbf{a}, \mathbf{b} \rangle = \int_{\Gamma} \mathbf{a} \cdot \mathbf{bds}$. The EFIE can also be discretized on the dual mesh using dual functions defined on the barycentric refinement. Both Buffa-Christiansen [38] and Chen-Wilton [39] elements can be used for this dual discretization. For the sake of brevity, we will omit the explicit definitions of the dual elements that will be denoted by $\{g_n\}_n$ in the following; the reader can refer to [7] and references therein for a more detailed treatment. We will also need the definition of the standard and dual Gram matrices whose entries are $[\mathbf{G}]_{mn} = \langle \mathbf{f}_m, \mathbf{f}_n \rangle$ and $[\mathbb{G}]_{mn} = \langle \mathbf{g}_m, \mathbf{g}_n \rangle$. While they are not required for the discretization of the EFIE itself, we introduce the patch and pyramids scalar basis functions sets $\{p_n\}$ and $\{\lambda_n\}$, respectively composed of $N_S$ and $N_L$ functions, that will be required for some of the following developments. These functions are defined as

$$
p_m(\mathbf{r}) = \begin{cases} A_m^{-1} & \text{if } \mathbf{r} \in c_m, \\ 0 & \text{otherwise,} \end{cases} \tag{6}
$$

and

$$
\lambda_m(\mathbf{r}) = \begin{cases} 1 & \mathbf{r} = \mathbf{v}_m, \\ 0 & \mathbf{r} = \mathbf{v}_n, n \neq m, \\ \text{linear} & \text{otherwise,} \end{cases} \quad (7)
$$

where {$v_n$}$_n$ are the vertices of the mesh. The number of
these functions can be deduced from the mesh properties:
$N_S$ is the number of mesh triangles and $N_L$ is the number
of mesh vertices. The Gram matrices corresponding to these
bases are $G_p$ for the patch functions and $G_\lambda$ for the pyramids
with $[G_p]_{mn} = \langle p_m, p_n \rangle$ and $[G_\lambda]_{mn} = \langle \lambda_m, \lambda_n \rangle$. The dual
of these functions, living in the barycentric refinement of the
original mesh will also be required, and their definitions, omit-
ted here for conciseness, can be found in [30]. The $N_L$ dual
patches will be designated as {$\tilde{p}_n$}$_n$, the $N_S$ dual pyramids as
{$\tilde{\lambda}_n$}$_n$, and the corresponding gram matrices as $\mathbb{G}_{\tilde{p}}$ and $\mathbb{G}_{\tilde{\lambda}}$
with $[\mathbb{G}_{\tilde{p}}]_{mn} = \langle \tilde{p}_m, \tilde{p}_n \rangle$ and $[\mathbb{G}_{\tilde{\lambda}}]_{mn} = \langle \tilde{\lambda}_m, \tilde{\lambda}_n \rangle$.

Because this contribution deals with discrete quasi-
Helmholtz decompositions we will recall some of their prop-
erties. The continuous solution **J** can be decomposed into a
solenoidal, irrotational, and (on non simply-connected mani-
folds) harmonic components

$$
\mathbf{J} = \nabla \times \hat{\mathbf{n}}\lambda + \nabla_s \phi + \mathbf{h}. \quad (8)
$$
<figure><img src="image_3.png" /><figcaption>Fig. 1: Convention used for the RWGs: each function is defined on the two triangles $c_i^+$ and $c_i^-$ that are formed with their common edge $e_i$ and the vertices $r_i^+$ and $r_i^-$, respectively.</figcaption></figure>

When $J$ is discretized by the approximate expansion in RWG functions, a discrete counterpart of (8) holds for the coefficient vector $j$

$$j = \Lambda l + \Sigma \sigma + \mathbf{H} h, \quad (9)$$

where $\Lambda \in \mathbb{R}^{N \times N_L}$ and $\Sigma \in \mathbb{R}^{N \times N_S}$ are the Loop-to-RWG and Star-to-RWG transformation matrices [9], [12], [40]–[42] defined, following the convention of Fig. 1 and the definition of the RWG in (4), as

$$[\Lambda]_{mn} = \begin{cases} 1 & \text{if node } n \text{ equals } v_m^+, \\ -1 & \text{if node } n \text{ equals } v_m^-, \\ 0 & \text{otherwise,} \end{cases} \quad (10)$$

and

$$[\Sigma]_{mn} = \begin{cases} 1 & \text{if the cell } n \text{ equals } c_m^+, \\ -1 & \text{if the cell } n \text{ equals } c_m^-, \\ 0 & \text{otherwise.} \end{cases} \quad (11)$$

With these definitions $\Lambda^T \Lambda$ and $\Sigma^T \Sigma$ are respectively the vertices- and the cells-based graph Laplacians [12]. The explicit use of the change of basis matrix $\mathbf{H}$ will not be required and we omit here its explicit definition, for the sake of conciseness, which could however be found in [7] and references therein.

With the definitions above, the standard quasi-Helmholtz projectors [12], [27] are defined as

$$\begin{aligned} \mathbf{P}^{\Sigma} &= \Sigma (\Sigma^{\mathrm{T}} \Sigma)^{+} \Sigma^{\mathrm{T}}, \\ \mathbf{P}^{\Lambda H} &= \mathbf{I} - \mathbf{P}^{\Sigma} \end{aligned} \quad (12)$$

for the primal ones,

$$\begin{aligned} \mathbb{P}^{\Lambda} &= \Lambda (\Lambda^{\mathrm{T}} \Lambda)^{+} \Lambda^{\mathrm{T}}, \\ \mathbb{P}^{\Sigma H} &= \mathbf{I} - \mathbb{P}^{\Lambda} \end{aligned} \quad (13)$$

for the dual ones, and

$$\mathbf{P}^{H} (= \mathbb{P}^{H}) = \mathbf{I} - \mathbf{P}^{\Sigma} - \mathbb{P}^{\Lambda} \quad (14)$$

for the projector to quasi-harmonic subspace, where + denotes the Moore-Penrose pseudo-inverse.

### III. LAPLACIAN FILTERED LOOP-STAR DECOMPOSITIONS

In this section, we will extend the notion of Loop-Star bases by introducing the concept of filtered (generating) functions. We will first treat graph-based decompositions (a direct generalization of the standard case) and we will then move on to their Gram matrix normalized counterparts that will be more effective in treating problems involving inhomogeneous meshes.

#### A. The Standard Case

Consider the singular value decomposition (SVD) [43] of a matrix $\mathbf{X} \in \mathbb{R}^{N \times N_x}$

$$\mathbf{X} = \mathbf{U}_{\mathbf{X}} \mathbf{S}_{\mathbf{X}} \mathbf{V}_{\mathbf{X}}^{\mathrm{T}} \quad (15)$$

where $\mathbf{X}$ is a placeholder for either $\Sigma$ or $\Lambda$, $\mathbf{U}_{\mathbf{X}} \in \mathbb{R}^{N \times N}$, $\mathbf{V}_{\mathbf{X}} \in \mathbb{R}^{N_x \times N_x}$, and $\mathbf{S}_{\mathbf{X}} \in \mathbb{R}^{N \times N_x}$. The matrices $\mathbf{U}_{\mathbf{X}}$ and $\mathbf{V}_{\mathbf{X}}$ are unitary and $\mathbf{S}_{\mathbf{X}}$ is a block diagonal matrix with the singular values $\sigma_{\mathbf{X},i}$ of $\mathbf{X}$ as entries (in decreasing order). Clearly the SVD of $(\mathbf{X}^{\mathrm{T}}\mathbf{X})$ is $\mathbf{V}_{\mathbf{X}}\mathbf{S}_{\mathbf{X}}^{\mathrm{T}}\mathbf{S}_{\mathbf{X}}\mathbf{V}_{\mathbf{X}}^{\mathrm{T}}$, and, by defining the diagonal matrix $\mathbf{L}_{\mathbf{X},n} \in \mathbb{R}^{N_x \times N_x}$, with $1 \le n \le N_x$, such that

$$[\mathbf{L}_{\mathbf{X},n}]_{ii} = \begin{cases} \sigma_{\mathbf{X},i} & \text{if } i > N_x - n, \\ 0 & \text{otherwise,} \end{cases} \quad (16)$$

we define the filtered graph Laplacians

$$(\mathbf{X}^{\mathrm{T}}\mathbf{X})_n := \mathbf{V}_{\mathbf{X}} \mathbf{L}_{\mathbf{X},n}^2 \mathbf{V}_{\mathbf{X}}^{\mathrm{T}}, \quad (17)$$

from which we introduce the filtered Loop-to-RWG and filtered Star-to-RWG matrices we propose in this work

$$\Sigma_n = \Sigma (\Sigma^T \Sigma)^+ (\Sigma^T \Sigma)_n, \quad (18)$$

$$\Lambda_n = \Lambda (\Lambda^T \Lambda)^+ (\Lambda^T \Lambda)_n, \quad (19)$$

These matrices contains the coefficients of sets of linearly dependent filtered Loop-Star functions.

**Properties:** We now study some properties of the filtered Loop-Star matrices. Because $\Sigma^T \Lambda = \mathbf{0}$ [7], we have $\forall n, m$

$$\begin{aligned} \Sigma_n^T \Lambda_m &= (\Sigma^T \Sigma)_n (\Sigma^T \Sigma)^+ \Sigma^T \Lambda (\Lambda^T \Lambda)^+ (\Lambda^T \Lambda)_m \\ &= \mathbf{0}. \end{aligned} \quad (20)$$

Otherwise said, the new filtered Loop-Star functions are coefficient-orthogonal ($l^2$-orthogonal) like their non-filtered, standard counterparts.

From the definition of $\mathbf{L}_{\mathbf{X},n}$ in (16), it follows that $\mathbf{L}_{\mathbf{X},n}\mathbf{L}_{\mathbf{X},m} = \mathbf{L}_{\mathbf{X},\min\{n,m\}}^2$. Thus from (17)

$$\begin{aligned} (\mathbf{X}^{\mathrm{T}}\mathbf{X})_n (\mathbf{X}^{\mathrm{T}}\mathbf{X})_m &= \mathbf{V}_{\mathbf{X}}\mathbf{L}_{\mathbf{X},n}\mathbf{V}_{\mathbf{X}}^{\mathrm{T}}\mathbf{V}_{\mathbf{X}}\mathbf{L}_{\mathbf{X},m}\mathbf{V}_{\mathbf{X}}^{\mathrm{T}} \\ &= \mathbf{V}_{\mathbf{X}}\mathbf{L}_{\mathbf{X},n}\mathbf{L}_{\mathbf{X},m}\mathbf{V}_{\mathbf{X}}^{\mathrm{T}} = \mathbf{V}_{\mathbf{X}}\mathbf{L}_{\mathbf{X},\min\{n,m\}}\mathbf{L}_{\mathbf{X},\min\{n,m\}}\mathbf{V}_{\mathbf{X}}^{\mathrm{T}} \\ &= \mathbf{V}_{\mathbf{X}}\mathbf{L}_{\mathbf{X},\min\{n,m\}}\mathbf{V}_{\mathbf{X}}^{\mathrm{T}}\mathbf{V}_{\mathbf{X}}\mathbf{L}_{\mathbf{X},\min\{n,m\}}\mathbf{V}_{\mathbf{X}}^{\mathrm{T}} \\ &= (\mathbf{X}^{\mathrm{T}}\mathbf{X})_{\min\{n,m\}}^2. \end{aligned} \quad (21)$$

We thus have

$$\begin{aligned} \Sigma_m^T \Sigma_n &= (\Sigma^T \Sigma)_m (\Sigma^T \Sigma)^+ \Sigma^T \Sigma (\Sigma^T \Sigma)^+ (\Sigma^T \Sigma)_n \\ &= (\Sigma^T \Sigma)^+ (\Sigma^T \Sigma)_m (\Sigma^T \Sigma)_n \\ &= (\Sigma^T \Sigma)^+ (\Sigma^T \Sigma)_{\min\{n,m\}} (\Sigma^T \Sigma)_{\min\{n,m\}} \quad (22) \\ &= \Sigma_{\min\{n,m\}}^T \Sigma_{\min\{n,m\}}. \end{aligned}$$

Similarly,

$$\Lambda_m^T \Lambda_n = \Lambda_{\min\{n,m\}}^T \Lambda_{\min\{n,m\}}. \quad (23)$$

Given integers such that $m < n < p < q$, the property

$$\begin{aligned} (\Sigma_m - \Sigma_n)^{\mathrm{T}} (\Sigma_p - \Sigma_q) &= \Sigma_m^{\mathrm{T}} (\Sigma_p - \Sigma_q) - \Sigma_n^{\mathrm{T}} (\Sigma_p - \Sigma_q) = \mathbf{0}, \end{aligned} \quad (24)$$
holds and, similarly,

$$ (\mathbf{\Lambda}_m - \mathbf{\Lambda}_n)^{\mathrm{T}} (\mathbf{\Lambda}_p - \mathbf{\Lambda}_q) = \mathbf{0}. \qquad (25) $$

Properties (24) and (25) show that non-intersecting differences of filtered Star or Loop bases are mutually orthogonal (and thus generate linearly independent spaces), a property that will be useful to build invertible changes of basis, as will be shown in Section VI-A.

## B. Generalization for Non-homogeneously Meshed Geometries

When the filtered Loop-Star decompositions are to be used on geometries with non-homogenous discretizations, both the standard discretizations of the EFIE and the graph Laplacian matrices may lead to suboptimal performance and a proper normalization with Gram matrices must be employed. In this context, we define the normalized EFIE electromagnetic operator matrices

$$ \tilde{\mathbf{T}} = \mathbf{G}^{-1/2} \mathbf{T} \mathbf{G}^{-1/2}, \qquad (26) $$

$$ \tilde{\mathbf{T}}_s = \mathbf{G}^{-1/2} \mathbf{T}_s \mathbf{G}^{-1/2}, \qquad (27) $$

$$ \tilde{\mathbf{T}}_h = \mathbf{G}^{-1/2} \mathbf{T}_h \mathbf{G}^{-1/2}, \qquad (28) $$

and the normalized Loop and Star matrices

$$ \tilde{\mathbf{\Sigma}} = \mathbf{G}^{-1/2} \mathbf{\Sigma} \mathbf{G}_{p}^{1/2}, \qquad (29) $$

$$ \tilde{\mathbf{\Lambda}} = \mathbf{G}^{1/2} \mathbf{\Lambda} \mathbf{G}_{\lambda}^{-1/2}. \qquad (30) $$

Following the same strategy as in (18) and (19), the normalized filtered Loop-Star matrices are consistently defined as

$$ \tilde{\mathbf{\Sigma}}_n = \tilde{\mathbf{\Sigma}} \left( \tilde{\mathbf{\Sigma}}^{\mathrm{T}} \tilde{\mathbf{\Sigma}} \right)^{+} \left( \tilde{\mathbf{\Sigma}}^{\mathrm{T}} \tilde{\mathbf{\Sigma}} \right)_n, \qquad (31) $$

$$ \tilde{\mathbf{\Lambda}}_n = \tilde{\mathbf{\Lambda}} \left( \tilde{\mathbf{\Lambda}}^{\mathrm{T}} \tilde{\mathbf{\Lambda}} \right)^{+} \left( \tilde{\mathbf{\Lambda}}^{\mathrm{T}} \tilde{\mathbf{\Lambda}} \right)_n. \qquad (32) $$

When dealing with dual Loop-Star decomposition matrices, the normalization is different from that of the primal ones, and the dually-normalized Loop and Star transformation matrices are defined as

$$ \tilde{\mathbf{\Sigma}} = \mathbb{G}^{1/2} \mathbf{\Sigma} \mathbb{G}_{\lambda}^{-1/2}, \qquad (33) $$

$$ \tilde{\mathbf{\Lambda}} = \mathbb{G}^{-1/2} \mathbf{\Lambda} \mathbb{G}_{p}^{1/2}, \qquad (34) $$

and the associated filtered decomposition matrices as

$$ \tilde{\mathbf{\Sigma}}_n = \tilde{\mathbf{\Sigma}} \left( \tilde{\mathbf{\Sigma}}^{\mathrm{T}} \tilde{\mathbf{\Sigma}} \right)^{+} \left( \tilde{\mathbf{\Sigma}}^{\mathrm{T}} \tilde{\mathbf{\Sigma}} \right)_n, \qquad (35) $$

$$ \tilde{\mathbf{\Lambda}}_n = \tilde{\mathbf{\Lambda}} \left( \tilde{\mathbf{\Lambda}}^{\mathrm{T}} \tilde{\mathbf{\Lambda}} \right)^{+} \left( \tilde{\mathbf{\Lambda}}^{\mathrm{T}} \tilde{\mathbf{\Lambda}} \right)_n. \qquad (36) $$

**Properties:** The primal and dual normalized Loop-Star bases keep satisfying the orthogonality properties

$$ \tilde{\mathbf{\Lambda}}^{\mathrm{T}} \tilde{\mathbf{\Sigma}} = \mathbf{G}_{\lambda}^{-1/2} \mathbf{\Lambda}^{\mathrm{T}} \mathbf{G}^{1/2} \mathbf{G}^{-1/2} \mathbf{\Sigma} \mathbf{G}_{p}^{1/2} = \mathbf{0}, \qquad (37) $$

$$ \tilde{\mathbf{\Sigma}}^{\mathrm{T}} \tilde{\mathbf{\Lambda}} = \mathbb{G}_{p}^{1/2} \mathbf{\Sigma}^{\mathrm{T}} \mathbb{G}^{-1/2} \mathbb{G}_{p}^{1/2} \mathbf{\Lambda} \mathbb{G}_{\lambda}^{-1/2} = \mathbf{0}, \qquad (38) $$

because $\mathbf{\Sigma}^{\mathrm{T}}\mathbf{\Lambda} = \mathbf{0}$. Moreover, because (21) holds, we obtain, similarly to Section III-A, that

$$ \tilde{\mathbf{\Sigma}}_{m}^{\mathrm{T}} \tilde{\mathbf{\Sigma}}_{n} = \tilde{\mathbf{\Sigma}}_{\mathrm{min}\{n,m\}}^{\mathrm{T}} \tilde{\mathbf{\Sigma}}_{\mathrm{min}\{n,m\}}, \qquad (39) $$

$$ \tilde{\mathbf{\Lambda}}_{m}^{\mathrm{T}} \tilde{\mathbf{\Lambda}}_{n} = \tilde{\mathbf{\Lambda}}_{\mathrm{min}\{n,m\}}^{\mathrm{T}} \tilde{\mathbf{\Lambda}}_{\mathrm{min}\{n,m\}}, \qquad (40) $$

$$ \tilde{\mathbf{\Lambda}}_{m}^{\mathrm{T}} \tilde{\mathbf{\Lambda}}_{n} = \tilde{\mathbf{\Lambda}}_{\mathrm{min}\{n,m\}}^{\mathrm{T}} \tilde{\mathbf{\Lambda}}_{\mathrm{min}\{n,m\}}, \qquad (41) $$

$$ \tilde{\mathbf{\Sigma}}_{m}^{\mathrm{T}} \tilde{\mathbf{\Sigma}}_{n} = \tilde{\mathbf{\Sigma}}_{\mathrm{min}\{n,m\}}^{\mathrm{T}} \tilde{\mathbf{\Sigma}}_{\mathrm{min}\{n,m\}}, \qquad (42) $$

Using these properties and the same reasoning as previously, the counterparts of the properties of the non-normalized filtered Loop-Star matrices can be obtained. In particular the counterparts of (24) and (25) can be obtained by replacing each matrix with its normalized (“tilde”) counterpart.

## IV. QUASI-HELMHOLTZ LAPLACIAN FILTERS

Although explicit quasi-Helmholtz decomposition bases are useful in applications in which a direct access to the Helmholtz components of the current is required, oftentimes, especially when the main target is preconditioning and regularization, implicit Helmholtz decompositions can be more efficient. An implicit Helmholtz decomposition was obtained in [27], where the concept of quasi-Helmholtz projector was introduced. Following a similar philosophy, and leveraging the filtered Loop-Star functions introduced above, we can now define quasi-Helmholtz Laplacian filters.

### A. The Standard Case

The idea behind the projectors was to obtain a basis-free quasi-Helmholtz decomposition that would not worsen the conditioning of the original equation. If we follow the definitions (12) and (13) by replacing the standard Star basis with the new filtered sets, we obtain

$$ \mathbf{\Sigma}_n \left( \mathbf{\Sigma}_n^\mathrm{T} \mathbf{\Sigma}_n \right)^{+} \mathbf{\Sigma}_n^\mathrm{T} = \mathbf{\Sigma} \left( \mathbf{\Sigma}^\mathrm{T} \mathbf{\Sigma} \right)^{+} \left( \mathbf{\Sigma}^\mathrm{T} \mathbf{\Sigma} \right)_n \left( \left( \mathbf{\Sigma}^\mathrm{T} \mathbf{\Sigma} \right)_n \left( \mathbf{\Sigma}^\mathrm{T} \mathbf{\Sigma} \right)^{+} \mathbf{\Sigma}^\mathrm{T} \left( \mathbf{\Sigma}^\mathrm{T} \mathbf{\Sigma} \right)^{+} \left( \mathbf{\Sigma}^\mathrm{T} \mathbf{\Sigma} \right)_n \right)^{+} \left( \mathbf{\Sigma}^\mathrm{T} \mathbf{\Sigma} \right)_n \left( \mathbf{\Sigma}^\mathrm{T} \mathbf{\Sigma} \right)^{+} \mathbf{\Sigma}^\mathrm{T} = \mathbf{\Sigma} \left( \mathbf{\Sigma}^\mathrm{T} \mathbf{\Sigma} \right)_n^{+} \mathbf{\Sigma}^\mathrm{T}, \qquad (43) $$

and, similarly,

$$ \mathbf{\Lambda}_n \left( \mathbf{\Lambda}_n^\mathrm{T} \mathbf{\Lambda}_n \right)^{+} \mathbf{\Lambda}_n^\mathrm{T} = \mathbf{\Lambda} \left( \mathbf{\Lambda}^\mathrm{T} \mathbf{\Lambda} \right)_n^{+} \mathbf{\Lambda}^\mathrm{T}. \qquad (44) $$

This justifies the following definitions of the new primal filters

$$ \mathbf{P}_n^{\Sigma} = \mathbf{\Sigma} \left( \mathbf{\Sigma}^{\mathrm{T}} \mathbf{\Sigma} \right)_n^{+} \mathbf{\Sigma}^{\mathrm{T}}, \qquad (45) $$

$$ \mathbf{P}_n^{\Lambda H} = \mathbf{\Lambda} \left( \mathbf{\Lambda}^{\mathrm{T}} \mathbf{\Lambda} \right)_n^{+} \mathbf{\Lambda}^{\mathrm{T}} + \mathbf{I} - \mathbf{P}^{\Sigma} - \mathbb{P}^{\Lambda}, \qquad (46) $$

and dual filters

$$ \mathbb{P}_n^{\Lambda} = \mathbf{\Lambda} \left( \mathbf{\Lambda}^{\mathrm{T}} \mathbf{\Lambda} \right)_n^{+} \mathbf{\Lambda}^{\mathrm{T}}, \qquad (47) $$

$$ \mathbb{P}_n^{\Sigma H} = \mathbf{\Sigma} \left( \mathbf{\Sigma}^{\mathrm{T}} \mathbf{\Sigma} \right)_n^{+} \mathbf{\Sigma}^{\mathrm{T}} + \mathbf{I} - \mathbb{P}^{\Lambda} - \mathbf{P}^{\Sigma}. \qquad (48) $$

The reader should note that, in the special case of simply connected geometries $\mathbf{P}_n^{\Lambda H} = \mathbf{\Lambda} \left( \mathbf{\Lambda}^{\mathrm{T}} \mathbf{\Lambda} \right)_n^{+} \mathbf{\Lambda}^{\mathrm{T}}$ and $\mathbb{P}_n^{\Sigma H} = \mathbf{\Sigma} \left( \mathbf{\Sigma}^{\mathrm{T}} \mathbf{\Sigma} \right)_n^{+} \mathbf{\Sigma}^{\mathrm{T}}$ since $\mathbb{P}_n^{\Lambda} + \mathbb{P}_n^{\Sigma} = \mathbf{I}$. Moreover, by construction,

$$ \mathbf{P}_{N_S}^{\Sigma} = \mathbf{P}^{\Sigma}, \qquad (49) $$

$$ \mathbb{P}_{N_L}^{\Lambda} = \mathbb{P}^{\Lambda}, \qquad (50) $$
and thus

$$ \mathbf{P}_{N_L}^{\Lambda H} = \mathbb{P}^{\Lambda} + \mathbf{I} - \mathbf{P}^{\Sigma} - \mathbb{P}^{\Lambda} = \mathbf{P}^{\Lambda H}, \quad (51) $$

$$ \mathbb{P}_{N_S}^{\Sigma H} = \mathbf{P}^{\Sigma} + \mathbf{I} - \mathbb{P}^{\Lambda} - \mathbf{P}^{\Sigma} = \mathbb{P}^{\Sigma H}, \quad (52) $$

which means that, with these definitions, the quasi-Helmholtz Laplacian filters converge to the standard quasi-Helmholtz projectors when the Laplacian is unfiltered ($n = N_X$).

**Properties:** From these definitions, a few useful properties of the quasi-Helmholtz Laplacian filters can be derived and will be summarized here. First, the filters still behave as projectors since

$$ \begin{aligned} \mathbf{P}_n^{\Sigma} \mathbf{P}_n^{\Sigma} &= \mathbf{\Sigma} (\mathbf{\Sigma}^T \mathbf{\Sigma})_n^+ \mathbf{\Sigma}^T \mathbf{\Sigma} (\mathbf{\Sigma}^T \mathbf{\Sigma})_n^+ \mathbf{\Sigma}^T \\ &= \mathbf{\Sigma} (\mathbf{\Sigma}^T \mathbf{\Sigma})_n^+ \mathbf{\Sigma}^T = \mathbf{P}_n^{\Sigma}, \end{aligned} \quad (53) $$

$$ \begin{aligned} \mathbb{P}_n^\Lambda \mathbb{P}_n^\Lambda &= \mathbf{\Lambda} (\mathbf{\Lambda}^T \mathbf{\Lambda})_n^+ \mathbf{\Lambda}^T \mathbf{\Lambda} (\mathbf{\Lambda}^T \mathbf{\Lambda})_n^+ \mathbf{\Lambda}^T \\ &= \mathbf{\Lambda} (\mathbf{\Lambda}^T \mathbf{\Lambda})_n^+ \mathbf{\Lambda}^T = \mathbb{P}_n^\Lambda, \end{aligned} \quad (54) $$

and, similarly,

$$ \mathbf{P}_n^{\Lambda H} \mathbf{P}_n^{\Lambda H} = \mathbf{P}_n^{\Lambda H}, \quad (55) $$

$$ \mathbb{P}_n^{\Sigma H} \mathbb{P}_n^{\Sigma H} = \mathbb{P}_n^{\Sigma H}. \quad (56) $$

Moreover, $\forall m, n$

$$ \mathbf{P}_m^{\Sigma} \mathbf{P}_n^{\Lambda H} = \mathbf{\Sigma} (\mathbf{\Sigma}^T \mathbf{\Sigma})_n^+ \mathbf{\Sigma}^T \mathbf{\Lambda} (\mathbf{\Lambda}^T \mathbf{\Lambda})_n^+ \mathbf{\Lambda}^T \\ + \mathbf{\Sigma} (\mathbf{\Sigma}^T \mathbf{\Sigma})_n^+ \mathbf{\Sigma}^T (\mathbf{I} - \mathbf{P}^{\Sigma} - \mathbb{P}^{\Lambda}) = \mathbf{0}, \quad (57) $$

where the properties $\mathbf{\Sigma}^T \mathbf{\Lambda} = \mathbf{0}$ and $\mathbf{\Sigma}^T (\mathbf{I} - \mathbf{P}^{\Sigma} - \mathbb{P}^{\Lambda}) = \mathbf{0}$ have been used. A similar property and proof hold for the dual projectors

$$ \mathbb{P}_m^\Lambda \mathbb{P}_n^{\Sigma H} = \mathbf{0}, \quad \forall m, n. \quad (58) $$

For integers $m < n < p < q$, we have the following orthogonality property

$$ \begin{aligned} & (\mathbf{P}_m^\Sigma - \mathbf{P}_n^\Sigma) (\mathbf{P}_p^\Sigma - \mathbf{P}_q^\Sigma) \\ &= (\mathbf{\Sigma} (\mathbf{\Sigma}^T \mathbf{\Sigma})_m^+ \mathbf{\Sigma}^T - \mathbf{\Sigma} (\mathbf{\Sigma}^T \mathbf{\Sigma})_n^+ \mathbf{\Sigma}^T) \\ & \quad (\mathbf{\Sigma} (\mathbf{\Sigma}^T \mathbf{\Sigma})_p^+ \mathbf{\Sigma}^T - \mathbf{\Sigma} (\mathbf{\Sigma}^T \mathbf{\Sigma})_q^+ \mathbf{\Sigma}^T) \\ &= \mathbf{\Sigma} (\mathbf{\Sigma}^T \mathbf{\Sigma})_m^+ \mathbf{\Sigma}^T - \mathbf{\Sigma} (\mathbf{\Sigma}^T \mathbf{\Sigma})_m^+ \mathbf{\Sigma}^T \\ & \quad + \mathbf{\Sigma} (\mathbf{\Sigma}^T \mathbf{\Sigma})_n^+ \mathbf{\Sigma}^T - \mathbf{\Sigma} (\mathbf{\Sigma}^T \mathbf{\Sigma})_n^+ \mathbf{\Sigma}^T = \mathbf{0}, \end{aligned} \quad (59) $$

where (21) has been used. In a similar way, one can prove that

$$ (\mathbb{P}_m^\Lambda - \mathbb{P}_n^\Lambda) (\mathbb{P}_p^\Lambda - \mathbb{P}_q^\Lambda) = \mathbf{0}. \quad (60) $$

Moreover, given that $\mathbf{P}_n^{\Lambda H} - \mathbf{P}_m^{\Lambda H} = \mathbb{P}_n^{\Lambda} - \mathbb{P}_m^{\Lambda}$ and $\mathbb{P}_n^{\Sigma H} - \mathbb{P}_m^{\Sigma H} = \mathbb{P}_n^{\Sigma} - \mathbb{P}_m^{\Sigma}$ $\forall n, m$—which can be deduced from (46) and (48)—the remaining properties

$$ (\mathbf{P}_m^{\Lambda H} - \mathbf{P}_n^{\Lambda H}) (\mathbf{P}_p^{\Lambda H} - \mathbf{P}_q^{\Lambda H}) = \mathbf{0}, \quad (61) $$

$$ (\mathbb{P}_m^{\Sigma H} - \mathbb{P}_n^{\Sigma H}) (\mathbb{P}_p^{\Sigma H} - \mathbb{P}_q^{\Sigma H}) = \mathbf{0} \quad (62) $$

follow. All the properties listed above, will be useful when building invertible transforms, similarly to their basis-based counterpart (24).

## B. Generalization for Non-homogeneously Meshed Geometries

The definitions of the normalized Loop and Star matrices in (29) and (30) suggest the following definition for the associated normalized quasi-Helmholtz projectors

$$ \tilde{\mathbf{P}}^{\Sigma} = \tilde{\mathbf{\Sigma}} (\tilde{\mathbf{\Sigma}}^{T} \tilde{\mathbf{\Sigma}})_{n}^{+} \tilde{\mathbf{\Sigma}}^{T}, \quad (63) $$

$$ \tilde{\mathbf{P}}^{\Lambda} = \tilde{\mathbf{\Lambda}} (\tilde{\mathbf{\Lambda}}^{T} \tilde{\mathbf{\Lambda}})_{n}^{+} \tilde{\mathbf{\Lambda}}^{T}. \quad (64) $$

Moreover, as is proved in Appendix A, the complementarity property

$$ \tilde{\mathbf{P}}^{\Sigma} = \mathbf{I} - \tilde{\mathbf{P}}^{\Lambda} \quad (65) $$

holds on simply connected geometries; on general geometries and together with definitions (31) and (32), this justifies the following definition for the normalized quasi-Helmholtz Laplacian filters

$$ \tilde{\mathbf{P}}_n^{\Sigma} = \tilde{\mathbf{\Sigma}} (\tilde{\mathbf{\Sigma}}^T \tilde{\mathbf{\Sigma}})_n^+ \tilde{\mathbf{\Sigma}}^T, \quad (66) $$

$$ \tilde{\mathbf{P}}_n^{\Lambda H} = \tilde{\mathbf{\Lambda}} (\tilde{\mathbf{\Lambda}}^T \tilde{\mathbf{\Lambda}})_n^+ \tilde{\mathbf{\Lambda}}^T + \mathbf{I} - \tilde{\mathbf{P}}^{\Sigma} - \tilde{\mathbf{P}}^{\Lambda}. \quad (67) $$

By analogy, we can define the normalized dual quasi-Helmholtz projectors as

$$ \tilde{\mathbf{P}}^{\Lambda} = \tilde{\mathbf{\Lambda}} (\tilde{\mathbf{\Lambda}}^{T} \tilde{\mathbf{\Lambda}})_{n}^{+} \tilde{\mathbf{\Lambda}}^{T}, \quad (68) $$

$$ \tilde{\mathbf{P}}^{\Sigma} = \tilde{\mathbf{\Sigma}} (\tilde{\mathbf{\Sigma}}^T \tilde{\mathbf{\Sigma}})_{n}^{+} \tilde{\mathbf{\Sigma}}^T, \quad (69) $$

with the property

$$ \tilde{\mathbf{P}}^{\Lambda} = \mathbf{I} - \tilde{\mathbf{P}}^{\Sigma} \quad (70) $$

holding on simply connected geometries (see Appendix A for the proof). Thus, dually to the primal case, we define, on general geometries,

$$ \tilde{\mathbf{P}}_n^\Lambda = \tilde{\mathbf{\Lambda}} (\tilde{\mathbf{\Lambda}}^T \tilde{\mathbf{\Lambda}})_n^+ \tilde{\mathbf{\Lambda}}^T, \quad (71) $$

$$ \tilde{\mathbf{P}}_n^{\Sigma H} = \tilde{\mathbf{\Sigma}} (\tilde{\mathbf{\Sigma}}^T \tilde{\mathbf{\Sigma}})_n^+ \tilde{\mathbf{\Sigma}}^T + \mathbf{I} - \tilde{\mathbf{P}}^{\Lambda} - \tilde{\mathbf{P}}^{\Sigma}. \quad (72) $$

**Properties:** Since the primal and dual normalized Loop-Star bases still satisfy the orthogonality properties (37) and (38) and because of the properties (39)-(42), the same reasoning yields all counterparts of the properties (53)-(62), after replacing each matrix with its normalized (“tilde”) counterpart.

## V. EFFICIENT FILTERING ALGORITHMS

The definitions of the filtered Loop and Star functions and of the quasi-Helmholtz Laplacian filters in Sections III and IV involve an SVD which, while ensuring a clear and compact theoretical treatment, is in general computationally inefficient. This section will be devoted to presenting algorithms allowing for SVD-free matrix-vector products for the filtered graph Laplacians $(\mathbf{\Sigma}^T \mathbf{\Sigma})_n$ and $(\mathbf{\Lambda}^T \mathbf{\Lambda})_n$, which are the two key operations on both approaches presented in the previous section. Moreover, while our treatment will deal with the graph matrices $\mathbf{\Sigma}$ and $\mathbf{\Lambda}$, it is intended that substantially the same strategies can be applied when replacing those matrices with their normalized counterparts $\tilde{\mathbf{\Sigma}}$ and $\tilde{\mathbf{\Lambda}}$ with
minor modifications. In fact the additional products with the inverse square roots of (well-conditioned) Gram matrices can be obtained efficiently by using matrix function strategies [44].

### A. Power Method Filtering

For filters with a filtering index that is independent on the total number of degrees of freedom, preconditioned inverse power methods [43] yield the last singular vectors and singular values of the $\Sigma^T\Sigma$ and $\Lambda^T\Lambda$ matrices at the price of a constant number of matrix-vector products. Given that the matrices involved are sparse, the resulting method is linear in complexity and the filtered projectors can be efficiently obtained. These schemes are well known and we do not provide extensive details here for the sake of brevity. We just mention that a special care should be put into using this schemes in the presence of degenerate spectra (arising from symmetries for example); under these conditions, the scheme presented in the following section should rather be preferred.

### B. Butterworth Matrix Filters

An alternative strategy for the previous scenario, i.e. when a filter is needed that has filtering index which is independent on the total number of degrees of freedom, is provided by a matrix function and filtering approach. Given a scalar (squared) Butterworth filter of positive order $m$ and cutoff parameter $x_c > 0$, characterized by

$$f_{m,x_c}(x) = (1 + (x/x_c)^m)^{-1}, \quad x \ge 0, \quad (73)$$

the spectrum of a symmetric positive matrix $A \in \mathbb{R}^{N \times N}$ composed of the set of singular values $\{\sigma_i(A)\}_i$ can be filtered by generalizing $f_{m,x_c}$ to matrix arguments and applying it to $A$, yielding the filtered matrix

$$A_{\text{filt}} := f_{m,x_c}(A) = (I + (A/x_c)^m)^{-1}, \quad (74)$$

with singular values $\{f_{m,x_c}(\sigma_i(A))\}_i$. The filtered matrix $(\Sigma^T\Sigma)_n$ can now be expressed as

$$(\Sigma^T\Sigma)_n = (\Sigma^T\Sigma) \lim_{m \to \infty} f_{m,\sigma_n(\Sigma^T\Sigma)} (\Sigma^T\Sigma) . \quad (75)$$

The presence of high exponents in (75) may render its computation unstable. Hence we propose to use the following factorization formula that leverages the roots of unity

$$(\Sigma^T\Sigma)_n = (\Sigma^T\Sigma) \lim_{m \to \infty} \prod_{k=1}^{m} \left( \frac{\Sigma^T\Sigma}{\sigma_n(\Sigma^T\Sigma)} - e^{(2k+1)i\pi/N} \right)^{-1} . \quad (76)$$

For practical purposes the infinite products in this expression can be truncated at the desired precision. Regarding the value of $\sigma_n(\Sigma^T\Sigma)$, an approximation can be obtained either with ad-hoc heuristics or by the approximation $\sigma_n(\Sigma^T\Sigma) \approx (N_s - n)/\| (\Sigma^T\Sigma)^+ \|$. Finally, when the filtering point is a constant with respect to the number of unknowns, a multigrid approach is effective in providing the inverse required by (76).

### C. Filter Approximation via Chebyshev Polynomials

When the filtering index is proportional to the number of unknowns, the computational burden of the two methods above can become high. In this regime we can leverage the ideas of polynomial preconditioning and graph wavelets [31]–[34] and adopt a method based on a polynomial expansion of the spectral filter.

Because we are interested in cases in which the filtering index is proportional to the number of degrees of freedom (for instance, $n = N_s/2$) we can leverage a polynomial approximation of $f_{m,x_c}$ on the interval $[0, \sigma_{N_s}(\Sigma^T\Sigma)]$; a natural basis for this approximation is that of the Chebyshev polynomials $\{T_n(x)\}_n$, defined by the recurrence relation

$$T_n(x) = \begin{cases} 1 & \text{if } n = 0 \\ x & \text{if } n = 1 \\ 2xT_{n-1}(x) - T_{n-2}(x) & \text{otherwise.} \end{cases} \quad (77)$$

The approximated filtered matrix now reads

$$(\Sigma^T\Sigma)_n \approx -\frac{c_0}{2}I + \sum_{k=1}^{n_c} c_k T_k \left( \frac{\Sigma^T\Sigma}{\sigma_n(\Sigma^T\Sigma)} \right), \quad (78)$$

where the $c_n$ are the expansion coefficients of $f_{m,\sigma_n(\Sigma^T\Sigma)}$ in the basis of the first $n_c+1$ Chebyshev polynomials. Algorithms for their computation can be found, among others, in [45]. Because the cutoff frequency of this filter is proportional to the number of unknowns and so is the domain size, the order of the polynomial that is required to obtain a given approximation of the Butterworth filter, will not need to be changed with increasing discretizations. In other words, the filters obtained by following this approach will require the same number of sparse matrix-vector multiplication for increasing discretization when the filtering index will be proportional to the number of degrees of freedom. It should be noted that in the transition region between the filters described in the previous two sections (constant filtering index) and the scenario described here (filtering index will be proportional to the number of degrees of freedom) the Chebyshev approach decreases in efficiency and further treatments may be required [34].

## VI. A FIRST APPLICATION CASE SCENARIO: LAPLACIAN FILTER BASED PRECONDITIONING

As a first application case scenario of the new filters introduced here, we will develop two families of preconditioners for the EFIE in (5). This equation is known to suffer from ill-conditioning both for decreasing frequency and average mesh length $h$ (phenomena known as the low-frequency and $h$-refinement breakdowns, respectively, see [7] and references therein). In the following we will cure both breakdowns by developing preconditioners based both on filtered functions decompositions and on quasi-Helmholtz Laplacian filters.

The reader should note that in this Section and in the subsequent ones, we will study the singular value spectrum of potentially singular matrices. When dealing with such matrices, the condition number will be defined as $\text{cond}(A) = \|A\|\|A^+\|$. Moreover, inverse powers of singular matrices in the following will always denote the corresponding positive power of the pseudoinverse of the matrix.
A. Filtered Bases Approach

The primal and dual Laplacians can be used to precondition
the single layer and the hypersingular operator [1], [18], [46],
[47], thus $V_\Lambda$, and $V_\Sigma$ followed by a diagonal precondition-
ing are valid bases for regularizing the vector and scalar poten-
tial parts of the EFIE. In particular, for $T_h$, this results from
the fact that an operator spectrally equivalent to the single layer
can be obtained from $T_h$. In fact, noticing that $T_h = \Sigma R \Sigma^T$
[48], where $R$ is the patch-function discretization of the single
layer operator, i.e. $[R]_{mn} = \langle p_m, S p_n \rangle$ with

$$
(Sp)(\mathbf{r}) := \int_{\Gamma} \frac{e^{ik\|\mathbf{r}-\mathbf{r}'\|}}{4\pi \|\mathbf{r}-\mathbf{r}'\|} p(\mathbf{r}') dS(\mathbf{r}'), \quad (79)
$$

and defining $\tilde{\mathbf{R}} := \mathbf{G}_p^{-1/2} \mathbf{R} \mathbf{G}_p^{-1/2}$, we obtain $\tilde{\mathbf{T}}_h = \tilde{\mathbf{\Sigma}} \tilde{\mathbf{R}} \tilde{\mathbf{\Sigma}}^\mathrm{T}$.
The equivalence between $(\tilde{\mathbf{\Sigma}}^\mathrm{T} \tilde{\mathbf{\Sigma}})^+ \tilde{\mathbf{\Sigma}}^\mathrm{T} \tilde{\mathbf{T}}_h \tilde{\mathbf{\Sigma}} (\tilde{\mathbf{\Sigma}}^\mathrm{T} \tilde{\mathbf{\Sigma}})^+ $ and
$\tilde{\mathbf{R}}$ thus follows. To conclude the reasoning, we note that,
because

$$
\begin{align}
& (\tilde{\Sigma}^{\mathrm{T}} \tilde{\Sigma})^{1/4} \tilde{\mathbf{R}} (\tilde{\Sigma}^{\mathrm{T}} \tilde{\Sigma})^{1/4} \nonumber \\
& \quad = \tilde{V}_{\tilde{\Sigma}} (\tilde{S}_{\tilde{\Sigma}}^{\mathrm{T}} \tilde{S}_{\tilde{\Sigma}})^{1/4} \tilde{V}_{\tilde{\Sigma}}^{\mathrm{T}} \tilde{R} \tilde{V}_{\tilde{\Sigma}} (\tilde{S}_{\tilde{\Sigma}}^{\mathrm{T}} \tilde{S}_{\tilde{\Sigma}})^{1/4} \tilde{V}_{\tilde{\Sigma}}^{\mathrm{T}}, \tag{80}
\end{align}
$$

is well conditioned for increasing discretization—as a consequence of the results proven in [47], since $\tilde{\Sigma}^T \tilde{\Sigma}$ is a valid discretization of a Laplacian matrix [49]—we have

$$
\begin{align}
& \operatorname{cond} \left( \tilde{V}_{\tilde{\Sigma}} (\tilde{S}_{\tilde{\Sigma}}^{\mathrm{T}} \tilde{S}_{\tilde{\Sigma}})^{1/4} \tilde{V}_{\tilde{\Sigma}}^{\mathrm{T}} (\tilde{\Sigma}^{\mathrm{T}} \tilde{\Sigma})^{+} \tilde{\Sigma}^{\mathrm{T}} \tilde{T}_{\mathrm{h}} \right. \nonumber \\
& \qquad \left. \tilde{\Sigma} (\tilde{\Sigma}^{\mathrm{T}} \tilde{\Sigma})^{+} \tilde{V}_{\tilde{\Sigma}} (\tilde{S}_{\tilde{\Sigma}}^{\mathrm{T}} \tilde{S}_{\tilde{\Sigma}})^{1/4} \tilde{V}_{\tilde{\Sigma}}^{\mathrm{T}} \right) = O(1), h \to 0. \tag{81}
\end{align}
$$

The reader should note that, since $\tilde{V}_{\tilde{\Sigma}}$ is unitary, we also have

$$
\begin{align}
& \operatorname{cond} \left( \tilde{V}_{\tilde{\Sigma}} (\tilde{S}_{\tilde{\Sigma}}^{\mathrm{T}} \tilde{S}_{\tilde{\Sigma}})^{1/4} \tilde{V}_{\tilde{\Sigma}}^{\mathrm{T}} (\tilde{\Sigma}^{\mathrm{T}} \tilde{\Sigma})^{+} \tilde{\Sigma}^{\mathrm{T}} \tilde{T}_{\mathrm{h}} \right. \nonumber \\
& \qquad \left. \tilde{\Sigma} (\tilde{\Sigma}^{\mathrm{T}} \tilde{\Sigma})^{+} \tilde{V}_{\tilde{\Sigma}} (\tilde{S}_{\tilde{\Sigma}}^{\mathrm{T}} \tilde{S}_{\tilde{\Sigma}})^{1/4} \tilde{V}_{\tilde{\Sigma}}^{\mathrm{T}} \right) = \nonumber \\
& \operatorname{cond} \left( (\tilde{S}_{\tilde{\Sigma}}^{\mathrm{T}} \tilde{S}_{\tilde{\Sigma}})^{1/4} \tilde{V}_{\tilde{\Sigma}}^{\mathrm{T}} (\tilde{\Sigma}^{\mathrm{T}} \tilde{\Sigma})^{+} \tilde{\Sigma}^{\mathrm{T}} \tilde{T}_{\mathrm{h}} \right. \nonumber \\
& \qquad \left. \tilde{\Sigma} (\tilde{\Sigma}^{\mathrm{T}} \tilde{\Sigma})^{+} \tilde{V}_{\tilde{\Sigma}} (\tilde{S}_{\tilde{\Sigma}}^{\mathrm{T}} \tilde{S}_{\tilde{\Sigma}})^{1/4} \right). \tag{82}
\end{align}
$$

Such an approach would require the computation of the matrix
$\tilde{V}_{\tilde{\Sigma}}$ and $\tilde{S}_{\tilde{\Sigma}}$ which are prohibitively expensive to obtain. A key
observation, however, is that we do not need to use the entire
diagonal of $\tilde{S}_{\tilde{\Sigma}}$, but a logarithmic sampling of it will suffice.
In other words, define $D_{\tilde{\Sigma}}$ the vector containing the entries of
the diagonal of $\tilde{S}_{\tilde{\Sigma}}^{\mathrm{T}} \tilde{S}_{\tilde{\Sigma}}$ and define the block diagonal matrix

$$
\tilde{D}_{\tilde{\Sigma}, \alpha} = \operatorname{diag} \left( [D_{\tilde{\Sigma}}]_{N_S - N_{S, \alpha} + 1} \mathbf{I}_{N_{S, \alpha}^{\text{rem}}}, \right. \\
\left. [D_{\tilde{\Sigma}}]_{N_S - \frac{N_{S, \alpha}}{\alpha} + 1} \mathbf{I}_{\frac{N_{S, \alpha}}{\alpha}}, \dots, [D_{\tilde{\Sigma}}]_{N_S} \mathbf{I}_1 \right), \quad (83)
$$

where $N_{S,\alpha} = \alpha^{\lfloor \log_{\alpha}(N_S) \rfloor}$, $N_{S,\alpha}^{\text{rem}} = N_S - (1 - N_{S,\alpha})(1-\alpha)^{-1}$, and $\mathbf{I}_n$ is the identity matrix of size $n$, or, more programmatically,

$$
[\tilde{D}_{\tilde{\Sigma}, \alpha}]_{ii} = [D_{\tilde{\Sigma}}]_{f_{\tilde{\Sigma}}(i)}, \quad (84)
$$

with $f_{\tilde{\Sigma}}(i) = N_S - \alpha^{\lfloor \log_{\alpha}(N_S - i + 1) \rfloor} + 1$. Note that the construction of this matrix only requires explicit knowledge of $\log_{\alpha}(N_S)$ terms of $D_{\tilde{\Sigma}}$. Few passages—omitted here—suffice to show that

$$
\begin{equation}
\begin{split}
& \operatorname{cond} \left( \tilde{\mathbf{D}}_{\tilde{\Sigma}, \alpha}^{1/4} \tilde{\mathbf{V}}_{\tilde{\Sigma}}^{\mathrm{T}} (\tilde{\mathbf{\Sigma}}^{\mathrm{T}} \tilde{\mathbf{\Sigma}})^{+} \tilde{\mathbf{\Sigma}}^{\mathrm{T}} \tilde{\mathbf{T}}_h \tilde{\mathbf{\Sigma}} (\tilde{\mathbf{\Sigma}}^{\mathrm{T}} \tilde{\mathbf{\Sigma}})^{+} \tilde{\mathbf{V}}_{\tilde{\mathbf{\Sigma}}} \tilde{\mathbf{D}}_{\tilde{\mathbf{\Sigma}}, \alpha}^{1/4} \right) \\
& \qquad = O(\alpha) = O(1), \quad h \to 0,
\end{split}
\tag{85}
\end{equation}
$$

which is reminiscent of hierarchical strategies (see [7] and
references therein). Because $\tilde{V}_{\tilde{\Sigma}}$ is unitary, we obtain equiv-
alently

$$
\begin{equation}
\begin{split}
\operatorname{cond} \left( \tilde{V}_{\tilde{\Sigma}} \tilde{D}_{\tilde{\Sigma}, \alpha}^{1/4} \tilde{V}_{\tilde{\Sigma}}^{\mathrm{T}} (\tilde{\Sigma}^{\mathrm{T}} \tilde{\Sigma})^{+} \tilde{\Sigma}^{\mathrm{T}} \tilde{T}_{\mathrm{h}} \right. \\
\left. \qquad \tilde{\Sigma} (\tilde{\Sigma}^{\mathrm{T}} \tilde{\Sigma})^{+} \tilde{V}_{\tilde{\Sigma}} \tilde{D}_{\tilde{\Sigma}, \alpha}^{1/4} \tilde{V}_{\tilde{\Sigma}}^{\mathrm{T}} \right) &= O(\alpha) = O(1).
\end{split}
\tag{86}
\end{equation}
$$

This preconditioning strategy can be slightly altered to
leverage the filtered basis presented in Section III by introduc-
ing an additional Laplacian in (81) and adjusting the exponent
of $\tilde{S}_{\tilde{\Sigma}}^{\mathrm{T}} \tilde{S}_{\tilde{\Sigma}}$ accordingly. In particular, we have

$$
\begin{equation}
\begin{split}
& \tilde{\Sigma} (\tilde{\Sigma}^\mathrm{T} \tilde{\Sigma})^{+} \tilde{V}_{\tilde{\Sigma}} (\tilde{S}_{\tilde{\Sigma}}^\mathrm{T} \tilde{S}_{\tilde{\Sigma}})^{1/4} \tilde{V}_{\tilde{\Sigma}}^\mathrm{T} = \\
& \qquad \tilde{\Sigma} (\tilde{\Sigma}^\mathrm{T} \tilde{\Sigma})^{+} (\tilde{\Sigma}^\mathrm{T} \tilde{\Sigma}) \tilde{V}_{\tilde{\Sigma}} (\tilde{S}_{\tilde{\Sigma}}^\mathrm{T} \tilde{S}_{\tilde{\Sigma}})^{-3/4} \tilde{V}_{\tilde{\Sigma}}^\mathrm{T},
\end{split}
\tag{87}
\end{equation}
$$

which, following the reasoning detailed above, means that

$$
\tilde{B}_{\tilde{\Sigma}} := \tilde{\Sigma} (\tilde{\Sigma}^{\mathrm{T}} \tilde{\Sigma})^{+} (\tilde{\Sigma}^{\mathrm{T}} \tilde{\Sigma}) \tilde{V}_{\tilde{\Sigma}} \tilde{D}_{\tilde{\Sigma}, \alpha}^{-3/4} \tilde{V}_{\tilde{\Sigma}}^{\mathrm{T}} \quad (88)
$$

is a valid left and right symmetric preconditioner for $\tilde{T}_h$.
Finally, thanks to the properties introduced in Section III, we
have

$$
\begin{equation}
\begin{split}
& \tilde{\Sigma} (\tilde{\Sigma}^\mathrm{T} \tilde{\Sigma})^{+} (\tilde{\Sigma}^\mathrm{T} \tilde{\Sigma}) \tilde{V}_{\tilde{\Sigma}} \tilde{D}_{\tilde{\Sigma}, \alpha}^{-3/4} \tilde{V}_{\tilde{\Sigma}}^{\mathrm{T}} = \\
& \qquad \sum_{l=2}^{N_{S,\alpha}} (\tilde{\Sigma}_{\alpha^{l-1}-1} - \tilde{\Sigma}_{\alpha^{l-1}-1}) [D_{\tilde{\Sigma}}]_{N_S-\alpha^{l-1}+1}^{-3/4} \\
& \qquad + (\tilde{\Sigma} - \tilde{\Sigma}_{\alpha^{N_{S,\alpha}-1}}) [D_{\tilde{\Sigma}}]_{N_S-N_{S,\alpha}+1}^{-3/4} =: \tilde{\Sigma}_{p,\alpha} \quad (89)
\end{split}
\end{equation}
$$

and thus from (87) and (89) it follows that

$$
\operatorname{cond} \left( \tilde{\Sigma}_{p, \alpha}^{\mathrm{T}} \tilde{T}_{h} \tilde{\Sigma}_{p, \alpha} \right) = O(1), h \to 0. \quad (90)
$$

A similar reasoning for $\tilde{T}_s$, following from the precondi-
tioning of the hypersingular operator, leads to

$$
\operatorname{cond} \left( \left( (\tilde{S}_{\tilde{\Lambda}}^{\mathrm{T}} \tilde{S}_{\tilde{\Lambda}})^{-1/4} \tilde{V}_{\tilde{\Lambda}}^{\mathrm{T}} \tilde{\Lambda}^{\mathrm{T}} \tilde{T}_{\mathrm{s}} \tilde{\Lambda} \tilde{V}_{\tilde{\Lambda}} (\tilde{S}_{\tilde{\Lambda}}^{\mathrm{T}} \tilde{S}_{\tilde{\Lambda}})^{-1/4} \right) \right) = O(1) \tag{91}
$$

and

$$
\begin{equation}
\begin{split}
& \operatorname{cond} \left( \tilde{V}_{\tilde{\Lambda}} (\tilde{S}_{\tilde{\Lambda}}^{\mathrm{T}} \tilde{S}_{\tilde{\Lambda}})^{-1/4} \tilde{V}_{\tilde{\Lambda}}^{\mathrm{T}} \tilde{\Lambda}^{\mathrm{T}} \tilde{T}_{\mathrm{s}} \right. \\
& \qquad \left. \tilde{\Lambda} \tilde{V}_{\tilde{\Lambda}} (\tilde{S}_{\tilde{\Lambda}}^{\mathrm{T}} \tilde{S}_{\tilde{\Lambda}})^{-1/4} \tilde{V}_{\tilde{\Lambda}}^{\mathrm{T}} \right) = O(1), \quad h \to 0.
\end{split}
\tag{92}
\end{equation}
$$

From this, following a dual reasoning as the one of the
previous section, we obtain

$$
\operatorname{cond} \left( \tilde{\Lambda}_{p, \alpha}^{\mathrm{T}} \tilde{T}_{\mathrm{s}} \tilde{\Lambda}_{p, \alpha} \right) = O(1), h \to 0. \quad (93)
$$
where

$$
\begin{equation}
\begin{split}
\tilde{\Lambda} (\tilde{\Lambda}^T \tilde{\Lambda})^+ (\tilde{\Lambda}^T \tilde{\Lambda}) \tilde{V}_{\tilde{\Lambda}} \tilde{D}_{\tilde{\Lambda},\alpha}^{-1/4} \tilde{V}_{\tilde{\Lambda}}^T = & \\
& \sum_{l=2}^{N_{L,\alpha}} (\tilde{\Lambda}_{\alpha^{l-1}} - \tilde{\Lambda}_{\alpha^{l-1}-1}) [D_{\tilde{\Lambda}}]_{N_L-\alpha^{l-1}+1}^{-1/4} \\
& + (\tilde{\Lambda} - \tilde{\Lambda}_{\alpha^{N_L,\alpha}-1}) [D_{\tilde{\Lambda}}]_{N_L-N_{L,\alpha}+1}^{-1/4} := \tilde{\Lambda}_{p,\alpha},
\end{split}
\tag{94}
\end{equation}
$$

and

$$
\tilde{D}_{\tilde{\Lambda}, \alpha} = \operatorname{diag} \left( [D_{\tilde{\Lambda}}]_{N_L - N_{L, \alpha} + 1} \mathbf{I}_{N_{L, \alpha}^{\text{rem}}}, \right. \\
\left. [D_{\tilde{\Lambda}}]_{N_L - \frac{N_{L, \alpha}}{\alpha} + 1} \frac{\mathbf{I}_{\frac{N_{L, \alpha}}{\alpha}}}{\alpha}, \dots, [D_{\tilde{\Lambda}}]_{N_L} \mathbf{I}_1 \right), \quad (95)
$$

with $N_{L,\alpha}^{\text{rem}} = N_L - (1 - N_{L,\alpha})(1 - \alpha)^{-1}$, $D_{\tilde{\Lambda}}$ the vector containing the elements of the diagonal of $\tilde{S}_{\tilde{\Lambda}}^{\text{T}}\tilde{S}_{\tilde{\Lambda}}$, and $N_{L,\alpha} = \alpha^{\lfloor \log_{\alpha}(N_L) \rfloor}$.

The previous preconditioners can then be combined to
obtain a complete regularization of the EFIE system, for both
low-frequency and h-refinement breakdowns, that reads

$$
\tilde{\mathbf{W}}^{\mathrm{T}} \tilde{\mathbf{T}} \tilde{\mathbf{W}} \tilde{j} = \tilde{\mathbf{W}}^{\mathrm{T}} \tilde{\mathbf{v}}, \quad (96)
$$

where $\tilde{\mathbf{v}} = \mathbf{G}^{-1/2} \mathbf{v}$, $j = \mathbf{G}^{-1/2} \tilde{\mathbf{W}} \tilde{j}$, $\tilde{\mathbf{W}} = [\sqrt{c_{\tilde{\Lambda}}^T \tilde{\Lambda}_{p,\alpha}} \sqrt{c_{\tilde{\Sigma}}^T \tilde{\Sigma}_{p,\alpha}}]$,
$c_{\tilde{\Sigma}} = \|\tilde{\Sigma}_{p,\alpha}^T \tilde{T}_h \tilde{\Sigma}_{p,\alpha}\|^{-1}$,
$c_{\tilde{\Lambda}} = \|\tilde{\Lambda}_{p,\alpha}^T \tilde{T}_s \tilde{\Lambda}_{p,\alpha}\|^{-1}$, and where we assume that the appropriate
number of columns have been removed from $\tilde{\Sigma}_{p,\alpha}$ and $\tilde{\Lambda}_{p,\alpha}$
(e.g. 1 column must be removed from each for a simply con-
nected, closed scatterer) to account for the linear dependence
in the underlying Loop and Star bases [8], as is done in stan-
dard Loop-Star preconditioning. The reader should note that,
as in the case of standard Loop-Star functions, this operations
will create a small number of isolated singular values, that
however will not impact the convergence properties of the
preconditioned equation. This effect will not be present instead
in the scheme of next Section. The h-refinement regularization
effect of this preconditioner can be deduced from the previous
derivations for each of the potentials [50]. The low frequency
regularization, can be demonstrated following the same rea-
soning as for standard Loop-Star approaches [7], since the new
filtered bases retain the crucial properties that made Loop-Star
so adapted low-frequency regularization in the first place—
$\tilde{\Lambda}_{p,\alpha}^T \tilde{T}_h = \mathbf{0}$, $\tilde{T}_h \tilde{\Lambda}_{p,\alpha} = \mathbf{0}$, and $\tilde{\Lambda}_{p,\alpha}^T \tilde{\Sigma}_{p,\alpha} = \mathbf{0}$. Finally, we
have

$$
\mathrm{cond} \left( \tilde{\mathbf{W}}^{\mathrm{T}} \tilde{\mathbf{T}} \tilde{\mathbf{W}} \right) = O(1), \quad \text{when } h \to 0, k \to 0. \quad (97)
$$

B. Quasi-Helmholtz Filters Approach

In several application scenarios, an explicit quasi-Helmholtz
decomposition, such as the Loop-Star decomposition, is not
necessary, and quasi-Helmholtz projectors [7] could be used
instead. Similarly, instead of using filtered Loop-Star pre-
conditioning approaches, basis-free approaches, based on the
quasi-Helmholtz Laplacian filters, will often be more effective.
This Section will explore this approach that, as an additional
advantage, will also have the avoidance of the burden of
global-loop detection for multiply connected scatterers.

Following the same philosophy as in Section VI-A, we will
form preconditioners for the solenoidal part of $\tilde{\mathbf{T}}_s$ and for $\tilde{\mathbf{T}}_h$

that will then be combined into a full EFIE preconditioner. We
can transition from a basis-based Helmholtz decomposition
to a projector based Helmholtz decomposition by leveraging
the correspondences between $\tilde{\Sigma}$ and $\tilde{\Lambda}$ and their respective
projectors $\tilde{\mathbf{P}}^{\Sigma}$ and $\tilde{\mathbf{P}}^{\Lambda}$. In particular, because $\tilde{\mathbf{B}}_{\tilde{\Sigma}}$ was a valid
preconditioner for $\tilde{\mathbf{T}}_h$ (equations (89) and (90)), $[\tilde{\mathbf{B}}_{\tilde{\Sigma}} \mathbf{0}]$,
once applied left and right to $\tilde{\mathbf{T}}_h$ will yield a block diag-
onal matrix which is well conditioned away from its large
nullspace. This, in turns, means that $[\tilde{\mathbf{C}}_{\tilde{\Sigma}} \mathbf{0}]$, with $\tilde{\mathbf{C}}_{\tilde{\Sigma}} =$
$\tilde{\Sigma} (\tilde{\Sigma}^{\mathrm{T}} \tilde{\Sigma})^{+} (\tilde{\Sigma}^{\mathrm{T}} \tilde{\Sigma}) \tilde{V}_{\tilde{\Sigma}} \tilde{D}_{\tilde{\Sigma}}^{-5/4} \tilde{V}_{\tilde{\Sigma}}^{\mathrm{T}} \tilde{V}_{\tilde{\Sigma}} \tilde{D}_{\tilde{\Sigma}}^{1/2}$, will also yield
a well-conditioned (up to its nullspace) matrix. Finally, be-
cause multiplications by unitary matrices do not compromise
conditioning properties, we can form the preconditioner

$$
\tilde{\Sigma} (\tilde{\Sigma}^T \tilde{\Sigma})^+ (\tilde{\Sigma}^T \tilde{\Sigma}) \tilde{V}_{\tilde{\Sigma}} \tilde{D}_{\tilde{\Sigma}}^{-5/4} \tilde{V}_{\tilde{\Sigma}}^T \tilde{V}_{\tilde{\Sigma}} [\tilde{D}_{\tilde{\Sigma}}^{1/2} \mathbf{0}] \tilde{U}_{\tilde{\Sigma}}^T = \\
\tilde{\Sigma} (\tilde{\Sigma}^T \tilde{\Sigma})^+ \tilde{V}_{\tilde{\Sigma}} \tilde{D}_{\tilde{\Sigma}}^{-1/4} \tilde{V}_{\tilde{\Sigma}}^T \tilde{\Sigma}^T. \quad (98)
$$

This allows us to form the preconditioner $\tilde{Q}_{p,\alpha}^{\tilde{\Sigma}}$, of additive
Schwarz type, based on quasi-Helmholtz filters

$$
\begin{equation}
\begin{split}
\tilde{\Sigma} (\tilde{\Sigma}^T \tilde{\Sigma})^+ \tilde{V}_{\tilde{\Sigma}} \tilde{D}_{\tilde{\Sigma}}^{-1/4} \tilde{V}_{\tilde{\Sigma}}^T \tilde{\Sigma}^T = & \\
& \sum_{l=2}^{N_{S,\alpha}} (\tilde{\mathbf{P}}_{\tilde{\alpha}^{l-1}}^{\tilde{\Sigma}} - \tilde{\mathbf{P}}_{\tilde{\alpha}^{l-1}-1}^{\tilde{\Sigma}}) [D_{\tilde{\Sigma}}]_{N_S-\alpha^{l-1}+1}^{-1/4} \\
& + (\tilde{\mathbf{P}}^{\tilde{\Sigma}} - \tilde{\mathbf{P}}_{\tilde{\alpha}^{N_{S,\alpha}-1}}^{\tilde{\Sigma}}) [D_{\tilde{\Sigma}}]_{N_S-N_{S,\alpha}+1}^{-1/4} := \tilde{Q}_{p,\alpha}^{\tilde{\Sigma}}
\end{split}
\tag{99}
\end{equation}
$$

for which

$$
\mathrm{cond} \left( \tilde{Q}_{p,\alpha}^{\tilde{\Sigma}} \tilde{T}_h \tilde{Q}_{p,\alpha}^{\tilde{\Sigma}} \right) = O(1), h \to 0. \quad (100)
$$

Similarly, a preconditioner for the solenoidal part of $\tilde{\mathbf{T}}_s$ is

$$
\tilde{Q}_{p,\alpha}^{\tilde{\Lambda}} := \sum_{l=2}^{N_{L,\alpha}} \left( \tilde{\mathbf{P}}_{\tilde{\alpha}^{l-1}}^{\tilde{\Lambda}} - \tilde{\mathbf{P}}_{\tilde{\alpha}^{l-1}-1}^{\tilde{\Lambda}} \right) [\tilde{D}_{\tilde{\Lambda}}]_{N_L-\alpha^{l-1}+1}^{1/4} \\
+ \left( \tilde{\mathbf{P}}^{\tilde{\Lambda}} - \tilde{\mathbf{P}}_{\tilde{\alpha}^{N_{L,\alpha}-1}}^{\tilde{\Lambda}} \right) [\tilde{D}_{\tilde{\Lambda}}]_{N_L-N_{L,\alpha}+1}^{1/4} \quad (101)
$$

for which

$$
\mathrm{cond} \left( \tilde{Q}_{p,\alpha}^{\tilde{\Lambda}} \tilde{T}_s \tilde{Q}_{p,\alpha}^{\tilde{\Lambda}} \right) = O(1), h \to 0. \quad (102)
$$

The full EFIE preconditioner is then an appropriate linear
combination of the solenoidal and non-solenoidal precondi-
tioners above to cure also the low-frequency breakdown. In
particular we define

$$
\tilde{Q} = \sqrt{b_{\tilde{\Lambda}} \tilde{Q}_{p,\alpha}^{\tilde{\Lambda}}} + i \sqrt{b_{\tilde{\Sigma}} \tilde{Q}_{p,\alpha}^{\tilde{\Sigma}}} + \sqrt{b_{\tilde{H}} \tilde{P}^H}, \quad (103)
$$

where $\tilde{\mathbf{P}}^H = \mathbf{I} - \tilde{\mathbf{P}}^\Sigma - \tilde{\mathbf{P}}^\Lambda$ and

$$
b_{\tilde{\Lambda}} = \| \tilde{Q}_{p,\alpha}^{\tilde{\Lambda}} \tilde{T}_s \tilde{Q}_{p,\alpha}^{\tilde{\Lambda}} \|^{-1}, \quad (104)
$$

$$
b_{\tilde{\Sigma}} = \| \tilde{Q}_{p,\alpha}^{\tilde{\Sigma}} \tilde{T}_h \tilde{Q}_{p,\alpha}^{\tilde{\Sigma}} \|^{-1}, \quad (105)
$$

$$
b_{\tilde{H}} = \| \tilde{\mathbf{P}}^{\tilde{H}} \tilde{T}_s \tilde{\mathbf{P}}^{\tilde{H}} \|^{-1}, \quad (106)
$$

account for the frequency-scaling of the operators and the
diameter of Γ. The preconditioned EFIE system is

$$
\tilde{Q} \tilde{T} \tilde{Q} \tilde{j}_{qH} = \tilde{Q} \tilde{v}, \quad (107)
$$

with $j = \mathbf{G}^{-1/2} \tilde{Q} \tilde{j}_{qH}$.
## VII. IMPLEMENTATION RELATED DETAILS AND FURTHER IMPROVEMENTS

In addition to the efficient filtering algorithms presented in Section V, obtaining a fast and efficient implementation of the proposed preconditioning scheme based on filtered projectors requires that particular attention be given to parts of their implementation. First, all the terms of the form $T_h Q_{p,\alpha}^{\Lambda}$, $Q_{p,\alpha}^{\Lambda} T_h$, $P^H T_h$, or $T_h P^H$ must be explicitly set to 0 to avoid numerical instabilities. Further treatments on the right hand side and on the solution vector, are required to ensure that the solution of the system remains accurate until arbitrarily low frequencies. These treatments are straightforward generalization of those required for standard quasi-Helmholtz preconditioning techniques that can be found in [7].

The condition numbers obtained when employing the schemes introduced in Section VI-A and Section VI-B, while stable, can be further brought down by slightly modifying the preconditioners. The diagonal preconditioning based on the theoretical Laplacian eigenvalues can be altered to instead employ matrix norms; the new preconditioners then become

$$Q_{p,\alpha}^{\Sigma} = \sum_{l=2}^{N_{S,\alpha}} (P_{\alpha^{l-1}}^{\Sigma} - P_{\alpha^{l-1}-1}^{\Sigma}) b_l + (P^{\Sigma} - P_{\alpha^{N_{S,\alpha}-1}}^{\Sigma}) b_{N_{S,\alpha}+1}, \quad (108)$$

where

$$b_l = \|(P_{\alpha^{l-1}}^{\Sigma} - P_{\alpha^{l-1}-1}^{\Sigma})^{\mathrm{T}} T_h (P_{\alpha^{l-1}}^{\Sigma} - P_{\alpha^{l-1}-1}^{\Sigma})\|^{-1/2}, \\ 2 \le l \le N_{S,\alpha}, \quad (109)$$

$$b_{N_{S,\alpha}+1} = \|(P^{\Sigma} - P_{\alpha^{N_{S,\alpha}-1}}^{\Sigma})^{\mathrm{T}} T_h (P^{\Sigma} - P_{\alpha^{N_{S,\alpha}-1}}^{\Sigma})\|^{-1/2} \quad (110)$$

The same modification can be performed for $Q_{p,\alpha}^{\Lambda}$ that becomes

$$Q_{p,\alpha}^{\Lambda} = \sum_{l=2}^{N_{L,\alpha}} (P_{\alpha^{l-1}}^{\Lambda} - P_{\alpha^{l-1}-1}^{\Lambda}) d_l + (P^{\Lambda} - P_{\alpha^{N_{L,\alpha}-1}}^{\Lambda}) d_{N_{S,\alpha}+1}, \quad (111)$$

with

$$d_l = \|(P_{\alpha^{l-1}}^{\Lambda} - P_{\alpha^{l-1}-1}^{\Lambda})^{\mathrm{T}} T_s (P_{\alpha^{l-1}}^{\Lambda} - P_{\alpha^{l-1}-1}^{\Lambda})\|^{-1/2}, \\ 2 \le l \le N_{L,\alpha}, \quad (112)$$

$$d_{N_{S,\alpha}+1} = \|(P^{\Lambda} - P_{\alpha^{N_{L,\alpha}-1}}^{\Lambda})^{\mathrm{T}} T_s (P^{\Lambda} - P_{\alpha^{N_{L,\alpha}-1}}^{\Lambda})\|^{-1/2} \quad (113)$$

To ensure that the overall complexity of the algorithm is not increased, the values of {$b_l$}$_l$ and {$d_l$}$_l$ can be efficiently computed using, for example, power methods. The reader should note that the preconditioning approach delineated above requires filter profiles with support both proportional to and independent from the number of unknowns, which can be

<figure><img src="image_4.png" /><figcaption>Fig. 2: Spectrum of the solenoidal part of the vector potential, its preconditioner, and its preconditioned counterpart. These spectra have been obtained for a smoothly-deformed sphere with a maximum diameter of 7.17 m (see insert), a frequency of $10^6$ Hz, and for two different average edge lengths 0.31 m and 0.20 m. The spectra have been normalized so that their first singular value is one, for readability. Perfect filters built out of SVD have been used in these results.</figcaption></figure>

efficiently obtained with the approaches described in Section V. As said in the previous Section, filters in the transition region could be less efficient to obtain, as the Chebyshev approach decreases in efficiency away from the middle of the spectrum [34]. All preconditioning real case scenarios presented here, however, are not impacted by this fact as shown in Section VIII.

## VIII. NUMERICAL RESULTS

All numerical results presented in this section have been obtained with non-normalized matrices ($\Lambda$, $\Sigma$) to illustrate that graph matrices are often enough for practical cases. Equally good or superior performance, however, can be obtained by using normalized matrices ($\hat{\Lambda}$, $\hat{\Sigma}$) instead. In the first set of examples we have leveraged perfect filters obtained by SVD before presenting results based on SVD-free approaches. The filtered Loop-Star preconditioning approach presented in Section VI-A leverages the spectral equivalences between the appropriately scaled filtered bases and $T_s$ and $T_h$. To numerically illustrate these equivalences, the spectra of these operators and their preconditioned counterparts are illustrated in Figures 2 and 3. These spectra correspond to a smoothly deformed sphere (see Fig. 2 and 3), and the ordering of the singular values is obtained by projection against the graph Laplacians' eigenvectors. The original spectrum of $T_s$ and $T_h$ show the expected $\xi^{-1/2}$ and $\xi^{1/2}$—with $\xi$ the spectral index—behaviors, predicted by pseudo-differential operator theory. Given the construction of the preconditioners, it is then not surprising that the preconditioned operators show a spectrum bounded (and away from zero) with the expected variations in the spectrum.

To illustrate that the preconditioning schemes based on filtered bases do regularize the EFIE, the condition number of
<figure><img src="image_5.png" /><figcaption>Fig. 3: Spectrum of the non-solenoidal part of the scalar potential, its preconditioner, and its preconditioned counterpart. These spectra have been obtained for a smoothly-deformed sphere with a maximum diameter of 7.17 m (see insert), a frequency of $10^6$ Hz, and for two different average edge lengths 0.31 m and 0.20 m. The spectra have been normalized so that their first singular value is one, for readability. Perfect filters built out of SVD have been used in these results.</figcaption></figure>

the original and preconditioned schemes will be compared for varying frequencies and discretizations. First, the conditioning of a filtered Loop-Star preconditioned EFIE for the NASA almond [51] is reported in Figure 4. The low frequency and dense discretization breakdowns of the original equations are apparent, while the preconditioned equation (corresponding to (96)) shows a constant conditioning. This is in contrast with the standard Loop-Star approach that does regularize the low frequency conditioning breakdown, but actually worsens the dense discretization behavior of the equation.

A similar study has been performed with the filtered projectors schemes. In Figures 5 and 6 the spectra of the dominant solenoidal and non-solenoidal parts of the EFIE operators are displayed alongside their preconditioners. The preconditioning performance on the overall EFIE system is illustrated in Figure 7 for a torus. The approach yields satisfactory conditioning that remains stable in both low frequency and dense discretization, which in turns shows that the scheme can also handle multiply-connected geometries.

Finally, a conditioning study of the NASA almond is reported in Figure 8 that has been obtained using Chebyshev-interpolated filters (78) corresponding to Butterworth filters of order 100, expanded into 200 Chebyshev polynomials. The coefficients of the filters are obtained via the norm estimates detailed in (108) and (111) and the cutting point of the filters is determined using the approximate Laplacian spectrum described below (76). The excellent stability of preconditioned scheme for a structure such as the NASA almond showcases the effectiveness of the scheme when using the fast techniques presented in this paper.

<figure><img src="image_6.png" /><figcaption>Fig. 4: Condition number of the EFIE (5), Loop-Star EFIE, and filtered Loop-Star EFIE (96) as a function of discretization for several frequencies. The condition number has been obtained after eliminating the isolated singular values, which have minimal impact on the convergence, arising from the deletion of one column from each of the preconditioning matrices. The solid lines correspond to a simulating frequency of $10^6$ Hz and the dotted lines to a frequency of $10^4$ Hz. The simulated structure is the NASA almond re-scaled to be enclosed in a bounding box of diameter 1.09 m. Perfect filters built out of SVD have been used in these results.</figcaption></figure>

<figure><img src="image_7.png" /><figcaption>Fig. 5: Spectrum of the solenoidal part of the vector potential, its preconditioner, and its preconditioned counterpart. These spectra have been obtained for a smoothly-deformed sphere with a maximum diameter of 7.17 m (see insert), a frequency of $10^6$ Hz, and for two different average edge lengths 0.31 m and 0.20 m. The spectra have been normalized so that their first singular value is one, for readability. Perfect filters built out of SVD have been used in these results.</figcaption></figure>
<figure><img src="image_8.png" /><figcaption>Fig. 6: Spectrum of the non-solenoidal part of the scalar potential, its preconditioner, and its preconditioned counterpart. These spectra have been obtained for a smoothly-deformed sphere with a maximum diameter of 7.17 m (see insert), a frequency of $10^6$ Hz, and for two different average edge lengths 0.31 m and 0.20 m. The spectra have been normalized so that their first singular value is one, for readability. Perfect filters built out of SVD have been used in these results.</figcaption></figure>

<figure><img src="image_9.png" /><figcaption>Fig. 7: Condition number of the EFIE (5), quasi-Helmholtz (qH) projector EFIE, and filtered qH projector EFIE (107) as a function of discretization for several frequencies. The solid lines correspond to a simulating frequency of $10^6$ Hz and the dotted lines to a frequency of $10^4$ Hz. The simulated structure is a torus with inner radius 0.9 m and outer radius 1.1 m. Perfect filters built out of SVD have been used in these results.</figcaption></figure>

<figure><img src="image_10.png" /><figcaption>Fig. 8: Condition number of the EFIE (5), quasi-Helmholtz (qH) projector EFIE, and filtered qH projector EFIE (107) as a function of discretization for several frequencies. The simulated structure is the NASA almond re-scaled to be enclosed in a bounding box of diameter 1.09 m. The preconditioner is built without using SVDs, but by leveraging Chebyshev-interpolated filters (78) corresponding to Butterworth filters of order 100, expanded into 200 Chebyshev polynomials.</figcaption></figure>

## IX. CONCLUSION

A new family of strategies has been introduced for performing filtered quasi-Helmholtz decompositions of electromagnetic integral equations: the filtered Loop-Star decompositions and the quasi-Helmholtz Laplacian filters. These new tools are capable of manipulating large parts of the operators' spectra to obtain new families of preconditioners and fast solvers. A first application to the case of frequency and h-refinement preconditioning of the electric field integral equation has been presented and numerical results have shown the practical effectiveness of the newly proposed tools.

## APPENDIX A
COMPLEMENTARITY OF THE PROJECTORS

In this appendix, we show that the properties $\tilde{\mathbf{P}}^A + \tilde{\mathbf{P}}^\Sigma = \mathbf{I}$ and $\tilde{\mathbf{P}}^A + \tilde{\mathbf{P}}^\Sigma = \mathbf{I}$ hold true on simply connected geometries. To this end, we first prove that the normalized coefficients $\tilde{j}$ of the RWG functions can be decomposed with $\tilde{\mathbf{\Lambda}}$ and $\tilde{\mathbf{\Sigma}}$, similarly as in (9) where we assume that the proper number of columns from the matrices have been removed as is standard to ensure a full column rank, such that

$$ \tilde{j} = \mathbf{G}^{-\frac{1}{2}} \tilde{j} = \tilde{\mathbf{\Lambda}} \tilde{i} + \tilde{\mathbf{\Sigma}} \tilde{s} \quad (114) $$

in which $\tilde{i}$ and $\tilde{s}$ are the coefficient vectors of the normalized Loop and Star parts in this decomposition. Since $\mathbf{G}$, $\mathbf{G}_p$, and $\mathbf{G}_\lambda$ are invertible matrices, we have $\mathrm{rank}(\tilde{\mathbf{\Sigma}}) = \mathrm{rank}(\mathbf{\Sigma})$ and $\mathrm{rank}(\tilde{\mathbf{\Lambda}}) = \mathrm{rank}(\mathbf{\Lambda})$. Moreover, since $\tilde{\mathbf{\Lambda}}^T \tilde{\mathbf{\Sigma}} = \mathbf{0}$, we also obtain that $\tilde{\mathbf{\Lambda}}$ and $\tilde{\mathbf{\Sigma}}$ have their column linearly independent, which yields $\mathrm{rank}([\tilde{\mathbf{\Lambda}} \tilde{\mathbf{\Sigma}}]) = \mathrm{rank}(\tilde{\mathbf{\Lambda}}) + \mathrm{rank}(\tilde{\mathbf{\Sigma}}) = N$, from which the existence and (unicity) of (114) follows.

Subsequently, using (114), we can form a new set of normalized projectors to retrieve $\tilde{\mathbf{\Lambda}}\tilde{i}$ and $\tilde{\mathbf{\Sigma}}\tilde{s}$ separately. The
first step is to apply $\tilde{\Lambda}^T$ and $\tilde{\Sigma}^T$ to (114) to express $\tilde{j}$ in the
two different bases

$$
\tilde{\Lambda}^T \tilde{j} = \tilde{\Lambda}^T \tilde{\Lambda} \tilde{l}, \qquad (115)
$$

$$
\tilde{\Sigma}^T \tilde{j} = \tilde{\Sigma}^T \tilde{\Sigma} \tilde{s}, \quad (116)
$$

since $\tilde{\Lambda}^T \tilde{\Sigma} = \mathbf{0}$ and $\tilde{\Sigma}^T \tilde{\Lambda} = \mathbf{0}$, given (20). Subsequently, we
express the coefficients of the normalized Loop and Star bases
as a function of $\tilde{j}$

$$
\tilde{l} = (\tilde{\Lambda}^T \tilde{\Lambda})^+ \tilde{\Lambda}^T \tilde{j}, \quad (117)
$$

$$
\tilde{s} = (\tilde{\Sigma}^T \tilde{\Sigma})^+ \tilde{\Sigma}^T \tilde{j}. \quad (118)
$$

Finally, we express $\tilde{\Lambda}\tilde{l}$ and $\tilde{\Sigma}\tilde{s}$ in terms of $\tilde{j}$ by applying $\tilde{\Lambda}$ and $\tilde{\Sigma}$ to (117) and (118)

$$
\tilde{\Lambda}\tilde{l} = \tilde{\Lambda} (\tilde{\Lambda}^T \tilde{\Lambda})^+ \tilde{\Lambda}^T \tilde{j} = \tilde{P}^A \tilde{j}, \quad (119)
$$

$$
\tilde{\Sigma}\tilde{s} = \tilde{\Sigma} (\tilde{\Sigma}^T \tilde{\Sigma})^+ \tilde{\Sigma}^T \tilde{j} = \tilde{P}^{\Sigma} \tilde{j}, \quad (120)
$$

and we obtain that $\tilde{P}^A + \tilde{P}^\Sigma = \mathbf{I}$ by leveraging (119), (120), and (114). Following the same procedure, except that now $\tilde{\Lambda}$ and $\tilde{\Sigma}$ are employed in the initial decomposition, we can show that the property $\tilde{P}^A + \tilde{P}^\Sigma = \mathbf{I}$ also holds true.

ACKNOWLEDGMENT

This work has been funded in part by the European Research
Council (ERC) under the European Union’s Horizon 2020
research and innovation program (ERC project 321, grant
No.724846) and in part by the ANR Labex CominLabs under
the project “CYCLE”.

REFERENCES

[1] J.-C. Nedelec, *Acoustic and Electromagnetic Equations: Integral Representations for Harmonic Problems*. Springer Science & Business Media, Mar. 2001.

[2] W. C. Gibson, *The Method of Moments in Electromagnetics* 2nd Ed. CRC press, 2014.

[3] J.-M. Jin, *Theory and Computation of Electromagnetic Fields*, 2nd ed. Piscataway, NJ: IEEE Press, 2015.

[4] W. C. Chew, J.-M. Jin, E. Michielssen, and J. M. Song, Eds., *Fast and Efficient Algorithms in Computational Electromagnetics*. Artech House, 2001.

[5] O. Axelsson, *Iterative Solution Methods*. Cambridge university press, 1996.

[6] D. L. Colton and R. Kress, *Integral Equation Methods in Scattering Theory*. Philadelphia: Society for Industrial and Applied Mathematics, 2013.

[7] S. B. Adrian, A. Dély, D. Consoli, A. Merlini, and F. P. Andriulli, "Electromagnetic Integral Equations: Insights in Conditioning and Preconditioning," *IEEE Open Journal of Antennas and Propagation*, vol. 2, pp. 1143-1174, 2021.

[8] D. R. Wilton, “Topological consideration in surface patch and volume cell modeling of electromagnetic scatterers,” in *Proc. URSI Int. Symp. Electromagn. Theory*, 1983, pp. 65–68.

[9] G. Vecchi, “Loop-star decomposition of basis functions in the discretization of the EFIE,” *IEEE Transactions on Antennas and Propagation*, vol. 47, no. 2, pp. 339–346, 1999.

[10] J.-F. Lee, R. Lee, and R. Burkholder, “Loop star basis functions and a robust preconditioner for EFIE scattering problems,” *IEEE Transactions on Antennas and Propagation*, vol. 51, no. 8, pp. 1855–1863, Aug. 2003.

[11] T. F. Eibert, “Iterative-solver convergence for loop-star and loop-tree decompositions in method-of-moments solutions of the electric-field integral equation,” *Antennas and Propagation Magazine*, IEEE, vol. 46, no. 3, pp. 80–85, 2004.

[12] F. P. Andriulli, “Loop-Star and Loop-Tree Decompositions: Analysis and Efficient Algorithms,” *IEEE Transactions on Antennas and Propagation*, vol. 60, no. 5, pp. 2347–2356, May 2012.

[13] F. Vipiana, P. Pirinoli, and G. Vecchi, “A multiresolution method of moments for triangular meshes,” *IEEE transactions on antennas and propagation*, vol. 53, no. 7, pp. 2247–2258, 2005.

[14] F. P. Andriulli, A. Tabacco, and G. Vecchi, “A multiresolution approach to the electric field integral equation in antenna problems,” *SIAM Journal on Scientific Computing*, vol. 29, no. 1, pp. 1–21, 2007.

[15] F. P. Andriulli, F. Vipiana, and G. Vecchi, “Hierarchical Bases for Nonhierarchic 3-D Triangular Meshes,” *IEEE Transactions on Antennas and Propagation*, vol. 56, no. 8, pp. 2288–2297, Aug. 2008.

[16] R.-S. Chen, J. Ding, D. Ding, Z. Fan, and D. Wang, "A multiresolution curvilinear rao-wilton-glisson basis function for fast analysis of electromagnetic scattering," *IEEE Transactions on Antennas and Propagation*, vol. 57, no. 10, pp. 3179-3188, 2009.

[17] J. E. O. Guzman, S. Adrian, R. Mitharwal, Y. Beghein, T. Eibert, K. Cools, and F. Andriulli, “On the Hierarchical Preconditioning of the PMCHWT Integral Equation on Simply and Multiply Connected Geometries,” *IEEE Antennas and Wireless Propagation Letters*, vol. PP, no. 99, pp. 1–1, 2016.

[18] S. Adrian, F. Andriulli, and T. Eibert, “A hierarchical preconditioner for the electric field integral equation on unstructured meshes based on primal and dual Haar bases,” *Journal of Computational Physics*, vol. 330, pp. 365–379, Feb. 2017.

[19] R. Adams and G. Brown, “Stabilisation procedure for electric field integral equation,” *Electronics Letters*, vol. 35, no. 23, pp. 2015–2016, Nov. 1999.

[20] H. Contopanagos, B. Dembart, M. Epton, J. Ottusch, V. Rokhlin, J. Visher, and S. Wandzura, “Well-conditioned boundary integral equations for three-dimensional electromagnetic scattering,” *IEEE Transactions on Antennas and Propagation*, vol. 50, no. 12, pp. 1824–1830, Dec. 2002.

[21] R. Adams, “Physical and Analytical Properties of a Stabilized Electric Field Integral Equation,” *IEEE Transactions on Antennas and Propagation*, vol. 52, no. 2, pp. 362–372, Feb. 2004.

[22] R. Adams and N. Champagne, “A Numerical Implementation of a Modified Form of the Electric Field Integral Equation,” *IEEE Transactions on Antennas and Propagation*, vol. 52, no. 9, pp. 2262–2266, Sep. 2004.

[23] S. Borel, D. Levadoux, and F. Alouges, “A new well-conditioned Integral formulation for Maxwell equations in three dimensions,” *IEEE Transactions on Antennas and Propagation*, vol. 53, no. 9, pp. 2995–3004, Sep. 2005.

[24] F. P. Andriulli, K. Cools, H. Bagci, F. Olyslager, A. Buffa, S. Christiansen, and E. Michielssen, “A Multiplicative Calderon Preconditioner for the Electric Field Integral Equation,” *IEEE Transactions on Antennas and Propagation*, vol. 56, no. 8, pp. 2398–2412, Aug. 2008.

[25] M. B. Stephanson and J.-F. Lee, “Preconditioned Electric Field Integral Equation Using Calderon Identities and Dual Loop/Star Basis Functions,” *IEEE Transactions on Antennas and Propagation*, vol. 57, no. 4, pp. 1274–1279, Apr. 2009.

[26] S. Borel, D. P. Levadoux, and F. Alouges, “A new well-conditioned integral formulation for maxwell equations in three dimensions,” *IEEE Transactions on antennas and propagation*, vol. 53, no. 9, pp. 2995–3004, 2005.

[27] F. P. Andriulli, K. Cools, I. Bogaert, and E. Michielssen, “On a Well-Conditioned Electric Field Integral Operator for Multiply Connected Geometries,” *IEEE Transactions on Antennas and Propagation*, vol. 61, no. 4, pp. 2077–2087, Apr. 2013.

[28] D. Dobbelaere, D. De Zutter, J. Van Hese, J. Sercu, T. Boonen, and H. Rogier, “A Calderón multiplicative preconditioner for the electromagnetic Poincaré–Steklov operator of a heterogeneous domain with scattering applications,” *Journal of Computational Physics*, vol. 303, pp. 355–371, Dec. 2015.

[29] A. Merlini, Y. Beghein, K. Cools, E. Michielssen, and F. P. Andriulli, “Magnetic and Combined Field Integral Equations Based on the Quasi-Helmholtz Projectors,” *IEEE Transactions on Antennas and Propagation*, vol. 68, no. 5, pp. 3834–3846, May 2020.

[30] S. Adrian, F. Andriulli, and T. Eibert, “On a Refinement-Free Calderón Multiplicative Preconditioner for the Electric Field Integral Equation,” *Journal of Computational Physics*, Oct. 2018.

[31] O. G. Johnson, C. A. Micchelli, and G. Paul, “Polynomial preconditioners for conjugate gradient calculations,” *SIAM Journal on Numerical Analysis*, vol. 20, no. 2, pp. 362–376, 1983.

[32] S. F. Ashby, T. A. Manteuffel, and J. S. Otto, “A comparison of adaptive chebyshev and least squares polynomial preconditioning for
hermitian positive definite linear systems,” *SIAM Journal on Scientific and Statistical Computing*, vol. 13, no. 1, pp. 1–29, 1992.

[33] D. K. Hammond, P. Vandergheynst, and R. Gribonval, “Wavelets on graphs via spectral graph theory,” *Applied and Computational Harmonic Analysis*, vol. 30, no. 2, pp. 129–150, Mar. 2011.

[34] R. Levie, F. Monti, X. Bresson, and M. M. Bronstein, “Cayleynets: Graph convolutional neural networks with complex rational spectral filters,” *IEEE Transactions on Signal Processing*, vol. 67, no. 1, pp. 97–109, 2018.

[35] L. Rahmouni and F. P. Andriulli, “A New Preconditioner for the EFIE Based on Primal and Dual Graph Laplacian Spectral Filters,” in *2019 International Conference on Electromagnetics in Advanced Applications (ICEAA)*. Granada, Spain: IEEE, Sep. 2019, pp. 1342–1344.

[36] A. Merlini, C. Henry, D. Consoli, L. Rahmouni, and F. P. Andriulli, “Laplacian filters for integral equations: Further developments and fast algorithms,” in *2022 IEEE International Symposium on Antennas and Propagation and USNC-URSI Radio Science Meeting (AP-S/URSI)*, 2022, pp. 1932–1933.

[37] S. Rao, D. Wilton, and A. Glisson, “Electromagnetic scattering by surfaces of arbitrary shape,” *IEEE Transactions on Antennas and Propagation*, vol. 30, no. 3, pp. 409–418, May 1982.

[38] A. Buffa and S. Christiansen, “A dual finite element complex on the barycentric refinement,” *Mathematics of Computation*, vol. 76, no. 260, pp. 1743–1769, 2007.

[39] Q. Chen and D. Wilton, “Electromagnetic scattering by three-dimensional arbitrary complex material/conducting bodies,” in *Antennas and Propagation Society International Symposium, 1990. AP-S. Merging Technologies for the 90's. Digest.*, May 1990, pp. 590–593 vol.2.

[40] J. Mautz and R. Harrington, “An E-field solution for a conducting surface small or comparable to the wavelength,” *IEEE Transactions on Antennas and Propagation*, vol. 32, no. 4, pp. 330–339, 1984.

[41] J. S. Lim, S. Rao, and D. R. Wilton, “A Novel Technique to Calculate the Electromagnetic Scattering by Surfaces of Arbitrary Shape,” in *1993 URSI Radio Science Meeting Digest*, 1993, p. 322.

[42] W.-L. Wu, A. W. Glisson, and D. Kajfez, “A study of two numerical solution procedures for the electric field integral equation at low frequency,” *Applied Computational Electromagnetics Society Journal*, vol. 10, no. 3, pp. 69–80, 1995.

[43] G. H. Golub and C. F. Van Loan, *Matrix Computations*. JHU Press, 2012, vol. 3.

[44] N. Higham, *Functions of Matrices*, ser. Other Titles in Applied Mathematics. Society for Industrial and Applied Mathematics, Jan. 2008.

[45] W. H. Press, *Numerical Recipes 3rd Edition: The Art of Scientific Computing*. Cambridge University Press, Sep. 2007.

[46] R. Mitharwal and F. Andriulli, “On the Multiplicative Regularization of Graph Laplacians on Closed and Open Structures With Applications to Spectral Partitioning,” *IEEE Access*, vol. 2, pp. 788–796, 2014.

[47] M. O'Neil, “Second-kind integral equations for the Laplace-Beltrami problem on surfaces in three dimensions,” *Advances in Computational Mathematics*, vol. 44, no. 5, pp. 1385–1409, Oct. 2018.

[48] J.-S. Zhao and W. C. Chew, “Integral equation solution of Maxwell’s equations from zero frequency to microwave frequencies,” *Antennas and Propagation, IEEE Transactions on*, vol. 48, no. 10, pp. 1635–1645, 2000.

[49] D. N. Arnold, P. B. Bochev, R. B. Lehoucq, R. A. Nicolaides, and M. Shashkov, *Compatible spatial discretizations*. Springer Science & Business Media, 2007, vol. 142.

[50] Y. Boubendir and C. Turc, “Well-conditioned boundary integral equation formulations for the solution of high-frequency electromagnetic scattering problems,” *Computers & Mathematics with Applications*, vol. 67, no. 10, pp. 1772–1805, Jun. 2014.

[51] A. Woo, H. Wang, M. Schuh, and M. Sanders, “EM programmer’s notebook-benchmark radar targets for the validation of computational electromagnetics programs,” *IEEE Antennas and Propagation Magazine*, vol. 35, no. 1, pp. 84–89, 1993.

Michael Shell Biography text here.

PLACE
PHOTO
HERE

John Doe Biography text here.

Jane Doe Biography text here.