Merge branch 'develop' into fix_custom_gradient

.codecov.yml

@@ -0,0 +1,9 @@
+coverage:
+  status:
+    project:
+      default:
+        target: 89%
+    patch:
+      default:
+        target: auto
+        threshold: 0%

.github/CHANGELOG.md

@@ -31,6 +31,9 @@
 * Loosen the upper-bound on `thewalrus` and upgrade it.
   [(#454)](https://github.com/XanaduAI/MrMustard/pull/454)
 
+* Update major version of `rich` dependency to version 13.
+  [(#512)](https://github.com/XanaduAI/MrMustard/pull/512)
+
 ### Bug fixes
 * Fix the bug in the order of indices of the triples for DsMap CircuitComponent.
   [(#385)](https://github.com/XanaduAI/MrMustard/pull/385)

Makefile

@@ -1,6 +1,6 @@
 PYTHON3 := $(shell which python3 2>/dev/null)
 TESTRUNNER := -m pytest tests -p no:warnings
-COVERAGE := --cov=mrmustard --cov-report=html:coverage_html_report --cov-append
+COVERAGE := --cov=mrmustard --cov-report=html --cov-append
 
 ifdef check
     CHECK := --check --diff

mrmustard/lab_dev/transformations/__init__.py

@@ -20,11 +20,12 @@
 from .attenuator import *
 from .base import *
 from .bsgate import *
+from .cft import *
 from .dgate import *
 from .fockdamping import *
 from .ggate import *
+from .gaussrandnoise import *
 from .identity import *
 from .rgate import *
 from .s2gate import *
 from .sgate import *
-from .cft import *

mrmustard/lab_dev/transformations/gaussrandnoise.py

@@ -0,0 +1,82 @@
+# Copyright 2024 Xanadu Quantum Technologies Inc.
+
+# Licensed under the Apache License, Version 2.0 (the "License");
+# you may not use this file except in compliance with the License.
+# You may obtain a copy of the License at
+
+#     http://www.apache.org/licenses/LICENSE-2.0
+
+# Unless required by applicable law or agreed to in writing, software
+# distributed under the License is distributed on an "AS IS" BASIS,
+# WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+# See the License for the specific language governing permissions and
+# limitations under the License.
+
+"""
+The class representing a Gaussian random noise channel.
+"""
+
+from __future__ import annotations
+from typing import Sequence
+from mrmustard import math, settings
+from mrmustard.utils.typing import RealMatrix
+from .base import Channel
+from ...physics.representations import Bargmann
+from ...physics import triples
+from ..utils import make_parameter
+
+__all__ = ["GaussRandNoise"]
+
+
+class GaussRandNoise(Channel):
+    r"""
+    The Gaussian random noise channel.
+
+    The number of modes must match half of the size of the Y matrix.
+    .. code-block ::
+
+        >>> import numpy as np
+        >>> from mrmustard.lab_dev import GaussRandNoise
+
+        >>> channel = GaussRandNoise(modes=[1, 2], Y = .2 * np.eye(4))
+        >>> assert channel.modes == [1, 2]
+        >>> assert np.allclose(channel.Y.value, .2 * np.eye(4))
+
+    Args:
+        modes: The modes the channel is applied to
+        Y: The Y matrix of the Gaussian random noise
+        Y_train: whether the Y matrix is a trainable variable
+
+    ..details::
+        The Bargmann representation of the channel is computed via the formulas provided in the paper:
+        https://arxiv.org/pdf/2209.06069
+
+        The channel maps an inout covariance matrix ``cov`` as
+        ..math::
+                cov \mapsto cov + Y.
+    """
+
+    short_name = "GRN"
+
+    def __init__(
+        self,
+        modes: Sequence[int],
+        Y: RealMatrix,
+        Y_trainable: bool = False,
+    ):
+
+        if Y.shape[-1] // 2 != len(modes):
+            raise ValueError(
+                f"The number of modes {len(modes)} does not match the dimension of the "
+                f"Y matrix {Y.shape[-1] // 2}."
+            )
+
+        if (math.real(math.eigvals(Y)) >= -settings.ATOL).min() == 0:
+            raise ValueError("The input Y matrix has negative eigen-values.")
+
+        super().__init__(modes_out=modes, modes_in=modes, name="GRN")
+        self._add_parameter(make_parameter(Y_trainable, value=Y, name="Y", bounds=(None, None)))
+
+        self._representation = Bargmann.from_function(
+            fn=triples.gaussian_random_noise_Abc, Y=self.Y.value
+        )

mrmustard/physics/triples.py

@@ -22,7 +22,7 @@
 import numpy as np
 
 from mrmustard import math, settings
-from mrmustard.utils.typing import Matrix, Vector, Scalar
+from mrmustard.utils.typing import Matrix, Vector, Scalar, RealMatrix
 from mrmustard.physics.gaussian_integrals import complex_gaussian_integral_2
 
 
@@ -605,6 +605,58 @@ def fock_damping_Abc(
     return A, b, c
 
 
+def gaussian_random_noise_Abc(Y: RealMatrix) -> Union[Matrix, Vector, Scalar]:
+    r"""
+    The triple (A, b, c) for the gaussian random noise channel.
+
+    Args:
+        Y: the Y matrix of the Gaussian random noise channel.
+
+    Returns:
+        The ``(A, b, c)`` triple of the Gaussian random noise channel.
+    """
+    m = Y.shape[-1] // 2
+    xi = math.eye(2 * m, dtype=math.complex128) + Y / settings.HBAR
+    xi_inv = math.inv(xi)
+    xi_inv_in_blocks = math.block(
+        [[math.eye(2 * m) - xi_inv, xi_inv], [xi_inv, math.eye(2 * m) - xi_inv]]
+    )
+    R = (
+        1
+        / math.sqrt(complex(2))
+        * math.block(
+            [
+                [
+                    math.eye(m, dtype=math.complex128),
+                    1j * math.eye(m, dtype=math.complex128),
+                    math.zeros((m, 2 * m), dtype=math.complex128),
+                ],
+                [
+                    math.zeros((m, 2 * m), dtype=math.complex128),
+                    math.eye(m, dtype=math.complex128),
+                    -1j * math.eye(m, dtype=math.complex128),
+                ],
+                [
+                    math.eye(m, dtype=math.complex128),
+                    -1j * math.eye(m, dtype=math.complex128),
+                    math.zeros((m, 2 * m), dtype=math.complex128),
+                ],
+                [
+                    math.zeros((m, 2 * m), dtype=math.complex128),
+                    math.eye(m, dtype=math.complex128),
+                    1j * math.eye(m, dtype=math.complex128),
+                ],
+            ]
+        )
+    )
+
+    A = math.Xmat(2 * m) @ R @ xi_inv_in_blocks @ math.conj(R).T
+    b = math.zeros(4 * m)
+    c = 1 / math.sqrt(math.det(xi))
+
+    return A, b, c
+
+
 def bargmann_to_quadrature_Abc(n_modes: int, phi: float) -> tuple[Matrix, Vector, Scalar]:
     r"""
     The ``(A, b, c)`` triple of the multi-mode kernel :math:`\langle \vec{p}|\vec{z} \rangle` between bargmann representation with ABC Ansatz form and quadrature representation with ABC Ansatz.