FEAT: implement BreitWigner expression classes

docs/_extend_docstrings.py

@@ -181,6 +181,15 @@ def extend_BreakupMomentumSquared() -> None:
     _append_latex_doit_definition(expr, deep=True)
 
 
+def extend_BreitWigner() -> None:
+    from ampform.dynamics import BreitWigner
+
+    s, m0, w0, m1, m2, d = sp.symbols("s m0 Gamma0 m1 m2 d", nonnegative=True)
+    L = sp.Symbol("L", integer=True, nonnegative=True)
+    expr = BreitWigner(s, m0, w0, m1, m2, angular_momentum=L, meson_radius=d)
+    _append_latex_doit_definition(expr)
+
+
 def extend_ComplexSqrt() -> None:
     from ampform.sympy.math import ComplexSqrt
 
@@ -358,6 +367,36 @@ def extend_is_within_phasespace() -> None:
     )
 
 
+def extend_MultichannelBreitWigner() -> None:
+    from ampform.dynamics import (
+        ChannelArguments,
+        MultichannelBreitWigner,
+        SimpleBreitWigner,
+    )
+
+    s, m0, w0 = sp.symbols("s m0 Gamma0", nonnegative=True)
+    channels = [
+        ChannelArguments(*sp.symbols("Gamma1 m_a1 m_b1 L1 d", nonnegative=True)),
+        ChannelArguments(*sp.symbols("Gamma2 m_a2 m_b2 L2 d", nonnegative=True)),
+    ]
+    expr = MultichannelBreitWigner(s, m0, channels=channels)
+    _append_latex_doit_definition(expr)
+    bw = SimpleBreitWigner(s, m0, w0)
+    _append_to_docstring(
+        MultichannelBreitWigner,
+        Rf"""
+    with :math:`{sp.latex(bw)}` defined by Equation :eq:`SimpleBreitWigner`, a
+    `SimpleBreitWigner`.
+    """,
+    )
+    _append_to_docstring(
+        ChannelArguments.formulate_width,
+        Rf"""
+    .. math:: {sp.latex(channels[0].formulate_width(s, m0))}
+    """,
+    )
+
+
 def extend_PhaseSpaceFactor() -> None:
     from ampform.dynamics.phasespace import PhaseSpaceFactor
 
@@ -473,6 +512,14 @@ def extend_RotationZMatrix() -> None:
     )
 
 
+def extend_SimpleBreitWigner() -> None:
+    from ampform.dynamics import SimpleBreitWigner
+
+    s, m0, w0 = sp.symbols("s m0 Gamma0", nonnegative=True)
+    expr = SimpleBreitWigner(s, m0, w0)
+    _append_latex_doit_definition(expr)
+
+
 def extend_SphericalHankel1() -> None:
     from ampform.dynamics.form_factor import SphericalHankel1
 

docs/usage/dynamics.ipynb

@@ -409,6 +409,29 @@
     "rel_bw"
    ]
   },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {
+    "editable": true,
+    "jupyter": {
+     "source_hidden": true
+    },
+    "slideshow": {
+     "slide_type": ""
+    },
+    "tags": [
+     "hide-input"
+    ]
+   },
+   "outputs": [],
+   "source": [
+    "from ampform.dynamics import SimpleBreitWigner\n",
+    "\n",
+    "bw = SimpleBreitWigner(s, m0, w0)\n",
+    "Math(aslatex({bw: bw.doit(deep=False)}))"
+   ]
+  },
   {
    "cell_type": "markdown",
    "metadata": {},
@@ -457,16 +480,12 @@
    "cell_type": "code",
    "execution_count": null,
    "metadata": {
-    "jupyter": {
-     "source_hidden": true
-    },
     "tags": []
    },
    "outputs": [],
    "source": [
     "from ampform.dynamics import EnergyDependentWidth\n",
     "\n",
-    "L = sp.Symbol(\"L\", integer=True)\n",
     "width = EnergyDependentWidth(\n",
     "    s=s,\n",
     "    mass0=m0,\n",
@@ -476,8 +495,34 @@
     "    angular_momentum=L,\n",
     "    meson_radius=1,\n",
     "    phsp_factor=PhaseSpaceFactorSWave,\n",
-    ")\n",
-    "Math(aslatex({width: width.evaluate()}))"
+    ")"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {
+    "editable": true,
+    "jupyter": {
+     "source_hidden": true
+    },
+    "slideshow": {
+     "slide_type": ""
+    },
+    "tags": [
+     "hide-input"
+    ]
+   },
+   "outputs": [],
+   "source": [
+    "from ampform.dynamics import BreitWigner\n",
+    "\n",
+    "exprs = [\n",
+    "    BreitWigner(s, m0, w0, m1, m2, angular_momentum=L, meson_radius=d),\n",
+    "    SimpleBreitWigner(s, m0, w0),\n",
+    "    width,\n",
+    "]\n",
+    "Math(aslatex({e: e.doit(deep=False) for e in exprs}))"
    ]
   },
   {

src/ampform/dynamics/__init__.py

@@ -9,6 +9,7 @@
 from warnings import warn
 
 import sympy as sp
+from attrs import asdict, frozen
 
 # pyright: reportUnusedImport=false
 from ampform.dynamics.form_factor import (
@@ -31,9 +32,66 @@
     from sympy.printing.latex import LatexPrinter
 
 
+@unevaluated
+class SimpleBreitWigner(sp.Expr):
+    """Simple, non-relativistic Breit-Wigner with :math:`1` in the nominator."""
+
+    s: Any
+    mass: Any
+    width: Any
+    _latex_repr_ = R"\mathcal{{R}}^\mathrm{{BW}}\left({s}; {mass}, {width}\right)"
+
+    def evaluate(self):
+        s, m0, w0 = self.args
+        return 1 / (m0**2 - s - m0 * w0 * sp.I)
+
+
+@unevaluated
+class BreitWigner(sp.Expr):
+    """Relativistic Breit-Wigner with :math:`1` in the nominator.
+
+    `SimpleBreitWigner` with `EnergyDependentWidth` as width (see Equations
+    :eq:`SimpleBreitWigner` and :eq:`EnergyDependentWidth`).
+    """
+
+    s: Any
+    mass: Any
+    width: Any
+    m1: Any = 0
+    m2: Any = 0
+    angular_momentum: Any = 0
+    meson_radius: Any = 1
+    phsp_factor: PhaseSpaceFactorProtocol = argument(
+        default=PhaseSpaceFactor, sympify=False
+    )
+
+    def evaluate(self):
+        width = self.energy_dependent_width()
+        expr = SimpleBreitWigner(self.s, self.mass, width)
+        if self.angular_momentum == 0 and self.m1 == 0 and self.m2 == 0:
+            return expr.evaluate()
+        return expr
+
+    def energy_dependent_width(self) -> sp.Expr:
+        s, m0, w0, m1, m2, ang_mom, d = self.args
+        if ang_mom == 0 and m1 == 0 and m2 == 0:
+            return w0  # type:ignore[return-value]
+        return EnergyDependentWidth(s, m0, w0, m1, m2, ang_mom, d, self.phsp_factor)
+
+    def _latex_repr_(self, printer: LatexPrinter, *args) -> str:
+        s = printer._print(self.s)
+        function_symbol = R"\mathcal{R}^\mathrm{BW}"
+        mass = printer._print(self.mass)
+        width = printer._print(self.width)
+        arg = Rf"\left({s}; {mass}, {width}\right)"
+        angular_momentum = printer._print(self.angular_momentum)
+        if isinstance(self.angular_momentum, sp.Integer):
+            return Rf"{function_symbol}_{{L={angular_momentum}}}{arg}"
+        return Rf"{function_symbol}_{{{angular_momentum}}}{arg}"
+
+
 @unevaluated
 class EnergyDependentWidth(sp.Expr):
-    # cspell:ignore asner
     r"""Mass-dependent width, coupled to the pole position of the resonance.
 
     See Equation (50.28) in :pdg-review:`2021; Resonances; p.9` and
@@ -75,13 +133,59 @@ def _latex_repr_(self, printer: LatexPrinter, *args) -> str:
         return Rf"{name}\left({s}\right)"
 
 
+@unevaluated
+class MultichannelBreitWigner(sp.Expr):
+    """`BreitWigner` for multiple channels."""
+
+    s: Any
+    mass: Any
+    channels: list[ChannelArguments] = argument(sympify=False)
+
+    def evaluate(self):
+        s = self.s
+        m0 = self.mass
+        width = sum(channel.formulate_width(s, m0) for channel in self.channels)
+        return SimpleBreitWigner(s, m0, width)
+
+    def _latex_repr_(self, printer: LatexPrinter, *args) -> str:
+        latex = R"\mathcal{R}^\mathrm{BW}_\mathrm{multi}\left("
+        latex += printer._print(self.s) + "; "
+        latex += ", ".join(printer._print(channel.width) for channel in self.channels)
+        latex += R"\right)"
+        return latex
+
+
+@frozen
+class ChannelArguments:
+    """Arguments for a channel in a `MultichannelBreitWigner`."""
+
+    width: Any
+    m1: Any = 0
+    m2: Any = 0
+    angular_momentum: Any = 0
+    meson_radius: Any = 1
+
+    def __attrs_post_init__(self) -> None:
+        for name, value in asdict(self).items():
+            object.__setattr__(self, name, sp.sympify(value))
+
+    def formulate_width(self, s: Any, m0: Any) -> sp.Expr:
+        w0 = self.width
+        m1 = self.m1
+        m2 = self.m2
+        angular_momentum = self.angular_momentum
+        meson_radius = self.meson_radius
+        ff2 = FormFactor(s, m1, m2, angular_momentum, meson_radius) ** 2
+        return ff2 * w0 * m0 / sp.sqrt(s)
+
+
 def relativistic_breit_wigner(s, mass0, gamma0) -> sp.Expr:
     """Relativistic Breit-Wigner lineshape.
 
     See :ref:`usage/dynamics:_Without_ form factor` and
     :cite:`ParticleDataGroup:2020ssz`.
     """
-    return gamma0 * mass0 / (mass0**2 - s - gamma0 * mass0 * sp.I)
+    return gamma0 * mass0 * SimpleBreitWigner(s, mass0, gamma0)
 
 
 def relativistic_breit_wigner_with_ff(  # noqa: PLR0917
@@ -99,13 +203,11 @@ def relativistic_breit_wigner_with_ff(  # noqa: PLR0917
     See :ref:`usage/dynamics:_With_ form factor` and :pdg-review:`2021; Resonances;
     p.9`.
     """
-    form_factor = FormFactor(s, m_a, m_b, angular_momentum, meson_radius)
-    energy_dependent_width = EnergyDependentWidth(
+    ff = FormFactor(s, m_a, m_b, angular_momentum, meson_radius)
+    bw = BreitWigner(
         s, mass0, gamma0, m_a, m_b, angular_momentum, meson_radius, phsp_factor
     )
-    return (mass0 * gamma0 * form_factor) / (
-        mass0**2 - s - energy_dependent_width * mass0 * sp.I
-    )
+    return mass0 * gamma0 * ff * bw
 
 
 def formulate_form_factor(s, m_a, m_b, angular_momentum, meson_radius) -> sp.Expr:

tests/dynamics/test_builder.py

@@ -2,7 +2,7 @@
 import sympy as sp
 from qrules.particle import Particle
 
-from ampform.dynamics import EnergyDependentWidth
+from ampform.dynamics import EnergyDependentWidth, SimpleBreitWigner
 from ampform.dynamics.builder import (
     RelativisticBreitWignerBuilder,
     TwoBodyKinematicVariableSet,
@@ -43,7 +43,8 @@ def test_simple_breit_wigner(
         s = variable_set.incoming_state_mass**2
         m0 = sp.Symbol("m_{N}", nonnegative=True)
         w0 = sp.Symbol(R"\Gamma_{N}", nonnegative=True)
-        assert bw == w0 * m0 / (-sp.I * w0 * m0 - s + m0**2)
+        simple_bw = SimpleBreitWigner(s, m0, w0)
+        assert bw == w0 * m0 * simple_bw
         assert set(parameters) == {m0, w0}
         assert parameters[m0] == particle.mass
         assert parameters[w0] == particle.width