C++ Integration
primate’s generic API is enabled through function templates specialized using C++20 concepts. In other words, a function F requiring concept C will compile with any type T so long as T respects the constraints imposed by C. For example, generically, any type T respecting the LinearOperatorconcept shown below can be passed to the Lanczos method:
using FP = std::floating_point;
template < typename T, FP float_t = typename T::value_type >
concept LinearOperator = requires(T A, const float_t* v, float_t* o) {
{ A.matvec(v, o) }; // o = A v
{ A.shape() } -> std::convertible_to< std::pair< size_t, size_t > >;
};In english, an instance A of type T is said to support the LinearOperator concept if it has:
- A method Av \mapsto o, with signature
A.matvec(const float_t* v, float_t* o) - A method yielding (\mathrm{card}(o), \mathrm{card}(v)), with signature
A.shape() -> pair< ... >
shape() should return a pair (n,m) representing the sizes of the output and input vectors, respectively. Note in the matrix setting this corresponds to the number of rows and columns.
Other Concepts
Depending on the problem at hand, the supplied operator may need to meet other constraints. Here’s a short list additional operator concepts:
| Concept | Supports | Signature | Requires |
|---|---|---|---|
LinearOperator |
A v \mapsto o | A.matvec(v, o) |
NA |
AdjointOperator |
A^T v \mapsto o | A.rmatvec(v, o) |
LinearOperator |
AdditiveOperator |
o \gets o + \alpha Av | A.matvec_add(v, alpha, o) |
LinearOperator |
AffineOperator |
Sets t s.t. A + tB | A.set_parameter(t) |
LinearOperator |
QuadOperator |
v^T A v | A.quad(v) |
NA |
The exported methods in primate only need the minimum constraints to be satisfied to compile: if you need access to the Lanczos method, then just supporting the LinearOperator concept is sufficient. On the other hand, adding support for other constraints can optimize the efficiency of certain methods; for example, the hutch method technically only requires a LinearOperator to do trace estimation (via matvec calls), but will also compile and prefer calling quad with a QuadOperator as input. In such a situaton, if your operator has an efficient quadratic form v \mapsto v^T A v, then implementing quad may improve the performance of hutch.