ButterflyFactorizations.jl

ButterflyFactorizations

A Julia module for constructing, manipulating, and applying Butterfly Factorizations.

This package provides hierarchical low-rank approximation techniques typically used for highly oscillatory integral operators. It leverages tree-based domain partitioning and includes functionality for:

• Structured matrix-vector and matrix-matrix products.
• Algebraic operations on butterfly structures (addition, multiplication, and recompression).
• Petrov-Galerkin block matrix assembly of kernel matrices, with extension support for boundary element frameworks.

Main API

Compressors:

• PartialQR: QR-based low-rank approximation

Kernel matrices:

• BEASTKernelMatrix: BEAST boundary integral operator matrices

Example

check the test folder/Readme.md for example usage of the API, including how to construct a PetrovGalerkinBF from a BEAST operator and apply it to a vector.

See also

• BEAST.jl for boundary integral operators
• H2Trees.jl for hierarchical clustering
• BlockSparseMatrices.jl or Sparse SparseArrays.jl for sparse block storage

PetrovGalerkinBF(operator, testspace, trialspace, tree, k; kwargs...)

Constructs the complete Butterfly Factorization for a Petrov-Galerkin boundary element problem using the dictionary-based format.

This function separates the interactions based on the admissibility condition α. The near-field interactions are evaluated directly and stored as a BlockSparseMatrix. The far-field interactions are compressed block-by-block using the dictionary-based parameterized subroutine_BF, resulting in a collection of highly memory-efficient Butterfly blocks. Keep in mind that as off right now this work is still restricted to balanced trees in the sense that all leaves need to be on leaf level.

Arguments:

• operator: The integral operator (e.g., Maxwell single-layer) to be evaluated.
• testspace: The observer (test) function space.
• trialspace: The source (trial) function space.
• tree: A coupled BlockTree representing the hierarchical clustering.
• k: The physical wavenumber of the problem.

Keyword Arguments:

• Compressor: The low-rank approximation strategy (default: PartialQR()).
• tol: The relative precision tolerance for compression (default: 1e-3).
• ntasks: Number of threads to use for parallel near-field assembly.
• α: The geometric admissibility parameter (default: 2.0).

Returns:

• A PetrovGalerkinBF struct containing the near-field matrix and the array of far-field BF blocks.

PetrovGalerkinBF_mats(operator, testspace, trialspace, tree, k; kwargs...)

Constructs the complete Butterfly Factorization for a Petrov-Galerkin boundary element problem using the sparse matrix-based format.

Similar to PetrovGalerkinBF, this isolates the near-field into a BlockSparseMatrix. However, the far-field interactions are compressed using subroutine_BF_mats, which assembles the Butterfly factors (Q, R, P) explicitly as sparse block-diagonal matrices. This format requires more memory but provides significantly faster Matrix-Vector multiplication.

Arguments:

• operator: The integral operator to be evaluated.
• testspace: The observer (test) function space.
• trialspace: The source (trial) function space.
• tree: A coupled BlockTree representing the hierarchical clustering.
• k: The physical wavenumber of the problem.

Keyword Arguments:

• Compressor: The low-rank approximation strategy (default: PartialQR()).
• tol: The relative precision tolerance for compression (default: 1e-3).
• ntasks: Number of threads to use for parallel near-field assembly.
• α: The geometric admissibility parameter (default: 2.0).

Returns:

• A PetrovGalerkinBF_mats struct containing the near-field matrix and the array of far-field BF_Mats blocks.

subroutine_BF(kernelmatrix, H2Blocktree, NO, NS, k, τ; Compressor=PartialQR())

Constructs the Butterfly Factorization for a given block in a dictionary format.

This subroutine traverses the H2 tree structure from the leaf level moving up to the root of the source tree, and to the leaves of the observer tree. It computes the low-rank approximations to build the Q, R, and P factors as dictionaries. It also maps the necessary index permutations to allow for correct, independent Matrix-Vector (MV) products without permanently storing the tree.

Arguments:

• kernelmatrix: Function computing matrix entries for specified row/column indices.
• H2Blocktree: The paired source-observer tree structure.
• NO, NS: The root IDs of the observer (test) and source (trial) spaces.
• k, τ: Wavenumber (crucial for rank estimation) and precision tolerance.
• Compressor: Compression scheme for low-rank blocks (default: PartialQR).

Why Dictionary format? It is extremely memory efficient because it avoids the overhead of saving nonzero entry indices required by sparse matrices. It is deeply intuitive, makes algebraic manipulation (like butterfly multiplication or summation) flexible, and preserves the semantic structure other than a matrix format. It is ideal when still working with the tree structure is desired, such as in the case of recompression or algebraic operations on the factors. However, it requires careful handling of index permutations to ensure correct MV products, which is managed through the PermQ and PermP dictionaries.

subroutine_BF_mats(kernelmatrix, H2Blocktree, NO, NS, k, τ; Compressor=PartialQR())

Constructs the Butterfly Factorization for a given block in a sparse matrix format.

Similar to subroutine_BF, it traverses the H2 tree structure to compute the Q, R, and P factors, but specifically pieces them together into sparse block-diagonal matrices. It also returns single continuous permutation vectors (PermP, PermQ) to map interactions across the entire space.

Why Matrix format? While slightly less memory efficient than the dictionary format (due to sparsity tracking overhead), it allows for dramatically faster direct Matrix-Vector applications using standard linear algebra methods. It also provides a clear visual and algebraic representation of the overall block structure, which is invaluable for debugging.

PartialQR <: Abstractcompressor

A type representing the Partial QR compression strategy for low-rank approximations.

BEASTKernelMatrix

A wrapper structure for BEAST boundary integral operator matrices.

mulBFs(BF_1::BF, BF_2::BF, τ::Float64)

Algebraically multiplies two Butterfly Factorization (BF) structures and recompresses the result using the truncation tolerance τ.

Both factorizations must have the same number of levels, and the source dimensions of BF_1 must match the observer dimensions of BF_2. Additionally, the resulting structure is purely algebraic and may lose its direct physical interpretation, similar to multiplying two dense matrices directly. The function constructs an intermediate messenger structure to hold the products of the factors, then iteratively multiplies and recompresses the factors level by level, ultimately returning a new BF that represents the product of the two input factorizations. The resulting BF maintains the same hierarchical structure but with potentially reduced ranks in the R factors, leading to improved efficiency in storage and matrix-vector products while preserving the overall accuracy within the specified tolerance. In terms of storage future work has to include more aggressive recompression strategies to prevent the intermediate factors from growing too large.

There are exactly two cases to consider when adding Butterfly factorizations (BFs). In both cases, the corresponding matrix blocks are of equal dimensions.

Case 1 (Same Source and Observer Clusters): The BFs map between the same observer and source clusters. Because the underlying compression scheme relies on a single tree, the physical degrees of freedom (DoF) match. Thus, the Q and P matrices can be appropriately combined without disturbing the underlying physics, applying the addition described in the literature.

Case 2 (Disjoint Source and Observer Clusters): The BFs represent disjoint source and observer clusters. A direct spatial concatenation cannot occur, so the BFs are simply joined into a new structure. Note that this resulting structure is purely algebraic and loses its direct physical interpretation, similar to adding the two dense matrices directly.

recompress_BF(Butterfly::BF, τ)

Recompresses a structural Butterfly Factorization (BF) by extracting its algebraic factors, recompressing them with tolerance τ, and restructuring the output back into a BF. This process involves two main steps: first, the right factors are recompressed using a QR-based approach, and then the left factors are recompressed by transposing the structure, applying the same right recompression, and transposing back. The resulting BF maintains the same hierarchical structure but with potentially reduced ranks in the R factors, leading to improved efficiency in storage and matrix-vector products while preserving the overall accuracy within the specified tolerance. Any of the algebraic operations is only supported for the Dictionary versions of the Butterflies, as the matrix-based format is not designed for algebraic manipulations and would require a complete restructuring of the underlying data representation to support such operations effectively.

apply_BF(Butterfly::BF, v::AbstractVector{ComplexF64})

Applies the sequence of Butterfly Factorization factors (Q, R, and P) to a vector v and returns the resulting vector. Do note that this function is optimized for the structure of BF and uses a "ping-pong" strategy to minimize memory allocations during the sequential application of the R factors. The indices of the dictionaries are matching the row and column indices of the skeleton associated with the underlying Matrix structure thus we can perform a block by block matrix vector product using the well known algebra for matrices, while retaining the semantical structure of the corresponding clusters.