![]() |
Eigen-Contrib
5.0.1
|
#include <contrib/Eigen/src/StructuredMatrices/Bccb.h>
A block circulant matrix with circulant blocks (BCCB), the matrix of a two-dimensional circular convolution, represented by its n2 x n1 generating array.
A BCCB matrix is the N x N matrix, N = n1*n2, that is circulant at two levels: it is an n1 x n1 block circulant whose n2 x n2 blocks are themselves circulant. With the generating array G (column k holds the first column of the k-th block), entry (i,j) with i = b1*n2 + i2, j = c1*n2 + j2 equals
\[ C_{i,j}=G_{(i_2-j_2)\bmod n_2,\,(b_1-c_1)\bmod n_1}. \]
On a column-major reshaped vector, \(C\,\operatorname{vec}(X) =\operatorname{vec}(G\mathbin{\circledast}X)\), where \(\circledast\) denotes 2-D circular convolution.
BCCB matrices are diagonalized by the 2-D discrete Fourier transform \( F_{n_1} \otimes F_{n_2} \) ([1], [2]): the operator's symbol – the 2-D DFT of G – holds the eigenvalues. Products reuse that symbol or, for an awkward transform size, an equivalent cached padded-embedding symbol. This yields O(N log N) products (operator*), an O(N log N) direct (pseudo-inverse) solve (solve), and closed-form factorizations: the eigendecomposition (eigenvalues, eigenvectors) and the SVD (singularValues, matrixU, matrixV) in the 2-D Fourier basis, plus rank, inverse and determinant. The class is closed under transpose, conjugate and adjoint, which reuse the cached symbols. BCCB operators are the workhorse of image deblurring with periodic boundary conditions, and the natural preconditioners for two-level Toeplitz (BTTB) systems [2].
The operator stores its own copy of the generating array and derives from EigenBase. Because operator* returns an Eigen product expression, a Bccb also drops into the matrix-free iterative solvers, and it can be assigned to a dense matrix when an explicit representation is needed. As with any matrix-free operator, the iterative solvers must be instantiated with IdentityPreconditioner (e.g. ConjugateGradient<Bccb<double>,Lower|Upper,IdentityPreconditioner>): the default preconditioners read individual coefficients through col() or InnerIterator, which the structured operators do not expose.
Spectral operations use FFTs of the exact sizes n1 and n2. Products use those sizes when both are 5-smooth; otherwise an equivalent per-axis circulant embedding pads each awkward dimension to a 5-smooth size, avoiding the default kissfft backend's slow generic butterfly for large prime factors.
| Scalar_ | the scalar type, real or complex. |
| BlockSize_ | the circulant block dimension n2 at compile time, or Dynamic (the default). |
| NumBlocks_ | the number of blocks n1 at compile time, or Dynamic (the default). |
Inheritance diagram for Eigen::Bccb< Scalar_, BlockSize_, NumBlocks_ >:Public Member Functions | |
| Bccb | adjoint () const |
| template<typename Derived> | |
| Bccb (const MatrixBase< Derived > &generator) | |
| Index | blockSize () const |
| Scalar | coeff (Index row, Index col) const |
| Bccb | conjugate () const |
| Scalar | determinant () const |
| ComplexVector | eigenvalues () const |
| ComplexMatrix | eigenvectors () const |
| const GeneratorType & | generator () const |
| Bccb | inverse () const |
| ComplexMatrix | matrixU () const |
| ComplexMatrix | matrixV () const |
| Index | numBlocks () const |
| template<typename Rhs> | |
| Product< Bccb, Rhs > | operator* (const MatrixBase< Rhs > &x) const |
| Index | rank () const |
| RealVector | singularValues () const |
| template<typename Rhs> | |
| Matrix< Scalar, RowsAtCompileTime, Rhs::ColsAtCompileTime > | solve (const MatrixBase< Rhs > &b) const |
| ComplexArray | symbol () const |
| Bccb | transpose () const |
|
inlineexplicit |
Builds a BCCB matrix from its generating array generator: column k is the first column of the k-th circulant block.
When the matrix is large enough for products to take the FFT path, the 2-D DFT of the array – the eigenvalues of the matrix – and any padded product symbol are computed here. Subsequent products and solves reuse the applicable cached transform.
|
inline |
*this, itself a Bccb operator. The cached symbols, when present, are reused: the symbols of the adjoint are the elementwise conjugates of the symbols (the eigenvalues conjugate while the 2-D Fourier eigenbasis stays fixed).
|
inline |
n2.
|
inline |
|
inline |
*this, itself a Bccb operator. The cached symbols, when present, are reused: the symbols of the conjugate are the conjugated two-dimensional index reversals of the symbols.
|
inline |
m * 2^e (the split fraction/exponent determinant convention of LINPACK's xGEDI [4]) – every factor and the running product are renormalized to unit magnitude with the power of two tracked separately – so the partial products can neither overflow nor underflow when the determinant itself is representable, whatever the ordering of large and small eigenvalues. For a real operator the product is real up to roundoff, and its real part is returned.
|
inline |
f1*n2 + f2 is symbol()(f2, f1), and its (unit-norm) eigenvector is the Kronecker product of the 1-D Fourier vectors of frequencies f1 and f2, i.e. column f1*n2 + f2 of eigenvectors. Every BCCB matrix is diagonalized by this same 2-D Fourier basis.
|
inline |
f1*n2 + f2 is the 2-D Fourier vector matching eigenvalues()[f1*n2 + f2]. N x N matrix; unlike the other methods of this class this costs O(N^2) storage.
|
inline |
|
inline |
|
inline |
U: column t is the 2-D Fourier vector of the t-th largest symbol entry, scaled by its phase (phase 1 for a zero entry). Dense N x N, see the note in eigenvectors.
|
inline |
V: column t is the 2-D Fourier vector of the t-th largest symbol entry. Dense N x N, see the note in eigenvectors.
|
inline |
n1 in each block row.
|
inline |
(*this) * x, evaluated through a fast 2-D-FFT-based matrix-vector product. The expression carries the default product tag, so assigning it behaves like any dense product: a temporary resolves aliasing between the destination and x, and .noalias() skips it.
|
inline |
N * epsilon * max|symbol|, clamped from below by the smallest normal number like SVDBase::rank(). This is the same threshold solve uses to decide which Fourier components to invert, and the comparison is strict like SVDBase's, so an entry sitting exactly on the threshold still counts as non-zero.
|
inline |
*this = U * singularValues().asDiagonal() * V^H.
|
inline |
(*this) * x = b, computed directly in the 2-D Fourier domain. Symbol entries whose modulus reaches the rank threshold (see rank) are inverted; the remaining ones are treated as exact zeros, so the result is the pseudo-inverse applied to b. For a non-singular operator this is the exact solution. Supports multiple right-hand sides.
|
inline |
|
inline |
*this, itself a Bccb operator: the one generated by the array index-reversed in both dimensions. The cached symbols, when present, are reused – the symbols of the transpose are their two-dimensional index reversals (embedding a generator commutes with index-reversing it, per axis, so the padded product symbol follows the same rule at its own grid) – so no FFT is recomputed.