torch_vn_algebra.jordan¶
Source: torch_vn_algebra/jordan.py. Signatures and docstrings as in the source.
Module documentation¶
Jordan structure of finite-dimensional von Neumann algebras: Jordan products, the quadratic representation,
Jordan subalgebras of Hermitian matrices and spin factors.
Definitions
-----------
A **Jordan algebra** is a real vector space :math:`J` with a commutative bilinear product :math:`x\circ y` that
satisfies the **Jordan identity** :math:`(x^2\circ y)\circ x=x^2\circ(y\circ x)`. The self-adjoint part
:math:`M_{sa}` of :math:`M=\bigoplus_c M_{k_c}` (and every :math:`\operatorname{Herm}_n`) is a Jordan algebra with
.. math:: x\circ y=\tfrac12(xy+yx).
A Jordan algebra isomorphic to a Jordan subalgebra of :math:`(A,\circ)` for an associative algebra :math:`A` is
**special**; every Jordan subalgebra of :math:`\operatorname{Herm}_n` is special by definition. On :math:`M_{sa}`:
* the multiplication operator :math:`L(y)x=y\circ x`;
* the **quadratic representation** :math:`P(y)=2L(y)^2-L(y^2)`, i.e. :math:`P(y)x=yxy` (expand
:math:`2y\circ(y\circ x)-y^2\circ x`);
* the **Jordan triple product**
:math:`\{x,y,z\}=\tfrac12(xyz+zyx)=(x\circ y)\circ z+(y\circ z)\circ x-(x\circ z)\circ y`, with
:math:`\{x,y,x\}=P(x)y`;
* Jordan powers :math:`x^{\circ1}=x`, :math:`x^{\circ(k+1)}=x\circ x^{\circ k}`. They coincide with the associative
powers, because :math:`x` commutes with :math:`x^k` (power associativity), so :func:`jordan_power` is a matrix power.
Identities checked in ``tests/test_jordan.py`` (they hold in every Jordan algebra; in a special one they reduce
to associativity): the Jordan identity; :math:`P(y)=2L(y)^2-L(y^2)`; the **fundamental formula**
:math:`P(P(y)x)=P(y)P(x)P(y)`; :math:`\operatorname{Tr}(P(y)x)=\operatorname{Tr}(y^2x)`; :math:`P(y)` is invertible
(as a linear map of :math:`M_{sa}`) iff :math:`y` is invertible, and then :math:`P(y)^{-1}=P(y^{-1})`
(on a block :math:`M_k` of rank-:math:`r` :math:`y`, :math:`P(y)` has rank :math:`r^2` on
:math:`\operatorname{Herm}_k(\mathbb C)` and :math:`r(r+1)/2` on :math:`\operatorname{Herm}_k(\mathbb R)`).
The **positive cone** of a Jordan algebra is the cone of squares :math:`\{x^2\}`; in :math:`M_{sa}` it is the cone of
positive semidefinite operators. :math:`P(y)` preserves it, but the Jordan product of two positive elements need not be
positive (:math:`p=\mathrm{diag}(1,0)`, :math:`q=\tfrac12\begin{pmatrix}1&1\\1&1\end{pmatrix}`:
:math:`\det(p\circ q)=-1/16`).
Theorems (stated for orientation; the code does not depend on them)
--------------------------------------------------------------------
* **Jordan-von Neumann-Wigner classification** [JvNW34]. A finite-dimensional *formally real* Jordan algebra
(:math:`\sum_ix_i^2=0\Rightarrow` all :math:`x_i=0`; equivalently *Euclidean*: it carries an inner product with
:math:`\langle x\circ y,z\rangle=\langle y,x\circ z\rangle`) is a direct sum of simple ones, and every simple one is
isomorphic to exactly one of
:math:`\mathbb R`; the spin factor :math:`V_m=\mathbb R\oplus\mathbb R^m` (:math:`m\ge2`, see :class:`SpinFactor`);
:math:`\operatorname{Herm}_r(\mathbb K)` for :math:`\mathbb K\in\{\mathbb R,\mathbb C,\mathbb H\}`, :math:`r\ge3`;
the Albert algebra :math:`\operatorname{Herm}_3(\mathbb O)` (27-dimensional). Rank 2 gives no new algebras:
:math:`\operatorname{Herm}_2(\mathbb R)\cong V_2`, :math:`\operatorname{Herm}_2(\mathbb C)\cong V_3`,
:math:`\operatorname{Herm}_2(\mathbb H)\cong V_5`, :math:`\operatorname{Herm}_2(\mathbb O)\cong V_9`.
All of them are special except the Albert algebra [Alb34].
* **Koecher-Vinberg theorem** [FK94]. The interiors of the cones of squares of Euclidean Jordan algebras are exactly
the *symmetric cones*: open convex cones in a Euclidean space that are self-dual and homogeneous (their linear
automorphism group acts transitively). For :math:`M_{sa}` this is the cone of positive definite operators, for
:math:`V_m` the (open) Lorentz cone :math:`a>|u|`.
* **Jordan isomorphisms** [Kad51]. A Jordan :math:`*`-isomorphism of von Neumann algebras is the sum of a
:math:`*`-isomorphism and a :math:`*`-anti-isomorphism (on complementary central summands). For
:math:`M=\bigoplus_cM_{k_c}` this means: every Jordan automorphism of :math:`M_{sa}` permutes sectors of equal size
and acts on each block as :math:`x\mapsto uxu^*` or as :math:`x\mapsto ux^{T}u^*` with a unitary :math:`u`. The
transpose (:func:`transpose_map`) preserves :math:`\circ` and :math:`P` but reverses products; it is positive but
not completely positive (its Choi matrix on :math:`M_k`, :math:`k\ge2`, is the swap, with eigenvalue :math:`-1`)
[Cho75].
Jordan subalgebras of Hermitian matrices
----------------------------------------
:func:`jordan_algebra_generated` computes the smallest real subspace :math:`J\subseteq\operatorname{Herm}_n` that
contains the generators (and :math:`1` if ``unital``) and is closed under :math:`\circ`. Given the complex
:math:`*`-algebra :math:`A\subseteq M_n(\mathbb C)` generated by the same set (as in
:mod:`~torch_vn_algebra.decompose`), always :math:`J\subseteq A_{sa}` and
:math:`\dim_{\mathbb R}A_{sa}=\dim_{\mathbb C}A`; equality :math:`J=A_{sa}` holds iff :math:`J` is closed under
:math:`(x,y)\mapsto\frac i2[x,y]` (then :math:`J+iJ` is a :math:`*`-algebra, since
:math:`xy=x\circ y-i\cdot\frac i2[x,y]`). :func:`star_algebra_self_adjoint_part` computes :math:`A_{sa}` as the closure
under both products. Examples (``tests/test_jordan.py``):
* real symmetric generators: :math:`J=\operatorname{Herm}_n(\mathbb R)` (generic pair), :math:`\dim J=n(n+1)/2`, but
:math:`A=M_n(\mathbb C)`, :math:`\dim A_{sa}=n^2`. (The *real* :math:`*`-algebra they generate is
:math:`M_n(\mathbb R)`, whose self-adjoint part is :math:`J`; the comparison here is always with the complex one.)
* anticommuting symmetries :math:`s_1,\dots,s_m` (:func:`anticommuting_symmetries`): :math:`s_i\circ s_j=\delta_{ij}1`,
so :math:`J=\operatorname{span}\{1,s_1,\dots,s_m\}\cong V_m`, :math:`\dim J=m+1`. In the representation on
:math:`\mathbb C^{2^{\lfloor m/2\rfloor}}` used here, :math:`A` is the full matrix algebra (:math:`m\ge2`), so
:math:`J=A_{sa}` only for :math:`m=3` (:math:`V_3=\operatorname{Herm}_2(\mathbb C)`, the Pauli matrices) and
:math:`m=1`.
* the self-adjoint part of a :math:`*`-algebra :math:`A=U(\bigoplus_cM_{k_c}\otimes1_{m_c})U^*`: three generic
Hermitian elements of A generate :math:`J=A_{sa}`, and so do two if no :math:`k_c=2`. Two Hermitian
:math:`2\times2` matrices are simultaneously unitarily similar to real symmetric ones (their Pauli vectors span a
plane), so in a block :math:`M_2` two elements generate only :math:`\operatorname{Herm}_2(\mathbb R)\cong V_2`
(dimension 3) as a Jordan algebra, although they generate :math:`M_2(\mathbb C)` as a :math:`*`-algebra.
* one Hermitian generator :math:`x`: :math:`J` = real polynomials in :math:`x` = span of its spectral projections;
:math:`\dim J` = number of distinct eigenvalues (``unital=True``) or of distinct non-zero eigenvalues
(``unital=False``: polynomials without constant term; if :math:`0\notin\operatorname{spec}x` they contain
:math:`1` anyway).
*Algorithm* (randomised closure, in the spirit of the randomised algorithm of decompose). Hermitian matrices are
identified isometrically with :math:`\mathbb R^N`, :math:`N=n^2` (:math:`n(n+1)/2` for real symmetric ones), by the
diagonal and :math:`\sqrt2` times the real and imaginary parts of the upper triangle, so that the Euclidean inner
product is :math:`\operatorname{Tr}(ab)` (Hilbert-Schmidt). Starting from an orthonormal basis :math:`Q` of the span
of the generators (and :math:`1`), every round draws a random unit element :math:`r` of the current span :math:`V`
and adds the part of :math:`L(r)V=\{r\circ v_j\}` orthogonal to :math:`V` (projection twice, then an SVD as rank
decision). If :math:`V` is not closed, :math:`\{r\in V: L(r)V\subseteq V\}` is a proper linear subspace of :math:`V`
(it misses some :math:`a` with :math:`a\circ b\notin V`), so a random :math:`r` enlarges :math:`V` with probability 1;
a round that adds nothing is therefore a closure certificate, which is confirmed with fresh random pairs
:math:`r_k\circ s_k` (``diagnostics['closure_defect']``). Typically :math:`\dim V` doubles per round (for one
generator, :math:`r` of degree :math:`d` gives degree :math:`2d`), so there are :math:`O(\log N)` rounds; in any
case every productive round adds a dimension, so there are at most :math:`N` of them.
*Cost.* A round costs :math:`d` products of :math:`n\times n` matrices, :math:`O(dn^3)`, and the orthogonalisation,
:math:`O(Nd\,k)` for :math:`k` candidates; the SVD runs on the candidates whose residual is not negligible, at most
:math:`O(Nk^2)`. With :math:`d,k\le N\le n^2` a round costs at most :math:`O(n^6)`, and as :math:`\dim V` grows
geometrically the last rounds dominate: :math:`O(n^6)` in total for a generic generating set, memory
:math:`O(n^4)` (the basis); :math:`J=\operatorname{Herm}_{64}(\mathbb C)` (:math:`d=4096`) takes seconds to a minute
on a CPU and about 0.3 GB.
Tolerances
----------
As in :mod:`~torch_vn_algebra.decompose`, with :math:`\varepsilon` the machine epsilon of the working precision
(double by default; single for float32/complex64 input) and ``tol`` defaulting to :math:`\varepsilon^{1/2}`:
* generators are normalised to unit Frobenius norm; Hermitian and anti-Hermitian parts (the latter enter as
:math:`(S-S^*)/2i`, so that a non-Hermitian generator :math:`S=H_1+iH_2` contributes :math:`H_1,H_2`, the same
convention as decompose) of relative size below ``tol`` are dropped as rounding noise;
* rank decision: residuals (after projection onto :math:`V`) of unit-norm candidates with singular values above
:math:`\theta=\max(\mathrm{tol},10^3n\varepsilon)` are new directions. The candidates :math:`r\circ v_j` have
:math:`\|r\|_F=\|v_j\|_F=1`, hence norm at most 1 and rounding errors of order :math:`n\varepsilon`; features of
relative size below :math:`\theta` (e.g. eigenvalues of a single generator closer than about :math:`\theta`) are
merged. A direction added with singular value :math:`\sigma` is accurate to about :math:`n\varepsilon/\sigma`
(Wedin's theorem), and its products with the basis have components of that size outside the true algebra;
:math:`\theta` is therefore raised to :math:`10\,n\varepsilon/\sigma_{\min}` when this is larger (without it, a
nearly degenerate input, e.g. eigenvalues :math:`1` and :math:`1+10^{-7}`, would cascade to spurious dimensions).
This does not model the compounding of errors over rounds, a worst case not observed in practice; the closure
certificate below checks the result. Singular values within a factor :math:`\varepsilon^{\mp1/8}` of the
current :math:`\theta` trigger a :class:`ClosureWarning` (the dimension is ill-determined);
* closure: :meth:`JordanAlgebra.closure_defect` is :math:`\max\|P_{V^\perp}(a\circ b)\|_F` over pairs of orthonormal
basis elements (exhaustive) or over random unit pairs (a randomised lower bound that is non-zero almost surely if
:math:`V` is not closed); a subspace counts as closed if the defect is at most :math:`\theta`.
Operator-level functions
------------------------
:func:`jordan_product`, :func:`quadratic_representation`, :func:`jordan_triple`, :func:`jordan_power`,
:func:`jordan_multiplication` and :func:`transpose_map` act on operators of a :class:`~torch_vn_algebra.TypeIAlgebra`
(tensors ``(batch, C, k_max, k_max)``, batch 1 broadcasting against batch B). They are lazy like the library's other
derived operators, differentiable, run as batched matrix products over all sectors at once (the zero padding stays
zero) and read nothing back to the host. Property flags are set only where they are mathematically guaranteed.
References
----------
* [JvNW34] P. Jordan, J. von Neumann, E. Wigner, *On an algebraic generalization of the quantum mechanical
formalism*, Ann. of Math. 35 (1934) 29-64.
* [Alb34] A. A. Albert, *On a certain algebra of quantum mechanics*, Ann. of Math. 35 (1934) 65-73.
* [Kad51] R. V. Kadison, *Isometries of operator algebras*, Ann. of Math. 54 (1951) 325-338.
* [Cho75] M.-D. Choi, *Completely positive linear maps on complex matrices*, Linear Algebra Appl. 10 (1975) 285-290.
* [FK94] J. Faraut, A. Korányi, *Analysis on Symmetric Cones*, Oxford University Press, 1994.
* K. McCrimmon, *A Taste of Jordan Algebras*, Springer, 2004 (quadratic representation, fundamental formula).
* P.-Å. Wedin, *Perturbation bounds in connection with singular value decomposition*, BIT 12 (1972) 99-111.
ClosureWarning¶
class ClosureWarning()
A rank decision of the subalgebra closure was close to its threshold (ill-conditioned input).
jordan_product¶
def jordan_product(A: Operator, B: Operator) -> Operator
Jordan product :math:`A\circ B=\frac12(AB+BA)` (lazy), batch by batch and sector by sector.
A, B: operators of the same algebra, matrices (batch, C, k_max, k_max) (batch 1 broadcasts).
Returns an operator of the same shape. Flags: self-adjoint (and normal) if A and B are; A∘A = A² is positive if A
is self-adjoint. The product of two positive operators is *not* positive in general (module docstring), so no
positivity flag is set otherwise.
jordan_multiplication¶
def jordan_multiplication(Y: Operator) -> Callable[[Operator], Operator]
The multiplication operator :math:`L(Y):X\mapsto Y\circ X` as a callable (lazy results, see jordan_product).
quadratic_representation¶
def quadratic_representation(Y: Operator, X: Optional[Operator]=None) -> Union[Operator, Callable]
Quadratic representation :math:`P(Y)X=YXY^*` (lazy); for self-adjoint Y this is the Jordan-algebra map
:math:`P(Y)X=YXY=2Y\circ(Y\circ X)-Y^2\circ X`. For a non-self-adjoint Y it is the congruence
:math:`X\mapsto YXY^*` (not the quadratic representation of the complex Jordan algebra :math:`(M,\circ)`, which is
:math:`X\mapsto YXY=` ``jordan_triple(Y, X, Y)``).
Y, X: operators of the same algebra, (batch, C, k_max, k_max). Without X, returns the map :math:`X\mapsto P(Y)X`.
Flags: self-adjoint if X is (for every Y: :math:`(YXY^*)^*=YX^*Y^*`); **positive if X is positive** (for every Y):
:math:`X=X^{1/2}X^{1/2}` with self-adjoint :math:`X^{1/2}`, hence
:math:`YXY^*=(YX^{1/2})(YX^{1/2})^*\ge0`, and for a vector :math:`\xi`,
:math:`\langle\xi,YXY^*\xi\rangle=\|X^{1/2}Y^*\xi\|^2\ge0`. So P(Y) maps the positive cone into itself;
invertible if X and Y are (and P(Y) is a bijection of the cone for invertible Y, with inverse
:math:`P(Y^{-1})`).
jordan_triple¶
def jordan_triple(X: Operator, Y: Operator, Z: Operator) -> Operator
Jordan triple product :math:`\{X,Y,Z\}=\frac12(XYZ+ZYX)` (lazy). For self-adjoint arguments
:math:`\{x,y,z\}=(x\circ y)\circ z+(y\circ z)\circ x-(x\circ z)\circ y` and :math:`\{x,y,x\}=P(x)y`. (The
:math:`JB^*`-triple product of :math:`M` is :math:`\frac12(XY^*Z+ZY^*X)`; it agrees for self-adjoint Y.)
X, Y, Z: operators of the same algebra, (batch, C, k_max, k_max). Flags: self-adjoint if all three are;
positive if X is Z, X is self-adjoint and Y positive (then :math:`\{X,Y,X\}=P(X)Y`).
jordan_power¶
def jordan_power(A: Operator, p: int) -> Operator
Jordan power :math:`A^{\circ p}` (:math:`A^{\circ0}=1`, :math:`A^{\circ(k+1)}=A\circ A^{\circ k}`) for an integer
:math:`p\ge0`, lazy. By power associativity (A commutes with :math:`A^k`) it equals the associative power
:math:`A^p`, computed by repeated squaring (:math:`O(\log p)` batched products). Shape (batch, C, k_max, k_max).
Negative powers: ``A.inverse()`` (the Jordan inverse of an invertible self-adjoint element is its inverse) or
``A.power(p)`` for positive A.
Flags: self-adjoint, normal, invertible and projection are inherited; positive if A is positive, or if A is
self-adjoint and p even.
transpose_map¶
def transpose_map(A: Operator) -> Operator
Blockwise transpose :math:`A=\bigoplus_cA_c\mapsto\bigoplus_cA_c^T` (lazy; transpose in the standard basis of every
block, not the adjoint). Shape (batch, C, k_max, k_max).
It is complex linear, unital and preserves adjoints, Jordan products, the quadratic representation and spectra
(hence positivity), but reverses products: :math:`(AB)^T=B^TA^T`. So it is a Jordan :math:`*`-automorphism and a
:math:`*`-anti-automorphism, not a :math:`*`-automorphism (for blocks of size :math:`\ge2`). It is positive but
not completely positive: its Choi matrix :math:`\sum_{ij}E_{ij}\otimes E_{ji}` on :math:`M_k` is the swap, with
eigenvalue :math:`-1` for :math:`k\ge2` [Cho75]. In a real algebra the transpose is the adjoint, so it is the
identity on self-adjoint operators.
Flags: all property flags (and a known spectrum) are inherited.
JordanAlgebra¶
class JordanAlgebra(basis: torch.Tensor, field: str, unital: bool, kind: str, threshold: float, diagnostics: Optional[Dict]=None, seed: Optional[int]=0)
A Jordan subalgebra :math:`J\subseteq\operatorname{Herm}_n(\mathbb K)` (:math:`\mathbb K=\mathbb R` or
:math:`\mathbb C`), as returned by :func:`jordan_algebra_generated` and :func:`star_algebra_self_adjoint_part`.
Attributes
----------
basis : (d, n, n) Hermitian matrices (real symmetric for field 'real'), orthonormal for
:math:`\langle a,b\rangle=\operatorname{Tr}(ab)`.
dim : d = :math:`\dim_{\mathbb R}J`. n : matrix size. field : 'real' or 'complex'.
unital : whether :math:`1` was adjoined. kind : 'jordan' (closure under ∘) or 'star' (closure under ∘ and
:math:`\frac i2[\cdot,\cdot]`: the self-adjoint part of a complex :math:`*`-algebra).
threshold : the final rank threshold :math:`\theta` (module docstring, 'Tolerances'); closure checks compare
with it.
diagnostics : dict (rounds, verifications, closure_defect, max_residual, min_kept, max_dropped, ...).
JordanAlgebra.coordinates¶
def coordinates(self, x: torch.Tensor) -> torch.Tensor
Real coefficients :math:`\operatorname{Re}\operatorname{Tr}(b_i^*x)`: x (..., n, n) -> (..., d).
JordanAlgebra.project¶
def project(self, x: torch.Tensor) -> torch.Tensor
Orthogonal projection (real Hilbert-Schmidt inner product) of x (..., n, n) onto J.
JordanAlgebra.residual¶
def residual(self, x: torch.Tensor) -> torch.Tensor
:math:`\|x-\mathrm{proj}_J(x)\|_F/\|x\|_F` (0 for x = 0), shape (...).
JordanAlgebra.contains¶
def contains(self, x: torch.Tensor, tol: Optional[float]=None) -> torch.Tensor
x in J up to the relative tolerance (default: the threshold of J). Bool tensor of shape (...).
JordanAlgebra.closure_defect¶
def closure_defect(self, product: str='jordan', exhaustive: Optional[bool]=None, samples: int=32, seed: Optional[int]=None) -> float
:math:`\max\|P_{J^\perp}(a\cdot b)\|_F` for the product 'jordan' (:math:`a\circ b`) or 'lie'
(:math:`\frac i2[a,b]`), over all pairs of basis elements (exhaustive, :math:`d(d+1)/2` products) or over
``samples`` random unit pairs of J (a randomised lower bound of the sup over unit pairs, non-zero almost
surely if J is not closed). Default: exhaustive when :math:`d(d+1)n^3/2\le10^{10}`.
JordanAlgebra.is_closed¶
def is_closed(self, exhaustive: Optional[bool]=None) -> bool
J ∘ J ⊆ J up to the threshold (closure_defect('jordan') <= threshold).
JordanAlgebra.is_self_adjoint_part¶
def is_self_adjoint_part(self, exhaustive: Optional[bool]=None) -> bool
Whether :math:`J=A_{sa}` for the complex :math:`*`-algebra :math:`A` generated by J, i.e. whether J is closed
under :math:`\frac i2[x,y]` (module docstring). False e.g. for :math:`\operatorname{Herm}_n(\mathbb R)`,
:math:`n\ge2`, and for spin factors :math:`V_m`, :math:`m\notin\{1,3\}`.
JordanAlgebra.star_algebra¶
def star_algebra(self, **kwargs) -> 'JordanAlgebra'
The self-adjoint part of the complex *-algebra generated by J (star_algebra_self_adjoint_part).
jordan_algebra_generated¶
def jordan_algebra_generated(generators: MatrixLike, *, unital: bool=True, field: Optional[str]=None, tol: Optional[float]=None, seed: Optional[int]=0, precision: Optional[str]=None, device=None, max_rounds: Optional[int]=None, samples: int=8) -> JordanAlgebra
The Jordan algebra generated by matrices: the smallest real subspace :math:`J\subseteq\operatorname{Herm}_n`
containing the Hermitian generators (and :math:`1` if ``unital``) with :math:`J\circ J\subseteq J`.
generators: (J, n, n) tensor or list of (n, n) matrices; a non-Hermitian :math:`S=H_1+iH_2` contributes
:math:`H_1` and :math:`H_2` (module docstring). Zero generators are allowed.
unital: adjoin the identity (default). Without it, e.g. one generator gives the polynomials without constant term.
field: 'real' (coordinates of real symmetric matrices, :math:`N=n(n+1)/2`; needs real symmetric generators),
'complex' (:math:`N=n^2`) or None (real iff the generators are real symmetric). J does not depend on it
(products of real symmetric matrices are real symmetric); only the cost and the dtype of the basis do.
tol: relative tolerance (default :math:`\sqrt\varepsilon`), seed: random elements (default 0: reproducible;
the result does not depend on it except for rounding), precision: 'double' (default) or 'single' (default for
float32/complex64 input), device: of the computation (default: of the input), max_rounds: bound on the
closure rounds (default :math:`2N+2`), samples: random pairs of the closure certificate.
Returns a :class:`JordanAlgebra` with an orthonormal basis (d, n, n). Algorithm, cost and tolerances: module
docstring. Issues a :class:`ClosureWarning` if a rank decision was close to the threshold.
star_algebra_self_adjoint_part¶
def star_algebra_self_adjoint_part(generators: MatrixLike, *, unital: bool=True, tol: Optional[float]=None, seed: Optional[int]=0, precision: Optional[str]=None, device=None, max_rounds: Optional[int]=None, samples: int=8) -> JordanAlgebra
The self-adjoint part :math:`A_{sa}` of the complex :math:`*`-algebra :math:`A` generated by matrices (with
:math:`1` if ``unital``: then :math:`A` is the algebra of :func:`~torch_vn_algebra.decompose.wedderburn` and
:math:`\dim_{\mathbb R}A_{sa}=\dim_{\mathbb C}A=\sum_ck_c^2`). Computed as the closure of the Hermitian generators
under :math:`\circ` and :math:`\frac i2[\cdot,\cdot]` (Lie-Jordan closure), always in complex coordinates. Proof
that this is :math:`A_{sa}`: the closure K lies in :math:`A_{sa}`; :math:`K+iK` is closed under products
(:math:`xy=x\circ y-i\cdot\frac i2[x,y]`) and adjoints, contains the generators, hence contains A, and its
self-adjoint part is K. Arguments as in :func:`jordan_algebra_generated`; returns a :class:`JordanAlgebra` with
kind 'star' and field 'complex'.
is_jordan_subalgebra¶
def is_jordan_subalgebra(basis: MatrixLike, *, tol: Optional[float]=None, exhaustive: Optional[bool]=None, precision: Optional[str]=None) -> bool
Whether the real span of Hermitian matrices (d, n, n) (any spanning set) is closed under :math:`\circ`. The span is
orthonormalised (rank threshold :math:`\theta` as in the module docstring) and
:meth:`JordanAlgebra.closure_defect` is compared with :math:`\theta` (exhaustive over basis pairs when
:math:`d(d+1)n^3/2\le10^{10}`, else randomised). Raises ValueError for non-Hermitian input.
compare_with_star_algebra¶
def compare_with_star_algebra(generators: MatrixLike, *, unital: bool=True, **kwargs) -> Dict
The Jordan algebra J and the self-adjoint part :math:`A_{sa}` of the complex :math:`*`-algebra generated by the
same matrices: ``{'jordan_dim', 'star_dim', 'equal', 'jordan', 'star'}``. Always :math:`J\subseteq A_{sa}`, so
'equal' (``jordan_dim == star_dim``) means :math:`J=A_{sa}`; otherwise :math:`J\subsetneq A_{sa}`. (Every Jordan
subalgebra of :math:`\operatorname{Herm}_n` is special; :math:`J\ne A_{sa}` is a statement about this embedding.)
anticommuting_symmetries¶
def anticommuting_symmetries(m: int, dtype: torch.dtype=torch.complex128, device=None) -> torch.Tensor
Hermitian :math:`s_1,\dots,s_m` with :math:`s_is_j+s_js_i=2\delta_{ij}1`, shape (m, D, D) with
:math:`D=2^{\lfloor m/2\rfloor}` (:math:`D=2` for :math:`m=1`): with :math:`p=\lfloor m/2\rfloor` qubits, the
Jordan-Wigner strings :math:`s_{2j-1}=Z^{\otimes(j-1)}\otimes X\otimes1`,
:math:`s_{2j}=Z^{\otimes(j-1)}\otimes Y\otimes1` (:math:`j\le p`) and, for odd m, :math:`s_m=Z^{\otimes p}`; for
:math:`m=1`, :math:`s_1=Z`. Examples: m = 3 gives the
Pauli matrices X, Y, Z on :math:`\mathbb C^2`; m = 5 acts on :math:`\mathbb C^4`. :math:`\operatorname{Tr}s_i=0` and
:math:`\operatorname{Tr}(s_is_j)=D\delta_{ij}`, so :math:`1,s_1,\dots,s_m` are linearly independent. For
:math:`m\ge2` the :math:`s_i` generate the full matrix algebra :math:`M_D(\mathbb C)` (the first :math:`2p` already
do; for odd m, :math:`s_m` is a multiple of their product).
SpinFactor¶
class SpinFactor(m: int, dtype: torch.dtype=torch.float64, device=None)
The spin factor :math:`V_m=\mathbb R\oplus\mathbb R^m` (:math:`m\ge1`), a Euclidean Jordan algebra of rank 2 and
dimension :math:`m+1` (simple for :math:`m\ge2`; :math:`V_1\cong\mathbb R\oplus\mathbb R`). Elements
:math:`x=(a,u)` are real tensors of shape (..., m+1) (any leading batch shape; everything is batched and
differentiable), with
* product :math:`(a,u)\circ(b,v)=(ab+u\cdot v,\;av+bu)`, unit :math:`(1,0)`;
* spectral decomposition :math:`x=\lambda_+e_++\lambda_-e_-` with :math:`\lambda_\pm=a\pm|u|` and the orthogonal
idempotents :math:`e_\pm=\frac12(1,\pm u/|u|)` (:math:`e_\pm^2=e_\pm`, :math:`e_+\circ e_-=0`,
:math:`e_++e_-=1`); for :math:`u=0` any unit vector replaces :math:`u/|u|` (here the first basis vector);
* functional calculus :math:`f(x)=f(\lambda_+)e_++f(\lambda_-)e_-=\big(\frac{f(\lambda_+)+f(\lambda_-)}2,\;
f[\lambda_+,\lambda_-]\,u\big)` with the divided difference :math:`f[\lambda_+,\lambda_-]=
(f(\lambda_+)-f(\lambda_-))/(2|u|)` (:math:`f'(a)` at :math:`u=0`): smooth in x for smooth f, also at u = 0;
* trace :math:`\operatorname{tr}x=\lambda_++\lambda_-=2a` (the canonical trace), determinant
:math:`\det x=\lambda_+\lambda_-=a^2-|u|^2`, inner product :math:`\operatorname{tr}(x\circ y)=2(ab+u\cdot v)`;
* inverse :math:`x^{-1}=(a,-u)/\det x` for :math:`\det x\ne0`;
* quadratic representation :math:`P(x)y=2x\circ(x\circ y)-x^2\circ y`;
* positive cone (cone of squares) = Lorentz (second-order) cone :math:`\{a\ge|u|\}`, i.e. :math:`\lambda_-\ge0`.
Representation (:meth:`to_matrix`): :math:`\phi(a,u)=a1+\sum_iu_is_i` with the anticommuting symmetries of
:func:`anticommuting_symmetries` on :math:`\mathbb C^D`. Since :math:`s_i\circ s_j=\delta_{ij}1`, :math:`\phi` is an
injective Jordan homomorphism; :math:`(u\cdot s)^2=|u|^21` and :math:`\operatorname{Tr}(u\cdot s)=0` give the
eigenvalues :math:`\lambda_\pm` of :math:`\phi(x)`, each with multiplicity :math:`D/2`; so
:math:`\phi(x)\ge0\iff a\ge|u|` and :math:`\operatorname{Tr}\phi(x)=\frac D2\operatorname{tr}x`. For m = 3,
:math:`\phi` is onto :math:`\operatorname{Herm}_2(\mathbb C)` (:math:`V_3\cong\operatorname{Herm}_2(\mathbb C)`).
SpinFactor.element¶
def element(self, a, u) -> torch.Tensor
(a, u) with a (...) and u (..., m) -> x (..., m+1).
SpinFactor.split¶
def split(self, x: torch.Tensor) -> Tuple[torch.Tensor, torch.Tensor]
x (..., m+1) -> (a (...), u (..., m)).
SpinFactor.unit¶
def unit(self, batch_shape: Tuple[int, ...]=()) -> torch.Tensor
SpinFactor.random¶
def random(self, batch_size: Optional[int]=None, generator: Optional[torch.Generator]=None, cone: bool=False) -> torch.Tensor
Gaussian element(s) (batch_size, m+1) (or (m+1,)); cone=True: a square x², an element of the cone.
SpinFactor.product¶
def product(self, x: torch.Tensor, y: torch.Tensor) -> torch.Tensor
x∘y = (ab + u·v, a v + b u); (..., m+1) x (..., m+1) -> (..., m+1) (leading shapes broadcast).
SpinFactor.square¶
def square(self, x: torch.Tensor) -> torch.Tensor
SpinFactor.power¶
def power(self, x: torch.Tensor, p: int) -> torch.Tensor
Jordan power x^p for an integer p >= 0 (power associative: repeated squaring).
SpinFactor.trace¶
def trace(self, x: torch.Tensor) -> torch.Tensor
tr x = λ+ + λ- = 2a, shape (...).
SpinFactor.determinant¶
def determinant(self, x: torch.Tensor) -> torch.Tensor
det x = λ+ λ- = a² - |u|², shape (...).
SpinFactor.inner¶
def inner(self, x: torch.Tensor, y: torch.Tensor) -> torch.Tensor
tr(x∘y) = 2(ab + u·v), shape (...): the trace inner product (associative: <x∘y, z> = <y, x∘z>).
SpinFactor.inverse¶
def inverse(self, x: torch.Tensor) -> torch.Tensor
x^{-1} = (a, -u) / det x (inf/nan where det x = 0).
SpinFactor.quadratic_representation¶
def quadratic_representation(self, x: torch.Tensor, y: Optional[torch.Tensor]=None)
P(x)y = 2x∘(x∘y) - x²∘y; without y, the map y -> P(x)y.
SpinFactor.triple¶
def triple(self, x: torch.Tensor, y: torch.Tensor, z: torch.Tensor) -> torch.Tensor
Jordan triple product {x,y,z} = (x∘y)∘z + (y∘z)∘x - (x∘z)∘y ({x,y,x} = P(x)y).
SpinFactor.eigenvalues¶
def eigenvalues(self, x: torch.Tensor) -> torch.Tensor
(λ-, λ+) = (a - |u|, a + |u|), ascending, shape (..., 2).
SpinFactor.idempotents¶
def idempotents(self, x: torch.Tensor) -> torch.Tensor
Jordan frame (e-, e+) = (½(1, -w), ½(1, w)), w = u/|u| (first basis vector if u = 0); (..., 2, m+1).
SpinFactor.spectral_decomposition¶
def spectral_decomposition(self, x: torch.Tensor) -> Tuple[torch.Tensor, torch.Tensor]
(eigenvalues (..., 2), idempotents (..., 2, m+1)) with x = Σ λ_i e_i.
SpinFactor.apply_function¶
def apply_function(self, x: torch.Tensor, f: Callable[[torch.Tensor], torch.Tensor]) -> torch.Tensor
f(x) = f(λ+) e+ + f(λ-) e- for an elementwise real f (called once on the eigenvalues (..., 2)), as
((f(λ+) + f(λ-))/2, f[λ+, λ-] u) with the divided difference, which is f' at (nearly) coincident eigenvalues
(|λ+ - λ-| <= eps^(1/3) max|λ±|, the convention of functional.divided_differences). Differentiable in x and
in parameters of f.
SpinFactor.sqrt¶
def sqrt(self, x: torch.Tensor) -> torch.Tensor
The square root in the cone (eigenvalues clamped at 0): the unique y >= 0 with y² = x for x >= 0.
SpinFactor.exp¶
def exp(self, x: torch.Tensor) -> torch.Tensor
SpinFactor.log¶
def log(self, x: torch.Tensor) -> torch.Tensor
log of an element of the interior of the cone (λ- > 0).
SpinFactor.in_cone¶
def in_cone(self, x: torch.Tensor, tol: Optional[float]=None, interior: bool=False) -> torch.Tensor
Membership in the Lorentz cone a >= |u| (interior: a > |u|), up to tol (default 8 eps) relative to |a| + |u|.
Bool tensor of shape (...); no host synchronisation.
SpinFactor.symmetries¶
def symmetries(self, dtype: Optional[torch.dtype]=None) -> torch.Tensor
(m, D, D) anticommuting symmetries of anticommuting_symmetries(m) (complex).
SpinFactor.representation_dim¶
@property
def representation_dim(self) -> int
SpinFactor.to_matrix¶
def to_matrix(self, x: torch.Tensor) -> torch.Tensor
φ(a, u) = a 1 + Σ u_i s_i: (..., m+1) -> Hermitian (..., D, D) (complex).
SpinFactor.from_matrix¶
def from_matrix(self, X: torch.Tensor) -> torch.Tensor
Coefficients of the orthogonal projection onto φ(V_m): a = Tr X / D, u_i = Re Tr(s_i X) / D; (..., m+1).