diff --git a/docs/_extend_docstrings.py b/docs/_extend_docstrings.py index 60fe5cb6..d3ce68bd 100644 --- a/docs/_extend_docstrings.py +++ b/docs/_extend_docstrings.py @@ -670,9 +670,9 @@ def extend_relativistic_breit_wigner_with_ff() -> None: .. math:: {sp.latex(rel_bw_with_ff)} :label: relativistic_breit_wigner_with_ff - where :math:`\Gamma(s)` is defined by :eq:`EnergyDependentWidth`, :math:`B_L^2` is - defined by :eq:`BlattWeisskopfSquared`, and :math:`q^2` is defined by - :eq:`BreakupMomentumSquared`. + where :math:`\Gamma(s)` is defined by :eq:`EnergyDependentWidth`, + :math:`\hat{{B}}_L^2` is defined by :eq:`BlattWeisskopfSquared`, and :math:`q^2` is + defined by :eq:`BreakupMomentumSquared`. """, ) diff --git a/docs/analyticity/integration-algorithms.ipynb b/docs/analyticity/integration-algorithms.ipynb index 083ef2eb..312eb38d 100644 --- a/docs/analyticity/integration-algorithms.ipynb +++ b/docs/analyticity/integration-algorithms.ipynb @@ -147,7 +147,7 @@ " \\frac\n", " {\\rho\\!\\left(s'\\right) n_\\ell^2\\!\\left(s'\\right) ds'}\n", " {\\left(s' - s_\\mathrm{thr}\\right) \\left(s'- s\\right)} \\\\\n", - "n_\\ell^2(s') \\;&=\\; \\mathcal{F}_\\ell^2\\!\\left(s', m_1, m_2\\right) \\\\\n", + "n_\\ell^2(s') \\;&=\\; \\hat{\\mathcal{F}}_\\ell^2\\!\\left(s', m_1, m_2\\right) \\\\\n", "s_\\mathrm{thr} \\;&=\\; (m_1 + m_2)^2\n", "\\end{aligned}\n", "$$\n", diff --git a/docs/conf.py b/docs/conf.py index f0143cee..17713690 100644 --- a/docs/conf.py +++ b/docs/conf.py @@ -192,7 +192,6 @@ def _get_dataclasses(module): "sphinx_comments", "sphinx_copybutton", "sphinx_design", - "sphinx_hep_pdgref", "sphinx_pybtex_etal_style", "sphinx_thebe", "sphinx_togglebutton", diff --git a/docs/dynamics.ipynb b/docs/dynamics.ipynb index f5fba495..659d7a29 100644 --- a/docs/dynamics.ipynb +++ b/docs/dynamics.ipynb @@ -80,7 +80,7 @@ "cell_type": "markdown", "metadata": {}, "source": [ - "AmpForm uses Blatt–Weisskopf functions $B_L$ as _barrier factors_ (also called _form factors_, see {class}`.BlattWeisskopfSquared` and **[TR-029](https://compwa.github.io/report/029)**):" + "AmpForm uses Blatt–Weisskopf functions $\\hat{B}_L$ as _barrier factors_ (also called _form factors_, see {class}`.BlattWeisskopfSquared` and **[TR-029](https://compwa.github.io/report/029)**). The hat indicates that they are normalized to $\\hat{B}_L(1)=1$:" ] }, { @@ -151,7 +151,7 @@ "cell_type": "markdown", "metadata": {}, "source": [ - "The Blatt–Weisskopf form factor is used to 'dampen' the breakup-momentum at threshold and when going to infinity. A usual choice for $z$ is therefore $z=q^2d^2$ with $q^2$ the {class}`.BreakupMomentumSquared` and $d$ the impact parameter (also called meson radius). The {class}`.FormFactor` expression class can be used for this:" + "The Blatt–Weisskopf form factor is used to 'dampen' the breakup-momentum at threshold and when going to infinity. A usual choice for $z$ is therefore $z=q^2d^2$ with $q^2$ the {class}`.BreakupMomentumSquared` and $d$ the impact parameter (also called meson radius). The {class}`.FormFactor` expression class can be used for this. Set `normalize=False` to omit the normalization constant $\\left|h_L^{(1)}(1)\\right|$, which gives the non-normalized vertex factor $n_L$ of Equation (50.33) in the [PDG review on resonances](https://pdg.lbl.gov/2026/reviews/rpp2026-rev-resonances.pdf#page=12):" ] }, { @@ -692,7 +692,7 @@ "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", - "version": "3.13.13" + "version": "3.13.15" } }, "nbformat": 4, diff --git a/docs/dynamics/k-matrix.ipynb b/docs/dynamics/k-matrix.ipynb index cc4da6b0..0f7f1f37 100644 --- a/docs/dynamics/k-matrix.ipynb +++ b/docs/dynamics/k-matrix.ipynb @@ -16,7 +16,7 @@ "cell_type": "markdown", "metadata": {}, "source": [ - "While {mod}`ampform` does not yet provide a generic way to formulate an amplitude model with $\\boldsymbol{K}$-matrix dynamics, the (experimental) {mod}`.kmatrix` module makes it fairly simple to produce a symbolic expression for a parameterized $\\boldsymbol{K}$-matrix with an arbitrary number of poles and channels and play around with it interactively. For more info on the $\\boldsymbol{K}$-matrix, see the classic paper by Chung {cite}`Chung:1995dx`, {pdg-review}`2021; Resonances`, or this instructive presentation {cite}`Meyer:2008-MatrixTutorial`.\n", + "While {mod}`ampform` does not yet provide a generic way to formulate an amplitude model with $\\boldsymbol{K}$-matrix dynamics, the (experimental) {mod}`.kmatrix` module makes it fairly simple to produce a symbolic expression for a parameterized $\\boldsymbol{K}$-matrix with an arbitrary number of poles and channels and play around with it interactively. For more info on the $\\boldsymbol{K}$-matrix, see the classic paper by Chung {cite}`Chung:1995dx`, [PDG2026, §Resonances, p.14](https://pdg.lbl.gov/2026/reviews/rpp2026-rev-resonances.pdf#page=14), or this instructive presentation {cite}`Meyer:2008-MatrixTutorial`.\n", "\n", "Section {ref}`dynamics/k-matrix:Physics` summarizes {cite}`Chung:1995dx`, so that the {mod}`.kmatrix` module can reference to the equations. It also points out some subtleties and deviations.\n", "\n", @@ -415,7 +415,11 @@ "\n", "with $\\gamma_{R,i}$ some _real_ constants and $\\Gamma^0_{R,i}$ the **partial width** of each pole. In the Lorentz-invariant form, the fixed width $\\Gamma^0$ is replaced by an \"energy dependent\" {class}`.EnergyDependentWidth` $\\Gamma(s)$.[^phase-space-factor-normalization] The **width** for each pole can be computed as $\\Gamma^0_R = \\sum_i\\Gamma^0_{R,i}$.\n", "\n", - "[^phase-space-factor-normalization]: Unlike Eq. (77) in {cite}`Chung:1995dx`, AmpForm defines {class}`.EnergyDependentWidth` as in {pdg-review}`2021; Resonances; p.6`, Eq. (50.28). The difference is that the phase space factor denoted by $\\rho_i$ in Eq. (77) in {cite}`Chung:1995dx` is divided by the phase space factor at the pole position $m_R$. So in AmpForm, the choice is $\\rho_i \\to \\frac{\\rho_i(s)}{\\rho_i(m_R)}$." + ":::{warning}\n", + "Eq. (50.28) no longer appears in [PDG2026, §Resonances, p.12](https://pdg.lbl.gov/2026/reviews/rpp2026-rev-resonances.pdf#page=12). The width is now defined in terms of bare couplings, $\\Gamma_b(s) = g_b^2\\rho_b(s)n_b^2(s)/m_\\mathrm{BW}$ (Eq. (50.32)), and Eq. (50.35) trades $g_b$ for the partial width $\\Gamma_{\\mathrm{BW},b}$. Combining the two gives the old Eq. (50.28), but the PDG stresses that this substitution is only valid for narrow resonances with all channel thresholds below $m_\\mathrm{BW}$.\n", + ":::\n", + "\n", + "[^phase-space-factor-normalization]: Unlike Eq. (77) in {cite}`Chung:1995dx`, AmpForm defines {class}`.EnergyDependentWidth` as in [PDG2021, Eq. (50.28)](https://pdg.lbl.gov/2021/reviews/rpp2021-rev-resonances.pdf#page=9). The difference is that the phase space factor denoted by $\\rho_i$ in Eq. (77) in {cite}`Chung:1995dx` is divided by the phase space factor at the pole position $m_R$. So in AmpForm, the choice is $\\rho_i \\to \\frac{\\rho_i(s)}{\\rho_i(m_R)}$." ] }, { @@ -441,7 +445,7 @@ "\n", "with $B_{R,i}(q(s))$ the **centrifugal damping factor** (see {class}`.FormFactor` and {class}`.BlattWeisskopfSquared`) for channel $i$ and $\\beta_R^0$ some (generally complex) constants that describe the production information of the decaying state $R$. Usually, these constants are rescaled just like the residue functions in {eq}`residue-function`:\n", "\n", - "[^damping-factor-P-parametrization]: Just as with [^phase-space-factor-normalization], we have smuggled a bit in the last equation in order to be able to reproduce Equation (50.23) in {pdg-review}`2021; Resonances; p.9` in the case $n=1,n_R=1$, on which {func}`.relativistic_breit_wigner_with_ff` is based.\n", + "[^damping-factor-P-parametrization]: Just as with [^phase-space-factor-normalization], we have smuggled a bit in the last equation in order to be able to reproduce [PDG2026, Eq. (50.37)](https://pdg.lbl.gov/2026/reviews/rpp2026-rev-resonances.pdf#page=13) in the case $n=1,n_R=1$, on which {func}`.relativistic_breit_wigner_with_ff` is based.\n", "\n", "```{margin}\n", "{cite}`Chung:1995dx` Eq. (121)\n", @@ -1332,7 +1336,7 @@ "cell_type": "markdown", "metadata": {}, "source": [ - "[^pole-vs-resonance]: See {pdg-review}`2021; Resonances`, Section 50.1, for a discussion about what poles and resonances are. See also the intro to Section 5 in {cite}`Chung:1995dx`." + "[^pole-vs-resonance]: See [PDG2026, §50.1](https://pdg.lbl.gov/2026/reviews/rpp2026-rev-resonances.pdf#page=1) for a discussion about what poles and resonances are. See also the intro to Section 5 in {cite}`Chung:1995dx`." ] }, { diff --git a/pyproject.toml b/pyproject.toml index 4a2c1bc5..f3b5a058 100644 --- a/pyproject.toml +++ b/pyproject.toml @@ -71,7 +71,6 @@ doc = [ "sphinx-comments", "sphinx-copybutton", "sphinx-design", - "sphinx-hep-pdgref", "sphinx-pybtex-etal-style", "sphinx-thebe", "sphinx-togglebutton", diff --git a/src/ampform/dynamics/__init__.py b/src/ampform/dynamics/__init__.py index 7b45c893..60bba6d8 100644 --- a/src/ampform/dynamics/__init__.py +++ b/src/ampform/dynamics/__init__.py @@ -5,7 +5,7 @@ from __future__ import annotations -from typing import TYPE_CHECKING, Any +from typing import TYPE_CHECKING, Any, Literal from warnings import warn import sympy as sp @@ -35,24 +35,47 @@ @unevaluated class SimpleBreitWigner(sp.Expr): - r"""Simple Breit–Wigner with :math:`m_0 \Gamma_0` in the numerator.""" + r"""Simple Breit–Wigner with a configurable numerator. + + With the default ``numerator="mass-width"``, the propagator is multiplied by + :math:`m_0 \Gamma_0`, so that + :math:`\left|\hat{\mathcal{R}}^\mathrm{BW}(m_0^2)\right| = 1`. Set + ``numerator="unity"`` for the dressed propagator of `PDG2026, Eq. (50.31) + `__, which is + also the convention of `MultichannelBreitWigner`. + """ s: Any mass: Any width: Any - _latex_repr_ = R"\mathcal{{R}}^\mathrm{{BW}}\left({s}; {mass}, {width}\right)" + numerator: Literal["mass-width", "unity"] = argument( + default="mass-width", kw_only=True, sympify=False + ) def evaluate(self): s, m0, w0 = self.args - return m0 * w0 * _formulate_breit_wigner(s, m0, w0) + numerator = _formulate_numerator(self.numerator, m0, w0) + return numerator * _formulate_breit_wigner(s, m0, w0) + + def _latex_repr_(self, printer: LatexPrinter, *args) -> str: + s, mass, width = map(printer._print, self.args) + function_symbol = _get_breit_wigner_symbol(self.numerator) + return Rf"{function_symbol}\left({s}; {mass}, {width}\right)" @unevaluated class BreitWigner(sp.Expr): - r"""Relativistic Breit–Wigner with :math:`m_0 \Gamma_0` in the numerator. - - Uses an `EnergyDependentWidth` in the denominator (see Equations - :eq:`BreitWigner` and :eq:`EnergyDependentWidth`). + r"""Relativistic Breit–Wigner with a configurable numerator. + + Uses an `EnergyDependentWidth` in the denominator (see Equations :eq:`BreitWigner` + and :eq:`EnergyDependentWidth`). With the default ``numerator="mass-width"``, the + propagator is multiplied by :math:`m_0 \Gamma_0`, so that + :math:`\left|\hat{\mathcal{R}}^\mathrm{BW}(m_0^2)\right| = 1`, because + :math:`\Gamma(m_0^2) = \Gamma_0`. Set ``numerator="unity"`` for the dressed + propagator of `PDG2026, Eq. (50.31) + `__, which is + also the convention of `MultichannelBreitWigner`. The numerator does not affect the + `.FormFactor` inside the `EnergyDependentWidth`, where its normalization cancels. """ s: Any @@ -65,12 +88,14 @@ class BreitWigner(sp.Expr): phsp_factor: PhaseSpaceFactorProtocol = argument( default=PhaseSpaceFactor, sympify=False ) # ty: ignore[invalid-assignment] + numerator: Literal["mass-width", "unity"] = argument( + default="mass-width", kw_only=True, sympify=False + ) def evaluate(self): width = self.energy_dependent_width() - return ( - self.mass * self.width * _formulate_breit_wigner(self.s, self.mass, width) - ) + numerator = _formulate_numerator(self.numerator, self.mass, self.width) + return numerator * _formulate_breit_wigner(self.s, self.mass, width) def energy_dependent_width(self) -> EnergyDependentWidth | sp.Basic: s, m0, w0, m1, m2, ang_mom, d = self.args @@ -80,7 +105,7 @@ def energy_dependent_width(self) -> EnergyDependentWidth | sp.Basic: def _latex_repr_(self, printer: LatexPrinter, *args) -> str: s = printer._print(self.s) - function_symbol = R"\mathcal{R}^\mathrm{BW}" + function_symbol = _get_breit_wigner_symbol(self.numerator) mass = printer._print(self.mass) width = printer._print(self.width) arg = Rf"\left({s}; {mass}, {width}\right)" @@ -94,10 +119,20 @@ def _latex_repr_(self, printer: LatexPrinter, *args) -> str: class EnergyDependentWidth(sp.Expr): r"""Mass-dependent width, coupled to the pole position of the resonance. - See Equation (50.28) in :pdg-review:`2021; Resonances; p.9` and + See `PDG2021, Eq. (50.28) + `__ and :cite:`ParticleDataGroup:2012pjm`, equation (6). Default value for :code:`phsp_factor` is `.PhaseSpaceFactor`. + .. warning:: Equation (50.28) no longer appears in + `PDG2026, §Resonances, p.12 `__. + The width is now defined in terms of bare couplings, :math:`\Gamma_b(s) = + g_b^2 \rho_b(s) n_b^2(s) / m_\mathrm{BW}` (Equation (50.32)), and Equation (50.35) + trades :math:`g_b` for the partial width :math:`\Gamma_{\mathrm{BW},b}`. + Combining the two gives the old Equation (50.28), but the PDG stresses that this + substitution is only valid for narrow resonances with all channel thresholds + below :math:`m_\mathrm{BW}`. + Note that the `.FormFactor` of AmpForm is normalized in the sense that equal powers of :math:`z` appear in the nominator and the denominator, while the definition in the PDG (as well as some other sources), always have :math:`1` in the nominator of @@ -144,7 +179,9 @@ class MultichannelBreitWigner(sp.Expr): where :math:`g_i^2` is the coupling squared, :math:`\rho_i` is a `.PhaseSpaceFactor`, and :math:`F_{L_i}` is a `.FormFactor`. Unlike an - `EnergyDependentWidth`, a channel term is not normalized at the pole position. + `EnergyDependentWidth`, a channel term is not normalized at the pole position. See + `PDG2026, Eqs. (50.31) and (50.32) + `__. """ s: Any @@ -174,6 +211,9 @@ class ChannelArguments(sp.Expr): .. math:: \Gamma_i^\text{ch}(s) = \frac{g_i^2}{m_0} \rho_i(s) F_{L_i}^2(s) + + See `PDG2026, Eq. (50.32) + `__. """ s: Any @@ -196,6 +236,27 @@ def _formulate_breit_wigner(s: Any, mass: Any, width: Any) -> sp.Expr: return 1 / (mass**2 - s - sp.I * mass * width) +def _formulate_numerator(numerator: str, mass: Any, width: Any) -> Any: + _check_numerator(numerator) + if numerator == "mass-width": + return mass * width + return sp.S.One + + +def _get_breit_wigner_symbol(numerator: str) -> str: + _check_numerator(numerator) + if numerator == "mass-width": + return R"\hat{\mathcal{R}}^\mathrm{BW}" + return R"\mathcal{R}^\mathrm{BW}" + + +def _check_numerator(numerator: str) -> None: + allowed = ("mass-width", "unity") + if numerator not in allowed: + msg = f"Unknown numerator {numerator!r}, expected one of {', '.join(map(repr, allowed))}" + raise ValueError(msg) + + def relativistic_breit_wigner(s, mass0, gamma0) -> sp.Expr: """Relativistic Breit–Wigner lineshape. @@ -220,7 +281,8 @@ def relativistic_breit_wigner_with_ff( # ruff: ignore[too-many-positional-argum ) -> sp.Expr: """Relativistic Breit–Wigner with `.FormFactor`. - See :ref:`dynamics:_With_ form factor` and :pdg-review:`2021; Resonances; p.9`. + See :ref:`dynamics:_With_ form factor` and `PDG2026, §Resonances, p.12 + `__. """ ff = FormFactor(s, m_a, m_b, angular_momentum, meson_radius) bw = BreitWigner( diff --git a/src/ampform/dynamics/builder.py b/src/ampform/dynamics/builder.py index 4b08568a..d3d40e22 100644 --- a/src/ampform/dynamics/builder.py +++ b/src/ampform/dynamics/builder.py @@ -101,8 +101,9 @@ class RelativisticBreitWignerBuilder: Args: form_factor: Formulate a relativistic Breit–Wigner function multiplied - by a Blatt–Weisskopf form factor (`.FormFactor`), like in Equation (50.26) - on :pdg-review:`2021; Resonances; p.9`. + by a Blatt–Weisskopf form factor (`.FormFactor`), like in `PDG2026, Eqs. + (50.33) and (50.37) + `__. energy_dependent_width: Use an `.EnergyDependentWidth` in the denominator of the Breit–Wigner. phsp_factor: A class that complies with the diff --git a/src/ampform/dynamics/form_factor.py b/src/ampform/dynamics/form_factor.py index d2c3aa8a..e3789d72 100644 --- a/src/ampform/dynamics/form_factor.py +++ b/src/ampform/dynamics/form_factor.py @@ -8,19 +8,28 @@ import sympy as sp from ampform.kinematics.phasespace import BreakupMomentumSquared -from ampform.sympy import unevaluated +from ampform.sympy import argument, unevaluated if TYPE_CHECKING: from collections.abc import Callable + from sympy.printing.latex import LatexPrinter + @unevaluated class FormFactor(sp.Expr): - """Formulate a Blatt–Weisskopf form factor. - - Returns the production process factor :math:`n_a` from Equation (50.26) in - :pdg-review:`2021; Resonances; p.9`, which features the - `~sympy.functions.elementary.miscellaneous.sqrt` of a `.BlattWeisskopfSquared`. + r"""Formulate a Blatt–Weisskopf form factor. + + Returns the `~sympy.functions.elementary.miscellaneous.sqrt` of a + `.BlattWeisskopfSquared` with :math:`z = q^2 d^2`, where :math:`q^2` is the + `.BreakupMomentumSquared` and :math:`d` is the meson radius. With + ``normalize=False``, this is the production process factor :math:`n_a` from + `PDG2026, Eq. (50.33) + `__, with + :math:`d = 1/q_0`. The default normalized form factor :math:`\hat{\mathcal{F}}_L` + differs from :math:`n_a` by the constant :math:`\left|h_L^{(1)}(1)\right|`. This + constant cancels in ratios, such as in `.EnergyDependentWidth`, but not when the + form factor is used as a vertex factor. """ s: Any @@ -28,32 +37,50 @@ class FormFactor(sp.Expr): m2: Any angular_momentum: Any meson_radius: Any = 1 - - _latex_repr_ = R"\mathcal{{F}}_{{{angular_momentum}}}\left({s}, {m1}, {m2}\right)" + normalize: bool = argument(default=True, kw_only=True, sympify=False) def evaluate(self): s, m1, m2, angular_momentum, meson_radius = self.args q2 = BreakupMomentumSquared(s, m1, m2) - ff_squared = BlattWeisskopfSquared(q2 * meson_radius**2, angular_momentum) + ff_squared = BlattWeisskopfSquared( + q2 * meson_radius**2, angular_momentum, normalize=self.normalize + ) return sp.sqrt(ff_squared) + def _latex_repr_(self, printer: LatexPrinter, *args) -> str: + s, m1, m2, angular_momentum = map(printer._print, self.args[:4]) + symbol = R"\hat{\mathcal{F}}" if self.normalize else R"\mathcal{F}" + return Rf"{symbol}_{{{angular_momentum}}}\left({s}, {m1}, {m2}\right)" + @unevaluated class BlattWeisskopfSquared(sp.Expr): - r"""Normalized Blatt–Weisskopf function :math:`B_L^2(z)`, with :math:`B_L^2(1)=1`. + r"""Normalized Blatt–Weisskopf function :math:`\hat{B}_L^2(z)`. Args: - z: Argument of the Blatt–Weisskopf function :math:`B_L^2(z)`. A usual - choice is :math:`z = (d q)^2` with :math:`d` the impact parameter and - :math:`q` the breakup-momentum (see `.BreakupMomentumSquared`). + z: Argument of the Blatt–Weisskopf function. A usual choice is :math:`z = (d + q)^2` with :math:`d` the impact parameter and :math:`q` the breakup-momentum + (see `.BreakupMomentumSquared`). angular_momentum: Angular momentum :math:`L` of the decaying particle. - Note that equal powers of :math:`z` appear in the nominator and the denominator, - while some sources define an *non-normalized* form factor :math:`F_L` with :math:`1` - in the nominator, instead of :math:`z^L`. See for instance Equation (50.27) in - :pdg-review:`2021; Resonances; p.9`. We normalize the form factor such that - :math:`B_L^2(1)=1` and that :math:`B_L^2` is unitless no matter what :math:`z` is. + normalize: Set to `False` to omit the normalization constant + :math:`\left|h_L^{(1)}(1)\right|^2`. The resulting :math:`B_L^2(z)` equals + :math:`z^L F_L^2(\sqrt{z})`, where :math:`F_L` is the non-normalized + Blatt–Weisskopf function of `PDG2026, Eq. (50.34) + `__. + This is the square of the factor :math:`n_L` of Equation (50.33). + + The hat indicates the normalization :math:`\hat{B}_L^2(1)=1`. Both variants have + equal powers of :math:`z` in the numerator and the denominator, so they are unitless + no matter what :math:`z` is. The PDG function :math:`F_L` instead has :math:`1` in + the numerator and carries the threshold factor :math:`z^L` separately. + + >>> z = sp.Symbol("z", nonnegative=True) + >>> BlattWeisskopfSquared(z, angular_momentum=2).doit() + 13*z**2/(z**2 + 3*z + 9) + >>> BlattWeisskopfSquared(z, angular_momentum=2, normalize=False).doit() + z**2/(z**2 + 3*z + 9) .. seealso:: :ref:`dynamics:Form factor`, :doc:`TR-029`, and :cite:`Chung:2015-FormulasAngularMomentumBarrier`. @@ -63,34 +90,40 @@ class BlattWeisskopfSquared(sp.Expr): z: Any angular_momentum: Any - _latex_repr_ = R"B_{{{angular_momentum}}}^2\left({z}\right)" + normalize: bool = argument(default=True, kw_only=True, sympify=False) def evaluate(self) -> sp.Expr: z, ell = self.args if ell.free_symbols: - return _formulate_blatt_weisskopf(ell, z) - expr = _get_polynomial_blatt_weisskopf(ell)(z) + return _formulate_blatt_weisskopf(ell, z, self.normalize) + expr = _get_polynomial_blatt_weisskopf(ell, self.normalize)(z) return sp.sympify(expr) + def _latex_repr_(self, printer: LatexPrinter, *args) -> str: + z, angular_momentum = map(printer._print, self.args) + symbol = R"\hat{B}" if self.normalize else "B" + return Rf"{symbol}_{{{angular_momentum}}}^2\left({z}\right)" + -@lru_cache(maxsize=20) -def _get_polynomial_blatt_weisskopf(ell: int | sp.Integer) -> Callable[[Any], Any]: +@lru_cache(maxsize=40) +def _get_polynomial_blatt_weisskopf( + ell: int | sp.Integer, normalize: bool = True +) -> Callable[[Any], Any]: """Get the Blatt–Weisskopf factor as a fraction of polynomials. See https://github.com/ComPWA/ampform/issues/426. """ z = sp.Symbol("z", nonnegative=True, real=True) - expr = _formulate_blatt_weisskopf(ell, z) + expr = _formulate_blatt_weisskopf(ell, z, normalize) expr = expr.doit().simplify() return sp.lambdify(z, expr, "math") -def _formulate_blatt_weisskopf(ell, z) -> sp.Expr: - return ( - sp.Abs(SphericalHankel1(ell, 1)) ** 2 - / sp.Abs(SphericalHankel1(ell, sp.sqrt(z))) ** 2 - / z - ) +def _formulate_blatt_weisskopf(ell, z, normalize: bool = True) -> sp.Expr: + expr = 1 / sp.Abs(SphericalHankel1(ell, sp.sqrt(z))) ** 2 / z + if normalize: + return sp.Abs(SphericalHankel1(ell, z=1)) ** 2 * expr + return expr @unevaluated diff --git a/src/ampform/dynamics/phasespace.py b/src/ampform/dynamics/phasespace.py index c4343c42..f650d2e2 100644 --- a/src/ampform/dynamics/phasespace.py +++ b/src/ampform/dynamics/phasespace.py @@ -1,9 +1,11 @@ """Different parametrizations of phase space factors. Phase space factors are computed by integrating over the phase space element given by -Equation (49.12) in :pdg-review:`2021; Kinematics; p.2`. See also Equation (50.9) on -:pdg-review:`2021; Resonances; p.6`. This integral is not always easy to solve, which -leads to different parametrizations. +`PDG2026, Eq. (49.12) +`__. See also +`PDG2026, Eq. (50.11) +`__. This integral +is not always easy to solve, which leads to different parametrizations. This module provides several parametrizations. They all comply with the `PhaseSpaceFactorProtocol`, so that they can be used in parametrizations like @@ -58,8 +60,10 @@ def __call__(self, s, m1, m2) -> sp.Expr: class PhaseSpaceFactor(sp.Expr): r"""Standard phase-space factor, using a definition consistent with `.BreakupMomentum`. - See :pdg-review:`2025; Resonances; p.6`, Equation (50.11). We ignore the factor - :math:`\frac{1}{16\pi}` as done in :cite:`Chung:1995-PrimerKmatrixFormalism`, p.5. + See `PDG2026, Eq. (50.11) + `__. We ignore + the factor :math:`\frac{1}{16\pi}` as done in + :cite:`Chung:1995-PrimerKmatrixFormalism`, p.5. Similarly to `.BreakupMomentum`, this class represents the numerator as a single square root for better numerical performance. This comes at the cost of a :ref:`more @@ -187,10 +191,14 @@ def _latex_repr_(self, printer: LatexPrinter, *args) -> str: @unevaluated class PhaseSpaceFactorSWave(sp.Expr): - """Phase space factor using :func:`ChewMandelstamSWave`. + r"""Phase space factor using :func:`ChewMandelstamSWave`. This `PhaseSpaceFactor` provides an analytic continuation for decay products with - both equal and unequal masses (compare `EqualMassPhaseSpaceFactor`). + both equal and unequal masses (compare `EqualMassPhaseSpaceFactor`). Following + `PDG2026, §50.3.3 + `__, the + Chew–Mandelstam function :math:`\Sigma(s)` replaces :math:`i\rho(s)`, so this class + returns :math:`-i\Sigma(s)`. """ s: Any @@ -211,10 +219,22 @@ def _latex_repr_(self, printer: LatexPrinter, *args) -> str: @unevaluated class ChewMandelstamSWave(sp.Expr): - """Chew–Mandelstam class for :math:`S`-waves (no angular momentum). + r"""Chew–Mandelstam class for :math:`S`-waves (no angular momentum). + + See `PDG2021, Eq. (50.40) + `__. As in + `PhaseSpaceFactor`, we ignore the factor :math:`\frac{1}{16\pi}`. As a trick, the square root in :math:`q` is defined with `.ComplexSqrt` so that this function has a well-defined behavior along the negative real axis. + + .. warning:: This function is given as `PDG2026, Eq. (50.46) + `__, which + contains two apparent typos: the denominator inside the first logarithm reads + :math:`2m_1m_1` instead of :math:`2m_1m_2`, and the last + term contains :math:`1/s_\mathrm{thr}^2` instead of + :math:`1/s_\mathrm{thr}=1/(m_1+m_2)^2`. This implementation follows the 2021 + version. """ s: Any @@ -308,6 +328,10 @@ def _latex_repr_(self, printer: LatexPrinter, *args) -> str: class ChewMandelstamIntegral(sp.Expr): """Dispersion integral for obtaining the analytic phase space factor for angular momenta L>0. + See `PDG2026, Eq. (50.45) + `__. The + integral is subtracted at the channel threshold. + Parameters: s: Mandelstam variable s. m1: Mass of particle 1. @@ -366,10 +390,18 @@ def _latex_repr_(self, printer: LatexPrinter, *args) -> str: class EqualMassPhaseSpaceFactor(sp.Expr): """Analytic continuation for the `PhaseSpaceFactor`. - See :pdg-review:`2018; Resonances; p.9` and :doc:`/analyticity/phasespace-factors`. + See `PDG2018, §Resonances, p.9 + `__ and + :doc:`/analyticity/phasespace-factors`. **Warning**: The PDG specifically derives this formula for a two-body decay *with equal masses*. + + .. warning:: This formula no longer appears in + `PDG2026, §Resonances, p.16 `__. + The PDG now gives the :math:`S`-wave Chew–Mandelstam function for arbitrary + masses instead (Equation (50.46)), which is subtracted at the channel threshold. + See `.ChewMandelstamSWave`. """ s: Any diff --git a/src/ampform/kinematics/phasespace.py b/src/ampform/kinematics/phasespace.py index eb54b9ab..5b94563b 100644 --- a/src/ampform/kinematics/phasespace.py +++ b/src/ampform/kinematics/phasespace.py @@ -22,7 +22,8 @@ class BreakupMomentum(sp.Expr): For a two-body decay :math:`R \to 12`, the *break-up momentum* is the absolute value of the momentum of both :math:`1` and :math:`2` in the rest frame of :math:`R`. See - Equation (50.7) on :pdg-review:`2024; Resonances; p.7`. + `PDG2026, Eq. (50.7) + `__. In AmpForm's standard implementation, the numerator is represented as a single square root. This results in :ref:`better computational performance @@ -59,11 +60,12 @@ class BreakupMomentumKallen(sp.Expr): """Two-body break-up momentum with a Källén function. This version of the `BreakupMomentum` represents the numerator using the `.Kallen` - function. This is common practice in literature (e.g. :pdg-review:`2024; Resonances; - p.7`), but results in a :ref:`more complicated cut - ` and :ref:`worse numerical - performance ` - than `BreakupMomentum`. + function. This is common practice in literature (e.g. `PDG2026, Eq. (50.7) + `__), but + results in a :ref:`more complicated cut ` and :ref:`worse numerical performance + ` than + `BreakupMomentum`. """ s: Any @@ -107,11 +109,14 @@ def _latex_repr_(self, printer: LatexPrinter, *args) -> str: @unevaluated class BreakupMomentumComplex(sp.Expr): - """Two-body break-up momentum with a square root that is defined on the real axis. + r"""Two-body break-up momentum with a square root that is defined on the real axis. In this version of the `BreakupMomentumSplitSqrt`, the square roots are replaced by `.ComplexSqrt`, which has a defined behavior for negative input values, so that it - can be evaluated on the entire real axis. + can be evaluated on the entire real axis. Between the pseudothreshold and the + threshold, this reproduces the analytic continuation :math:`q = i\sqrt{-q^2}` of + `PDG2026, Eq. (50.36) + `__. """ s: Any diff --git a/tests/dynamics/test_dynamics.py b/tests/dynamics/test_dynamics.py index 9d2a1e0c..406e6c66 100644 --- a/tests/dynamics/test_dynamics.py +++ b/tests/dynamics/test_dynamics.py @@ -109,6 +109,31 @@ def it_commutes_evaluation_with_parameter_substitution(method: str): def describe_BreitWigner(): + def it_has_unit_modulus_at_the_pole_with_a_mass_width_numerator(): + m0, w0, m1, m2 = sp.symbols("m0 Gamma0 m1 m2", positive=True) + breit_wigner = BreitWigner(m0**2, m0, w0, m1, m2, angular_momentum=1) + assert sp.simplify(breit_wigner.doit()) == sp.I + + def it_can_have_a_unity_numerator(): + s, m0, w0 = sp.symbols("s m0 Gamma0", nonnegative=True) + breit_wigner = BreitWigner(s, m0, w0, numerator="unity") + expected = 1 / (m0**2 - s - sp.I * m0 * w0) + assert breit_wigner.doit() == expected + + def it_marks_the_mass_width_numerator_with_a_hat(): + s, m0, w0 = sp.symbols("s m0 Gamma0", nonnegative=True) + arguments = R"_{L=0}\left(s; m_{0}, \Gamma_{0}\right)" + mass_width = BreitWigner(s, m0, w0) + unity = BreitWigner(s, m0, w0, numerator="unity") + assert sp.latex(mass_width) == R"\hat{\mathcal{R}}^\mathrm{BW}" + arguments + assert sp.latex(unity) == R"\mathcal{R}^\mathrm{BW}" + arguments + + def it_rejects_unknown_numerators(): + s, m0, w0 = sp.symbols("s m0 Gamma0", nonnegative=True) + breit_wigner = BreitWigner(s, m0, w0, numerator="mass") + with pytest.raises(ValueError, match=r"'mass-width', 'unity'"): + breit_wigner.doit() + def it_reduces_to_simple_breit_wigner(): s, m0, w0 = sp.symbols("s m0 Gamma0", nonnegative=True) breit_wigner = BreitWigner(s, m0, w0) @@ -163,6 +188,22 @@ def it_matches_serialized_l1405(): assert actual == pytest.approx(expected) +def describe_SimpleBreitWigner(): + def it_can_have_a_unity_numerator(): + s, m0, w0 = sp.symbols("s m0 Gamma0", nonnegative=True) + breit_wigner = SimpleBreitWigner(s, m0, w0, numerator="unity") + expected = 1 / (m0**2 - s - sp.I * m0 * w0) + assert breit_wigner.doit() == expected + + def it_marks_the_mass_width_numerator_with_a_hat(): + s, m0, w0 = sp.symbols("s m0 Gamma0", nonnegative=True) + arguments = R"\left(s; m_{0}, \Gamma_{0}\right)" + mass_width = SimpleBreitWigner(s, m0, w0) + unity = SimpleBreitWigner(s, m0, w0, numerator="unity") + assert sp.latex(mass_width) == R"\hat{\mathcal{R}}^\mathrm{BW}" + arguments + assert sp.latex(unity) == R"\mathcal{R}^\mathrm{BW}" + arguments + + def _subs(obj: sp.Basic, replacements: dict, method) -> sp.Expr: return getattr(obj, method)(replacements) diff --git a/tests/dynamics/test_form_factor.py b/tests/dynamics/test_form_factor.py index 8c89a362..855b520b 100644 --- a/tests/dynamics/test_form_factor.py +++ b/tests/dynamics/test_form_factor.py @@ -1,11 +1,76 @@ import pytest import sympy as sp -from ampform.dynamics.form_factor import _get_polynomial_blatt_weisskopf +from ampform.dynamics.form_factor import ( + BlattWeisskopfSquared, + FormFactor, + _get_polynomial_blatt_weisskopf, +) +from ampform.helicity import ParameterValues +from ampform.kinematics.phasespace import BreakupMomentumSquared z = sp.Symbol("z", nonnegative=True, real=True) +def describe_BlattWeisskopfSquared(): + @pytest.mark.parametrize( + ("ell", "expected"), + [ + (0, 1), + (1, z / (1 + z)), + (2, z**2 / (9 + 3 * z + z**2)), + ], + ) + def it_matches_the_pdg_convention_without_normalization( + ell: int, expected: sp.Expr + ): + expr = BlattWeisskopfSquared(z, angular_momentum=ell, normalize=False) + assert sp.simplify(expr.doit() - expected) == 0 + + def it_omits_the_normalization_for_symbolic_angular_momentum(): + ell = sp.Symbol("L", integer=True, nonnegative=True) + expr = BlattWeisskopfSquared(z, angular_momentum=ell, normalize=False).doit() + expected = z**2 / (9 + 3 * z + z**2) + assert sp.simplify(expr.subs(ell, 2).doit() - expected) == 0 + + def it_marks_the_normalized_function_with_a_hat(): + assert ( + sp.latex(BlattWeisskopfSquared(z, angular_momentum=1)) + == R"\hat{B}_{1}^2\left(z\right)" + ) + expr = BlattWeisskopfSquared(z, angular_momentum=1, normalize=False) + assert sp.latex(expr) == R"B_{1}^2\left(z\right)" + + +def describe_FormFactor(): + def it_forwards_the_normalization_flag(): + s, m1, m2, d = sp.symbols("s m1 m2 d", nonnegative=True) + form_factor = FormFactor( + s, m1, m2, angular_momentum=1, meson_radius=d, normalize=False + ) + q2 = BreakupMomentumSquared(s, m1, m2) + expected = BlattWeisskopfSquared(q2 * d**2, angular_momentum=1, normalize=False) + assert form_factor.evaluate() == sp.sqrt(expected) + + def it_marks_the_normalized_form_factor_with_a_hat(): + s, m1, m2 = sp.symbols("s m1 m2", nonnegative=True) + arguments = R"_{1}\left(s, m_{1}, m_{2}\right)" + normalized = FormFactor(s, m1, m2, angular_momentum=1) + non_normalized = FormFactor(s, m1, m2, angular_momentum=1, normalize=False) + assert sp.latex(normalized) == R"\hat{\mathcal{F}}" + arguments + assert sp.latex(non_normalized) == R"\mathcal{F}" + arguments + + def it_keeps_the_normalization_flag_when_substituting_parameter_values(): + s, m1, m2, d = sp.symbols("s m1 m2 d", nonnegative=True) + form_factor = FormFactor( + s, m1, m2, angular_momentum=2, meson_radius=d, normalize=False + ) + parameters = ParameterValues({d: 5.0, m1: 0.938, m2: 0.493}) + substituted = form_factor.xreplace(parameters) + assert substituted.normalize is False + assert substituted.meson_radius == 5.0 + + @pytest.mark.parametrize( ("ell", "expected"), [ diff --git a/uv.lock b/uv.lock index 10decf76..a3e7f595 100644 --- a/uv.lock +++ b/uv.lock @@ -92,7 +92,6 @@ dev = [ { name = "sphinx-copybutton" }, { name = "sphinx-design", version = "0.6.1", source = { registry = "https://pypi.org/simple" }, marker = "python_full_version < '3.11'" }, { name = "sphinx-design", version = "0.7.0", source = { registry = "https://pypi.org/simple" }, marker = "python_full_version >= '3.11'" }, - { name = "sphinx-hep-pdgref" }, { name = "sphinx-pybtex-etal-style" }, { name = "sphinx-thebe" }, { name = "sphinx-togglebutton" }, @@ -124,7 +123,6 @@ doc = [ { name = "sphinx-copybutton" }, { name = "sphinx-design", version = "0.6.1", source = { registry = "https://pypi.org/simple" }, marker = "python_full_version < '3.11'" }, { name = "sphinx-design", version = "0.7.0", source = { registry = "https://pypi.org/simple" }, marker = "python_full_version >= '3.11'" }, - { name = "sphinx-hep-pdgref" }, { name = "sphinx-pybtex-etal-style" }, { name = "sphinx-thebe" }, { name = "sphinx-togglebutton" }, @@ -255,7 +253,6 @@ dev = [ { name = "sphinx-comments" }, { name = "sphinx-copybutton" }, { name = "sphinx-design" }, - { name = "sphinx-hep-pdgref" }, { name = "sphinx-pybtex-etal-style" }, { name = "sphinx-thebe" }, { name = "sphinx-togglebutton" }, @@ -280,7 +277,6 @@ doc = [ { name = "sphinx-comments" }, { name = "sphinx-copybutton" }, { name = "sphinx-design" }, - { name = "sphinx-hep-pdgref" }, { name = "sphinx-pybtex-etal-style" }, { name = "sphinx-thebe" }, { name = "sphinx-togglebutton" }, @@ -5012,23 +5008,6 @@ wheels = [ { url = "https://files.pythonhosted.org/packages/30/cf/45dd359f6ca0c3762ce0490f681da242f0530c49c81050c035c016bfdd3a/sphinx_design-0.7.0-py3-none-any.whl", hash = "sha256:f82bf179951d58f55dca78ab3706aeafa496b741a91b1911d371441127d64282", size = 2220350, upload-time = "2026-01-19T13:12:51.077Z" }, ] -[[package]] -name = "sphinx-hep-pdgref" -version = "0.2.2" -source = { registry = "https://pypi.org/simple" } -dependencies = [ - { name = "attrs" }, - { name = "docutils", version = "0.21.2", source = { registry = "https://pypi.org/simple" }, marker = "python_full_version < '3.11'" }, - { name = "docutils", version = "0.22.4", source = { registry = "https://pypi.org/simple" }, marker = "python_full_version >= '3.11'" }, - { name = "sphinx", version = "8.1.3", source = { registry = "https://pypi.org/simple" }, marker = "python_full_version < '3.11'" }, - { name = "sphinx", version = "9.0.4", source = { registry = "https://pypi.org/simple" }, marker = "python_full_version == '3.11.*'" }, - { name = "sphinx", version = "9.1.0", source = { registry = "https://pypi.org/simple" }, marker = "python_full_version >= '3.12'" }, -] -sdist = { url = "https://files.pythonhosted.org/packages/2e/37/6a0e669bec2089db62e78c51900138fa1ea11cc51c8dc1ce433e1cb261af/sphinx_hep_pdgref-0.2.2.tar.gz", hash = "sha256:1257b50628bca25219112fca2ec3a5ab881c82b0cdc50be73eb56b2400c5ec2e", size = 17339, upload-time = "2026-02-06T12:03:29.602Z" } -wheels = [ - { url = "https://files.pythonhosted.org/packages/72/1e/fbdfe0b6131ac0d2b62052f2f44a48c8f02434ce16cef4383fcaad4fe5cf/sphinx_hep_pdgref-0.2.2-py3-none-any.whl", hash = "sha256:3b990e16ac449f967c2e54cbbeffdc9fda380e4112e9de9aac0fdc0aea285915", size = 7468, upload-time = "2026-02-06T12:03:28.512Z" }, -] - [[package]] name = "sphinx-pybtex-etal-style" version = "0.0.4"