Eigen-Contrib  5.0.1
 
Loading...
Searching...
No Matches
Eigen::Bccb< Scalar_, BlockSize_, NumBlocks_ > Class Template Reference

#include <contrib/Eigen/src/StructuredMatrices/Bccb.h>

Detailed Description

template<typename Scalar_, int BlockSize_, int NumBlocks_>
class Eigen::Bccb< Scalar_, BlockSize_, NumBlocks_ >

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.

Template Parameters
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).
See also
class Circulant, makeBccb()
+ 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
 

Constructor & Destructor Documentation

◆ Bccb()

template<typename Scalar_, int BlockSize_, int NumBlocks_>
template<typename Derived>
Eigen::Bccb< Scalar_, BlockSize_, NumBlocks_ >::Bccb ( const MatrixBase< Derived > & generator)
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.

Member Function Documentation

◆ adjoint()

template<typename Scalar_, int BlockSize_, int NumBlocks_>
Bccb Eigen::Bccb< Scalar_, BlockSize_, NumBlocks_ >::adjoint ( ) const
inline
Returns
the adjoint of *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).

◆ blockSize()

template<typename Scalar_, int BlockSize_, int NumBlocks_>
Index Eigen::Bccb< Scalar_, BlockSize_, NumBlocks_ >::blockSize ( ) const
inline
Returns
the circulant block dimension n2.

◆ coeff()

template<typename Scalar_, int BlockSize_, int NumBlocks_>
Scalar Eigen::Bccb< Scalar_, BlockSize_, NumBlocks_ >::coeff ( Index row,
Index col ) const
inline
Returns
the coefficient at row row and column col.

◆ conjugate()

template<typename Scalar_, int BlockSize_, int NumBlocks_>
Bccb Eigen::Bccb< Scalar_, BlockSize_, NumBlocks_ >::conjugate ( ) const
inline
Returns
the complex conjugate of *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.

◆ determinant()

template<typename Scalar_, int BlockSize_, int NumBlocks_>
Scalar Eigen::Bccb< Scalar_, BlockSize_, NumBlocks_ >::determinant ( ) const
inline
Returns
the determinant, i.e. the product of the eigenvalues (the symbol entries). The product is accumulated in the balanced form 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.

◆ eigenvalues()

template<typename Scalar_, int BlockSize_, int NumBlocks_>
ComplexVector Eigen::Bccb< Scalar_, BlockSize_, NumBlocks_ >::eigenvalues ( ) const
inline
Returns
the eigenvalues as the column-major flattening of the symbol: eigenvalue 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.

◆ eigenvectors()

template<typename Scalar_, int BlockSize_, int NumBlocks_>
ComplexMatrix Eigen::Bccb< Scalar_, BlockSize_, NumBlocks_ >::eigenvectors ( ) const
inline
Returns
the unitary matrix of eigenvectors: column f1*n2 + f2 is the 2-D Fourier vector matching eigenvalues()[f1*n2 + f2].
Note
The eigenvector matrix is materialized as a dense N x N matrix; unlike the other methods of this class this costs O(N^2) storage.

◆ generator()

template<typename Scalar_, int BlockSize_, int NumBlocks_>
const GeneratorType & Eigen::Bccb< Scalar_, BlockSize_, NumBlocks_ >::generator ( ) const
inline
Returns
the generating array.

◆ inverse()

template<typename Scalar_, int BlockSize_, int NumBlocks_>
Bccb Eigen::Bccb< Scalar_, BlockSize_, NumBlocks_ >::inverse ( ) const
inline
Returns
the inverse of *this, itself a Bccb operator: the one generated by the 2-D inverse DFT of the entrywise-inverted symbol.
Warning
The operator must be non-singular; use solve for a pseudo-inverse solve of a rank-deficient operator.

◆ matrixU()

template<typename Scalar_, int BlockSize_, int NumBlocks_>
ComplexMatrix Eigen::Bccb< Scalar_, BlockSize_, NumBlocks_ >::matrixU ( ) const
inline
Returns
the matrix of left singular vectors 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.

◆ matrixV()

template<typename Scalar_, int BlockSize_, int NumBlocks_>
ComplexMatrix Eigen::Bccb< Scalar_, BlockSize_, NumBlocks_ >::matrixV ( ) const
inline
Returns
the matrix of right singular vectors 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.

◆ numBlocks()

template<typename Scalar_, int BlockSize_, int NumBlocks_>
Index Eigen::Bccb< Scalar_, BlockSize_, NumBlocks_ >::numBlocks ( ) const
inline
Returns
the number of blocks n1 in each block row.

◆ operator*()

template<typename Scalar_, int BlockSize_, int NumBlocks_>
template<typename Rhs>
Product< Bccb, Rhs > Eigen::Bccb< Scalar_, BlockSize_, NumBlocks_ >::operator* ( const MatrixBase< Rhs > & x) const
inline
Returns
the product expression (*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.

◆ rank()

template<typename Scalar_, int BlockSize_, int NumBlocks_>
Index Eigen::Bccb< Scalar_, BlockSize_, NumBlocks_ >::rank ( ) const
inline
Returns
the numerical rank: the number of symbol entries whose modulus is no smaller than the threshold 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.

◆ singularValues()

template<typename Scalar_, int BlockSize_, int NumBlocks_>
RealVector Eigen::Bccb< Scalar_, BlockSize_, NumBlocks_ >::singularValues ( ) const
inline
Returns
the singular values, sorted in decreasing order: the moduli of the symbol entries. The ordering is shared with matrixU and matrixV, so together they form the SVD *this = U * singularValues().asDiagonal() * V^H.

◆ solve()

template<typename Scalar_, int BlockSize_, int NumBlocks_>
template<typename Rhs>
Matrix< Scalar, RowsAtCompileTime, Rhs::ColsAtCompileTime > Eigen::Bccb< Scalar_, BlockSize_, NumBlocks_ >::solve ( const MatrixBase< Rhs > & b) const
inline
Returns
the minimum-norm least-squares solution of (*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.

◆ symbol()

template<typename Scalar_, int BlockSize_, int NumBlocks_>
ComplexArray Eigen::Bccb< Scalar_, BlockSize_, NumBlocks_ >::symbol ( ) const
inline
Returns
the symbol of the operator: the 2-D DFT of the generating array, an n2 x n1 complex array whose entries are the eigenvalues of the matrix (see eigenvalues for the ordering). Cached when the operator is large enough for products to take the FFT path, computed on the fly for small operators.

◆ transpose()

template<typename Scalar_, int BlockSize_, int NumBlocks_>
Bccb Eigen::Bccb< Scalar_, BlockSize_, NumBlocks_ >::transpose ( ) const
inline
Returns
the transpose of *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.

The documentation for this class was generated from the following file: