From 7927904d4a08178b29e3dfabf6adc23f473cf481 Mon Sep 17 00:00:00 2001
From: Athan Reines <kgryte@gmail.com>
Date: Thu, 1 Dec 2022 02:43:36 -0800
Subject: [PATCH 1/4] Add complex number support to `eigh` and `eigvalsh`
---
spec/API_specification/array_api/linalg.py | 61 +++++++++++++++-------
1 file changed, 41 insertions(+), 20 deletions(-)
diff --git a/spec/API_specification/array_api/linalg.py b/spec/API_specification/array_api/linalg.py
index b3595e1fa..3fc699733 100644
--- a/spec/API_specification/array_api/linalg.py
+++ b/spec/API_specification/array_api/linalg.py
@@ -110,59 +110,80 @@ def diagonal(x: array, /, *, offset: int = 0) -> array:
"""
def eigh(x: array, /) -> Tuple[array]:
- """
- Returns an eigendecomposition x = QLQᵀ of a symmetric matrix (or a stack of symmetric matrices) ``x``, where ``Q`` is an orthogonal matrix (or a stack of matrices) and ``L`` is a vector (or a stack of vectors).
+ r"""
+ Returns an eigenvalue decomposition of a complex Hermitian or real symmetric matrix (or a stack of matrices) ``x``.
+
+ If ``x`` is real-valued, let :math:`\mathbb{K}` be the set of real numbers :math:`\mathbb{R}`, and, if ``x`` is complex-valued, let :math:`\mathbb{K}` be the set of complex numbers :math:`\mathbb{C}`.
+
+ The **eigenvalue decomposition** of a complex Hermitian or real symmetric matrix :math:`x \in\ \mathbb{K}^{n \times n}` is defined as
+
+ .. math::
+ x = Q \Lambda Q^H
+
+ where :math:`Q^H` is the conjugate transpose when :math:`Q` is complex and the transpose when :math:`Q` is real-valued and :math:`\Lambda` is a diagonal matrix whose diagonal elements are the corresponding eigenvalues. When ``x`` is real-valued, :math:`Q` is orthogonal, and, when ``x`` is complex, :math:`Q` is unitary.
.. note::
- The function ``eig`` will be added in a future version of the specification, as it requires complex number support.
+ The eigenvalues of a complex Hermitian or real symmetric matrix are always real.
- ..
- NOTE: once complex numbers are supported, each square matrix must be Hermitian.
+ .. warning:
+ The eigenvectors of a symmetric matrix are not unique and are not continuous with respect to ``x``. Because eigenvectors are not unique, different hardware and software may compute different eigenvectors.
+
+ Non-uniqueness stems from the fact that multiplying an eigenvector by :math:`-1` when ``x`` is real-valued and by :math:`e^{\phi j}` (:math:`\phi \in \mathbb{R}`) when ``x`` is complex produces another set of valid eigenvectors.
.. note::
- Whether an array library explicitly checks whether an input array is a symmetric matrix (or a stack of symmetric matrices) is implementation-defined.
+ Whether an array library explicitly checks whether an input array is Hermitian or a symmetric matrix (or a stack of matrices) is implementation-defined.
+
+ .. note::
+ The function ``eig`` will be added in a future version of the specification.
Parameters
----------
x: array
- input array having shape ``(..., M, M)`` and whose innermost two dimensions form square matrices. Must have a real-valued floating-point data type.
+ input array having shape ``(..., M, M)`` and whose innermost two dimensions form square matrices. Must have a floating-point data type.
Returns
-------
out: Tuple[array]
a namedtuple (``eigenvalues``, ``eigenvectors``) whose
- - first element must have the field name ``eigenvalues`` (corresponding to ``L`` above) and must be an array consisting of computed eigenvalues. The array containing the eigenvalues must have shape ``(..., M)``.
- - second element have have the field name ``eigenvectors`` (corresponding to ``Q`` above) and must be an array where the columns of the inner most matrices contain the computed eigenvectors. These matrices must be orthogonal. The array containing the eigenvectors must have shape ``(..., M, M)``.
-
- Each returned array must have the same real-valued floating-point data type as ``x``.
+ - first element must have the field name ``eigenvalues`` (corresponding to :math:`\operatorname{diag}\Lambda` above) and must be an array consisting of computed eigenvalues. The array containing the eigenvalues must have shape ``(..., M)`` and must have a real-valued floating-point data type whose precision matches the precision of ``x`` (e.g., if ``x`` is ``complex128``, then ``eigenvalues`` must be ``float64``).
+ - second element have have the field name ``eigenvectors`` (corresponding to :math:`Q` above) and must be an array where the columns of the inner most matrices contain the computed eigenvectors. These matrices must be orthogonal. The array containing the eigenvectors must have shape ``(..., M, M)`` and must have the same data type as ``x``.
.. note::
Eigenvalue sort order is left unspecified and is thus implementation-dependent.
"""
def eigvalsh(x: array, /) -> array:
- """
- Returns the eigenvalues of a symmetric matrix (or a stack of symmetric matrices) ``x``.
+ r"""
+ Returns the eigenvalues of a complex Hermitian or real symmetric matrix (or a stack of matrices) ``x``.
- .. note::
- The function ``eigvals`` will be added in a future version of the specification, as it requires complex number support.
+ If ``x`` is real-valued, let :math:`\mathbb{K}` be the set of real numbers :math:`\mathbb{R}`, and, if ``x`` is complex-valued, let :math:`\mathbb{K}` be the set of complex numbers :math:`\mathbb{C}`.
- ..
- NOTE: once complex numbers are supported, each square matrix must be Hermitian.
+ The **eigenvalues** of a complex Hermitian or real symmetric matrix :math:`x \in\ \mathbb{K}^{n \times n}` are defined as the roots (counted with multiplicity) of the polynomial :math:`p` of degree :math:`n` given by
+
+ .. math::
+ p(\lambda) = \operatorname{det}(x - \lambda I_n)
+
+ where :math:`\lambda \in \mathbb{R}` and where :math:`I_n` is the *n*-dimensional identity matrix.
+
+ .. note:;
+ The eigenvalues of a complex Hermitian or real symmetric matrix are always real.
+
+ .. note::
+ Whether an array library explicitly checks whether an input array is Hermitian or a symmetric matrix (or a stack of matrices) is implementation-defined.
.. note::
- Whether an array library explicitly checks whether an input array is a symmetric matrix (or a stack of symmetric matrices) is implementation-defined.
+ The function ``eigvals`` will be added in a future version of the specification.
Parameters
----------
x: array
- input array having shape ``(..., M, M)`` and whose innermost two dimensions form square matrices. Must have a real-valued floating-point data type.
+ input array having shape ``(..., M, M)`` and whose innermost two dimensions form square matrices. Must have a floating-point data type.
Returns
-------
out: array
- an array containing the computed eigenvalues. The returned array must have shape ``(..., M)`` and have the same data type as ``x``.
+ an array containing the computed eigenvalues. The returned array must have shape ``(..., M)`` and have a real-valued floating-point data type whose precision matches the precision of ``x`` (e.g., if ``x`` is ``complex128``, then must have a ``float64`` data type).
.. note::
Eigenvalue sort order is left unspecified and is thus implementation-dependent.
From a74013f7d9c5cf87bd128596a3b410062f5659ab Mon Sep 17 00:00:00 2001
From: Athan Reines <kgryte@gmail.com>
Date: Thu, 1 Dec 2022 02:52:22 -0800
Subject: [PATCH 2/4] Update copy
---
spec/API_specification/array_api/linalg.py | 2 +-
1 file changed, 1 insertion(+), 1 deletion(-)
diff --git a/spec/API_specification/array_api/linalg.py b/spec/API_specification/array_api/linalg.py
index 3fc699733..b13fbb38b 100644
--- a/spec/API_specification/array_api/linalg.py
+++ b/spec/API_specification/array_api/linalg.py
@@ -120,7 +120,7 @@ def eigh(x: array, /) -> Tuple[array]:
.. math::
x = Q \Lambda Q^H
- where :math:`Q^H` is the conjugate transpose when :math:`Q` is complex and the transpose when :math:`Q` is real-valued and :math:`\Lambda` is a diagonal matrix whose diagonal elements are the corresponding eigenvalues. When ``x`` is real-valued, :math:`Q` is orthogonal, and, when ``x`` is complex, :math:`Q` is unitary.
+ with :math:`Q \in \mathbb{K}^{n \times n}` and :math:`\Lambda \in \mathbb{R}^n` and where :math:`Q^H` is the conjugate transpose when :math:`Q` is complex and the transpose when :math:`Q` is real-valued and :math:`\Lambda` is a diagonal matrix whose diagonal elements are the corresponding eigenvalues. When ``x`` is real-valued, :math:`Q` is orthogonal, and, when ``x`` is complex, :math:`Q` is unitary.
.. note::
The eigenvalues of a complex Hermitian or real symmetric matrix are always real.
From 800fe5f2e78ad0484322a735399ce3ca36a6af30 Mon Sep 17 00:00:00 2001
From: Athan Reines <kgryte@gmail.com>
Date: Thu, 1 Dec 2022 02:54:38 -0800
Subject: [PATCH 3/4] Fix spacing
---
spec/API_specification/array_api/linalg.py | 2 ++
1 file changed, 2 insertions(+)
diff --git a/spec/API_specification/array_api/linalg.py b/spec/API_specification/array_api/linalg.py
index b13fbb38b..3cb73b5a9 100644
--- a/spec/API_specification/array_api/linalg.py
+++ b/spec/API_specification/array_api/linalg.py
@@ -149,6 +149,7 @@ def eigh(x: array, /) -> Tuple[array]:
- first element must have the field name ``eigenvalues`` (corresponding to :math:`\operatorname{diag}\Lambda` above) and must be an array consisting of computed eigenvalues. The array containing the eigenvalues must have shape ``(..., M)`` and must have a real-valued floating-point data type whose precision matches the precision of ``x`` (e.g., if ``x`` is ``complex128``, then ``eigenvalues`` must be ``float64``).
- second element have have the field name ``eigenvectors`` (corresponding to :math:`Q` above) and must be an array where the columns of the inner most matrices contain the computed eigenvectors. These matrices must be orthogonal. The array containing the eigenvectors must have shape ``(..., M, M)`` and must have the same data type as ``x``.
+
.. note::
Eigenvalue sort order is left unspecified and is thus implementation-dependent.
"""
@@ -185,6 +186,7 @@ def eigvalsh(x: array, /) -> array:
out: array
an array containing the computed eigenvalues. The returned array must have shape ``(..., M)`` and have a real-valued floating-point data type whose precision matches the precision of ``x`` (e.g., if ``x`` is ``complex128``, then must have a ``float64`` data type).
+
.. note::
Eigenvalue sort order is left unspecified and is thus implementation-dependent.
"""
From a1c9e7c963f169f92ebd3222d51d11dd03f74e5e Mon Sep 17 00:00:00 2001
From: Athan Reines <kgryte@gmail.com>
Date: Thu, 1 Dec 2022 02:55:45 -0800
Subject: [PATCH 4/4] Fix directive
---
spec/API_specification/array_api/linalg.py | 2 +-
1 file changed, 1 insertion(+), 1 deletion(-)
diff --git a/spec/API_specification/array_api/linalg.py b/spec/API_specification/array_api/linalg.py
index 3cb73b5a9..f93cdeca8 100644
--- a/spec/API_specification/array_api/linalg.py
+++ b/spec/API_specification/array_api/linalg.py
@@ -125,7 +125,7 @@ def eigh(x: array, /) -> Tuple[array]:
.. note::
The eigenvalues of a complex Hermitian or real symmetric matrix are always real.
- .. warning:
+ .. warning::
The eigenvectors of a symmetric matrix are not unique and are not continuous with respect to ``x``. Because eigenvectors are not unique, different hardware and software may compute different eigenvectors.
Non-uniqueness stems from the fact that multiplying an eigenvector by :math:`-1` when ``x`` is real-valued and by :math:`e^{\phi j}` (:math:`\phi \in \mathbb{R}`) when ``x`` is complex produces another set of valid eigenvectors.