# Fitting ultracold resonances without a fit

A. Simoni  Univ Rennes, CNRS, IPR (Institut de Physique de Rennes)-UMR 6251,
F-35000 Rennes, France

###### Abstract

We present a numerical procedure allowing one to extract Feshbach resonance
parameters from numerical calculations without relying on approximate fitting
procedures. Our approach is based on a simple decomposition of the reactance
matrix in terms of poles and residual background contribution, and can be
applied to the general situation of inelastic overlapping resonances. A simple
lineshape for overlapping inelastic resonances, equivalent to known results in
the particular cases of isolated and overlapping elastic features, is also
rigorously derived.

##  I Introduction

Feshbach resonances are an ubiquitous phenomenon taking place when the energy
of colliding particles is nearly degenerate with a quasi-bound level of the
system. Resonances are traditionally observed as a function of collision
energy or also by varying applied electromagnetic fields. In ultracold gases,
it has been long realized that resonances occurring at extremely low energies
can be induced by tuning an external magnetic field Tiesinga  _et al._ (1993),
allowing one to vary in a controlled way the effective interatomic
interaction. The possibility to control the interaction using resonances is at
the hearth of a series of groundbreaking experiments in the field of ultracold
gases; See e.g. Bloch  _et al._ (2008).

Feshbach resonances have been studied both theoretically and experimentally in
a number of ultracold atomic systems Chin  _et al._ (2010). Theory usually
conveys information extracted from numerical models in a small number of
lineshape parameters, namely position, width and the non-resonant (or
background) contribution, extracted by fitting the numerical result to an
analytical form. A celebrated example is the case of an isolated elastic
resonance Tiesinga  _et al._ (1993) where the $s$-wave scattering length
diverges both positively and negatively as a function of the magnetic field
$B$ across the center $B_{0}$ :

$a(B)=a_{\rm bg}\left(1-\dfrac{\Delta}{B-B_{0}}\right).$ (1)

Here $\Delta$ represents the magnetic width and the background scattering
length $a_{\rm bg}$ is the off-resonance value of $a(B)$.

A more general lineshape describing a pair and an arbitrary number of elastic
overlapping resonances has been derived using formal scattering theory and
quantum defect analysis in the Supplemental Material of Roy  _et al._ (2013)
and in Jachymski and Juilenne (2013), respectively. In this case, the
scattering length is still a purely real quantity showing a pair or a series
of nearby divergences.

Isolated resonances in the presence of inelastic processes have been discussed
in van Kempen  _et al._ (2002) and studied in greater detail in Hutson and
Soldán (2007); Naik  _et al._ (2011). The main influence of inelastic loss
channels is to make the scattering length a complex quantity and to suppress
the pole in its real part, that no longer diverges Ciurylo  _et al._ (2005);
Hutson and Soldán (2007). In the presence of inelastic processes both the
lineshape and the extraction of related resonance parameters become more
involved. A relatively elaborate fitting strategy based on a three-points
iterative procedure has been developed in Frye and Hutson (2017). However,
like any other approach relying on fitting the lineshape to an analytical
expression, the latter procedure presents few intrinsic drawbacks.

First, the precise value of the fitting parameters depends on the chosen
fitting interval and bears an error bar. Second, the analytical lineshape can
be a poor approximation to the fitted data, even to the extent that the
fitting procedures may not converge Frye and Hutson (2017). Even if the least
square iterations do converge, the physical meaning of the fitting parameters
will become unclear if the analytical lineshape poorly describes the data. For
instance, if the background scattering length presents an unusually rapid
variation, defining such quantity as the off-resonance value of $a(B)$ is not
a well defined prescription. Finally, the analytical lineshape for overlapping
resonances in the presence of inelastic processes is not yet available. This
situation which indeed occurs frequently for collisions in excited states Frye
_et al._ (2019) is still missing a systematic procedure allowing one to
summarize the numerical data in few theoretical parameters. A simple approach
consisting in presenting graphics of the numerical results without
accompanying tables of resonance parameters is not fully satisfactory, since
reading quantitative information from a plot without having access to the
rough data can be difficult or impossible.

In this paper, we present a formalism designed to overcome all of the above
difficulties. Our approach allows one to define and extract in an
algorithmically well defined way resonance parameters even in cases where the
lineshape cannot be modeled by a simple expression. For the sake of better
physical insight and clarity, our discussion is based on Feshbach resonance
theory Feshbach (1967). However, at a fundamental level the present strategy
is based on a general pole approximation that could be taken as starting point
in a similar fashion as scattering resonances can be defined as scattering-
matrix poles Taylor (2006). Special attention is granted to stable ways to
implement such pole decomposition, a numerically non trivial task in the
presence of naturally diverging quantities such as the scattering length on
resonance. As a by-product of the formalism, we also demonstrate that the
complex scattering length can always be parametrized in a simple and intuitive
way in terms of complex magnetic field locations and magnetic widths.

The paper is organized as follows. In section II we introduce the formalism
building on Feshbach resonance theory.

Section III proves that the present approach reduces under suitable conditions
to results previously appeared in the literature.

In section IV we present an application to two situations of physical interest
for 39K, the case of overlapping and inelastic $s$-wave resonances and of
finite-energy $p$-wave resonances.

A short conclusion ends this work.

##  II Formalism

###  II.1 Basic ideas for single channel scattering

The main idea can be illustrated with reference to the simple single channel
resonant expression of the $s$-wave scattering length in the presence of an
isolated resonance. Here we consider the most common case that the scattering
length is tuned via an external magnetic field $B$, but the formalism applies
equally well to a general external fields or variation of any physical
parameter in the Hamiltonian capable of inducing resonances. In the case of a
magnetic field, recalling (1), one can write :

$a(B)=a_{\rm bg}\left(1-\dfrac{\Delta}{B-B_{0}}\right)\equiv a_{\rm
bg}+\dfrac{p}{B-B_{0}}$ (2)

where we have defined $p=-a_{\rm bg}\Delta$ the field-independent pole
strength. Any variation of $a$ with $B$ beyond the pole contribution is
included in the background term, which will therefore be considered a function
of $B$.

The pole strength can be extracted from a theoretical calculation of $a$ as

$p={\lim}_{B\to B_{0}}\left[a(B)(B-B_{0})\right].$ (3)

It will be shown below how this limit can be performed in a numerically stable
way; See equation (16). Once the $p$ parameter is determined, one can compute
the on-resonance value of the background scattering length according to

$a_{\rm bg}(B_{0})={\lim}_{B\to B_{0}}\left[a(B)-\dfrac{p}{B-B_{0}}\right].$
(4)

Again, this limit can be computed in a stable way; See (19). We stress that
the usual interpretation of the background scattering length as the limit
value of $a$ far away from resonance may become ambiguous in cases where
$a_{\rm bg}$ is not strictly constant or, even worse, is rapidly varying. On
the converse, the quantity (4) is well defined in any case and denotes simply
the value of $a(B)$ at the resonance location after subtraction of the
divergent part. We now generalize these ideas to multichannel scattering.

###  II.2 Pole decomposition for multichannel scattering

Les us consider a general scattering problem with $N$ channel, $N_{c}$ closed
channels, and $N_{op}=(N-N_{c})$ open channels. We will further partition the
$N_{op}$ open channels into $N_{o}$ degenerate (also termed in scattering
theory super-elastic) channels with energy equal to the incoming channel, and
$N_{\bar{o}}$ inelastic channels with energy release (or absorption). We use
Feshbach theory, that formally partitions channel space into open and closed
channel subspaces $P$ and $Q$ Feshbach (1958). Scattering in the $P$ subspace
is supposed to be known, and the influence of the $Q$ subspace is taken into
account through the introduction of an effective Hamiltonian. The nature of
$P$ and $Q$ subspaces depends on the specific problem at hand and for general
formal manipulations does not even need to be specified Domcke (1983). A
resonance occurs when the energy of bound states in $Q$ subspace is quasi-
degenerate with the energy of the colliding particles. Under such conditions
the influence of the $Q$ subspace on the total scattering amplitude becomes
dramatic.

Feshbach formalism is most often applied to the scattering matrix Feshbach
(1958), but also provides a partition of the reactance matrix Feshbach (1967).
Working with the reactance matrix is convenient since the resulting scattering
matrix is automatically unitary, algebra is real, and the number of fitting
parameters is minimal Feshbach (1967). According to Feshbach theory, the
reactance matrix at total energy $E$ can be written

$\displaystyle{\mathbf{K}}$ $\displaystyle=$
$\displaystyle\left[\sin({\bm{\xi}}_{\rm bg})+\cos({\bm{\xi}}_{\rm
bg}){\mathbf{K}}_{\rm res}\right]\left[\cos({\bm{\xi}}_{\rm
bg})-\sin({\bm{\xi}}_{\rm bg}){\mathbf{K}}_{\rm res}\right]^{-1}$ (5)
$\displaystyle=$ $\displaystyle\cos({\bm{\xi}_{\rm
bg}})\left[\tan({\bm{\xi}}_{\rm bg})+{\mathbf{K}}_{\rm
res}\right]\left[{\bm{1}}_{P}-\tan({\bm{\xi}}_{\rm bg}){\mathbf{K}}_{\rm
res}\right]^{-1}\cos({\bm{\xi}_{\rm bg}})^{-1}$ (6)

where ${\bm{1}}_{P}$ is the identity in $P$ space and the resonant part of the
reactance matrix has the form

${\mathbf{K}}_{\rm
res}={\mathbf{y}}\left(E{\bm{1}}_{Q}-{\mathbf{E}}_{b}-{\mathbf{V}}\right)^{-1}{\mathbf{y}}^{t}.$
(7)

Here $\mathbf{y}$ is the $N_{op}\times N_{c}$ matrix formed by half
predissociation amplitudes and the operator
$\left(E{\bm{1}}_{Q}-{\mathbf{E}}_{b}-{\mathbf{V}}\right)$ acts in the closed
channel subspace, the $Q$-space of Feshbach theory. Zero energy is taken at
the energy of the separated atoms at the given value of the applied field. The
diagonal matrix ${\mathbf{E}}_{b}$ represents the energy levels of the bare
molecular states defining the $Q$-space. We will make the usual assumption
that the position of such levels depends linearly on the applied field
${\mathbf{E}}_{b}=-{\bm{\mu}}(B{\bm{1}}_{Q}-{\mathbf{B}}_{r})$ with $\bm{\mu}$
a diagonal matrix whose entries $\mu_{\alpha}$ are the differential magnetic
moment of state $\alpha$ with respect to the one of two separated atoms in the
initial channel. The entries $B_{r,\alpha}$ of the diagonal matrix
${\mathbf{B}}_{r}$ represent the bare resonance magnetic field, i.e the field
at which level $\alpha$ crosses the incoming channel dissociation threshold.
The coupling $\mathbf{V}$ contains in general an energy-independent direct
coupling $\mathbf{U}$ between molecular states as well as indirect coupling
through the continuum via the level shift operator ${\bm{\delta}}^{E}$ at
energy $E$ Cohen-Tannoudji  _et al._ (1992); Simoni  _et al._ (2002). At this
stage it is not necessary to distinguish the two contributions, only the
Hermitic nature of $\mathbf{V}$ being of importance.

To proceed further with the derivation, inspired by the single-channel case,
we factor out the magnetic moment in the denominator Tiesinga  _et al._
(1993). To this aim, we define
$\tilde{\mathbf{y}}={\bm{\mu}}^{-1/2}{\mathbf{y}}$ and
${\bm{\mathcal{B}}}={\mathbf{B}}_{r}+{\bm{\mu}}^{-1/2}({\mathbf{V}}-E{\bm{1}}_{Q}){\bm{\mu}}^{-1/2}$.
The resonant reactance matrix is expressed in terms of these redefined
quantities as :

${\mathbf{K}}_{\rm
res}=\tilde{\mathbf{y}}\left(B{\bm{1}}_{Q}-{\bm{\mathcal{B}}}\right)^{-1}\tilde{\mathbf{y}}^{t}.$
(8)

If the magnetic moments all share the same sign the matrix ${\bm{\mu}}^{1/2}$
will be either purely real or purely imaginary. In this case, by its
definition $\bm{\mathcal{B}}$ will be real and symmetric, thus with real
eigenvalues. If the magnetic moments present both negative or positive values
the matrix $\bm{\mathcal{B}}$ will in general be complex, yet still symmetric.

We now insert (8) into (6). After some algebra, see proof in Appendix A, one
can recast the $K$ matrix in the form

${\mathbf{K}}=\tan({\bm{\xi}}_{\rm bg})+\cos({\bm{\xi}}_{\rm
bg})^{-1}\tilde{\mathbf{y}}(B{\bm{1}}_{Q}-\hat{\bm{\mathcal{B}}})^{-1}\tilde{\mathbf{y}}^{t}\cos({\bm{\xi}}_{\rm
bg})^{-1}$ (9)

where the shifted resonance position matrix $\hat{\bm{\mathcal{B}}}$ is
defined as
$\hat{\bm{\mathcal{B}}}=\bm{\mathcal{B}}+\tilde{\mathbf{y}}^{t}\tan({\bm{\xi}}_{\rm
bg})\tilde{\mathbf{y}}$. If there are no inelastic channels and if threshold
collisions are considered ($E\to 0$), all elements of $\tan({\bm{\xi}}_{\rm
bg})$ as well as of the predissociation amplitude $\mathbf{y}$ tend to zero,
hence $\hat{\bm{\mathcal{B}}}\to\bm{\mathcal{B}}$. The
$\hat{\bm{\mathcal{B}}}$ operator is still complex symmetric, is real if
${\bm{\mathcal{B}}}$ is real, but unlike ${\bm{\mathcal{B}}}$ in general it
varies with magnetic field because of the background phase ${\bm{\xi}}_{\rm
bg}$ appearing in its definition. It can therefore be brought a diagonal form
$\hat{\mathbf{b}}$ through a $B$-dependent complex orthogonal transformation
$\mathbf{C}$ such that
$\hat{\bm{\mathcal{B}}}={\mathbf{C}}\hat{\mathbf{b}}{\mathbf{C}}^{t}$ and
${\mathbf{C}}^{t}{\mathbf{C}}={\mathbf{C}}{\mathbf{C}}^{t}={\mathbf{1}}$
Wilkinson (1965). With this transformation, equation (9) can be rewritten :

${\mathbf{K}}=\tan({\bm{\xi}}_{\rm
bg})+\hat{\mathbf{y}}(B{\bm{1}}_{Q}-\hat{\mathbf{b}})^{-1}\hat{\mathbf{y}}^{t}=\tan({\bm{\xi}}_{\rm
bg})+\sum_{\alpha}\frac{\hat{\mathbf{y}}_{\alpha}\hat{\mathbf{y}}^{t}_{\alpha}}{B-{\hat{b}}_{\alpha}}$
(10)

in terms of a redefined predissociation amplitude
$\hat{\mathbf{y}}=\cos({\bm{\xi}}_{\rm
bg})^{-1}\tilde{\mathbf{y}}{\mathbf{C}}^{t}$, which also depends on $B$.

The reactance matrix clearly diverges at the solutions, say $b_{\alpha}$, of
the nonlinear equation $B-{\hat{b}}_{\alpha}(B)=0$. The following
manipulations will turn out to be simpler if one isolates pole contributions
with a field-independent pole residue. To this aim, let us define the constant
amplitude vector
$\hat{\mathbf{y}}^{(0)}_{\alpha}=\hat{\mathbf{y}}_{\alpha}(b_{\alpha})$ and
rewrite (10)

${\mathbf{K}}=\tan(\tilde{\bm{\xi}}_{\rm
bg})+\hat{\mathbf{y}}^{(0)}(B{\bm{1}}_{Q}-{\mathbf{b}})^{-1}\hat{\mathbf{y}}^{(0)t}=\tan(\tilde{\bm{\xi}}_{\rm
bg})+\sum_{\alpha}\frac{\hat{\mathbf{y}}^{(0)}_{\alpha}\hat{\mathbf{y}}^{(0)t}_{\alpha}}{B-b_{\alpha}}.$
(11)

where the modified background

$\tan(\tilde{\bm{\xi}}_{\rm bg})=\tan({\bm{\xi}}_{\rm
bg})-\hat{\mathbf{y}}^{(0)}(B{\bm{1}}_{Q}-{\mathbf{b}})^{-1}\hat{\mathbf{y}}^{(0)t}+\hat{\mathbf{y}}(B{\bm{1}}_{Q}-\hat{\mathbf{b}})^{-1}\hat{\mathbf{y}}^{t}.$
(12)

is still regular at the poles and reduces to $\tan({\bm{\xi}}_{\rm bg})$ if
${\bm{\xi}}_{\rm bg}$ is strictly constant.

According to the remark below (8) only real values of $b_{\alpha}$ are
possible if all entries in matrix ${\bm{\mu}}$ have the same sign. If this is
not the case, it can be shown that solutions $b_{\alpha}$ are either real or
occur in pairs of complex conjugate quantities. On a physical ground, such
complex solutions arise when a pair of molecular levels with opposite slopes
avoid each other exactly at the collision threshold, such that neither level
crosses the latter. For purely elastic scattering, this unlikely case will
result in a scattering length that may have a rapid variation, but will not
diverge as $B$ varies on the real axis. We will not consider this situation as
a resonance condition and will henceforth assume that $b_{\alpha}$ is always
real, with the implicit assumption that the contribution from complex
solutions $b_{\alpha}$, if any, is included in the background.

Since we are mainly interested in scattering near a collision threshold, it is
convenient to factor out the low-energy limit behavior of the reactance matrix
for channels with orbital angular momentum $\ell$ as follows Dalitz (1961):

$\bar{\mathbf{K}}={\mathbf{q}}^{-(\ell+1/2)}{\mathbf{K}}{\mathbf{q}}^{-(\ell+1/2)}.$
(13)

The diagonal normalization matrix ${\mathbf{q}}^{-(\ell+1/2)}$ has by
definition elements $q_{i}^{-(\ell_{i}+1/2)}$ for all channels degenerate with
the incoming one, and to 1 for inelastic ones. The resulting normalized
$K$-matrix is thus finite at threshold and is the main quantity needed to
compute zero-energy scattering lengths. Performing such channel normalization
and introducing the background scattering length matrix
${\mathbf{A}}={\mathbf{q}}^{-(\ell+1/2)}\tan(\tilde{\bm{\xi}}_{\rm
bg}){\mathbf{q}}^{-(\ell+1/2)}$ as well as the resonance amplitude
$\bar{\mathbf{y}}={\mathbf{q}}^{-(\ell+1/2)}\hat{\mathbf{y}}^{(0)}$, equation
(11) transforms into

$\bar{\mathbf{K}}={\mathbf{A}}+\bar{\mathbf{y}}(B{\bm{1}}_{Q}-{\mathbf{b}})^{-1}\bar{\mathbf{y}}^{t}={\mathbf{A}}+\sum_{\alpha}\frac{\bar{\mathbf{y}}_{\alpha}\bar{\mathbf{y}}^{t}_{\alpha}}{B-b_{\alpha}}.$
(14)

Note that according to the Wigner laws Wigner1948, by construction both
${\mathbf{A}}$ and $\bar{\mathbf{y}}$ are finite at threshold.

With the assumption above, the matrix $\bar{\mathbf{K}}$ diverges at the real
magnetic field locations $b_{\alpha}$. This should be contrasted with the
behavior of the complex scattering length, that even in the presence of a pole
in $\bar{\mathbf{K}}$ on the real axis can in some cases present barely no
structure. In order to locate the poles numerically, we find it more
convenient to consider the inverse reactance matrix
${\mathbf{M}}={\bar{\mathbf{K}}}^{-1}$ and to search for the zeros in the
determinant $\det{\mathbf{M}}$. We first bracket zeros by a coarse search,
then quickly refine their position using Ridder’s method. Such nonlinear root-
finding algorithm has as main advantages a superlinear convergence and that of
keeping the root bracketed Press  _et al._ (1996).

We now show that the residue as well as the background contribution at the
pole position (recall that ${\mathbf{A}}$ will, in general, have a dependence
on $B$) can be accurately determined. Concerning the residue, by definition :

$\bar{\mathbf{y}}_{\alpha}\bar{\mathbf{y}}^{t}_{\alpha}={\lim}_{B\to
b_{\alpha}}(B-b_{\alpha})\bar{\mathbf{K}}(B).$ (15)

In order to evaluate the indeterminate $0\cdot\infty$ product, we write
$\bar{\mathbf{K}}={\mathbf{M}}^{-1}=\text{adj}{~{}\mathbf{M}}/\det{\mathbf{M}}$
and use l’Hôpital rule. Keeping into account continuity of the adjugate, one
promptly obtains

$\bar{\mathbf{y}}_{\alpha}\bar{\mathbf{y}}^{t}_{\alpha}=\dfrac{\text{adj}{~{}\mathbf{M}(b_{\alpha})}}{\det^{\prime}{\mathbf{M}}(b_{\alpha})}.$
(16)

Equation (16) gives a numerically stable way to evaluate the residue in terms
of $\mathbf{M}$. The needed derivatives of the determinant can be evaluated
using Jacobi’s formula :

${\det}^{\prime}{\mathbf{M}}(B)=\text{tr}\left(\text{adj}({\mathbf{M}}(B)){\mathbf{M}}^{\prime}\right)$
(17)

The local background at the pole is defined by subtracting the divergent
contribution from the global reactance matrix :

${\mathbf{A}}_{\text{loc}}(b_{\alpha})={\lim}_{B\to
b_{\alpha}}\left(\bar{\mathbf{K}}(B)-\frac{\bar{\mathbf{y}}_{\alpha}\bar{\mathbf{y}}^{t}_{\alpha}}{B-b_{\alpha}}\right).$
(18)

Note that ${\mathbf{A}}_{\text{loc}}$ contains contributions from the global
background ${\mathbf{A}}$ as well as from all resonances but the $\alpha$-th
one. The lhs presents itself in a $\infty-\infty$ form. In order to evaluate
the limit, we write again $\bar{\mathbf{K}}$ in terms of the adjugate and
determinant of $\mathbf{M}$. After bringing to common denominator and Taylor
expanding to second order one obtains

${\mathbf{A}}_{\text{loc}}(b_{\alpha})=\dfrac{2\det^{\prime}({\mathbf{M}}(b_{\alpha}))\text{adj}^{\prime}({\mathbf{M}(b_{\alpha})})-\det^{\prime\prime}{\mathbf{M}}(b_{\alpha})\text{adj}({\mathbf{M}(b_{\alpha})}}{2\left[\det^{\prime}{\mathbf{M}}(b_{\alpha})\right]^{2}}.$
(19)

One proceeds similarly for the first derivative, which will also be needed.
One starts from the definition

${\mathbf{A}}_{\text{loc}}^{\prime}(b_{\alpha})={\lim}_{b\to
b_{\alpha}}\left(\bar{\mathbf{K}}^{\prime}(B)+\frac{\bar{\mathbf{y}}_{\alpha}\bar{\mathbf{y}}^{t}_{\alpha}}{(B-b_{\alpha})^{2}}\right)$
(20)

and evaluates the indeterminate limit using Taylor expansion to get

${\mathbf{A}}_{\text{loc}}^{\prime}(b_{\alpha})=\dfrac{\begin{multlined}2(\det^{\prime\prime}{\mathbf{M}}(b_{\alpha}))^{2}\text{adj}({\mathbf{M}(b_{\alpha})}-2\det^{\prime}({\mathbf{M}}(b_{\alpha}))\det^{\prime\prime}({\mathbf{M}}(b_{\alpha}))\text{adj}^{\prime}({\mathbf{M}(b_{\alpha})})+2(\det^{\prime}{\mathbf{M}}(b_{\alpha}))^{2}\text{adj}^{\prime\prime}({\mathbf{M}(b_{\alpha})})\\\
-\det^{\prime}({\mathbf{M}}(b_{\alpha}))\det^{\prime\prime\prime}({\mathbf{M}}(b_{\alpha}))\text{adj}({\mathbf{M}(b_{\alpha})})-\det^{\prime}({\mathbf{M}}(b_{\alpha}))\det^{\prime\prime}({\mathbf{M}}(b_{\alpha}))\text{adj}^{\prime\prime}({\mathbf{M}(b_{\alpha})})\end{multlined}2(\det^{\prime\prime}{\mathbf{M}}(b_{\alpha}))^{2}\text{adj}({\mathbf{M}(b_{\alpha})}-2\det^{\prime}({\mathbf{M}}(b_{\alpha}))\det^{\prime\prime}({\mathbf{M}}(b_{\alpha}))\text{adj}^{\prime}({\mathbf{M}(b_{\alpha})})+2(\det^{\prime}{\mathbf{M}}(b_{\alpha}))^{2}\text{adj}^{\prime\prime}({\mathbf{M}(b_{\alpha})})\\\
-\det^{\prime}({\mathbf{M}}(b_{\alpha}))\det^{\prime\prime\prime}({\mathbf{M}}(b_{\alpha}))\text{adj}({\mathbf{M}(b_{\alpha})})-\det^{\prime}({\mathbf{M}}(b_{\alpha}))\det^{\prime\prime}({\mathbf{M}}(b_{\alpha}))\text{adj}^{\prime\prime}({\mathbf{M}(b_{\alpha})})}{2(\det^{\prime}{\mathbf{M}}(b_{\alpha}))^{3}}.$
(21)

The required derivatives of the adjugate can be evaluated keeping into account
that elements of the adjugate are in turn determinants of minors, and using
again (17). Higher order derivatives of a determinant (and thus of the
adjugate) follow by successive differentiation of (17). For instance, the
second derivative can be calculated as

${\det}^{\prime\prime}{\mathbf{M}}(B)=\text{tr}\left(\text{adj}^{\prime}({\mathbf{M}}(B)){\mathbf{M}}^{\prime}+\text{adj}({\mathbf{M}}(B)){\mathbf{M}}^{\prime\prime}\right).$
(22)

Once the local background matrix (which includes contributions from nearby
resonances) has been determined, the usually smooth global background and its
derivative are respectively given by

${\mathbf{A}}(b_{\alpha})={\mathbf{A}}_{\text{loc}}(b_{\alpha})-\sum_{\beta\neq\alpha}\frac{\bar{\mathbf{y}}_{\beta}\bar{\mathbf{y}}^{t}_{\beta}}{b_{\alpha}-b_{\beta}}.$
(23)

and

${\mathbf{A}}^{\prime}(b_{\alpha})={\mathbf{A}}_{\text{loc}}^{\prime}(b_{\alpha})+\sum_{\beta\neq\alpha}\frac{\bar{\mathbf{y}}_{\beta}\bar{\mathbf{y}}^{t}_{\beta}}{(b_{\alpha}-b_{\beta})^{2}}.$
(24)

In practice, it is clear that in these expansions one will not need to include
all poles of $\mathbf{K}$. One rather focuses on resonances occurring in a
predefined magnetic field region and includes all poles such that the addition
of additional ones has negligible influence on the parameters of the
resonances of interest. In simpler terms, only resonances significantly
overlapping with the ones under investigation must be taken into account. Is
is also clear that without a stable and accurate algorithm to compute the
$\mathbf{M}$ matrix and its derivatives the present formalism would be of
little practical use. Our numerical approach for computing such derivatives is
the subject of the following subsection.

###  II.3 Computing the derivatives

On a general ground, a basis representation of the wavefunction is
particularly suited to compute derivatives of any order, since the latter can
then be obtained in an essentially analytical fashion. In fact, while
propagation algorithms also can be adapted for obtaining derivatives (see e.g.
Aquilanti  _et al._ (2005)), algorithmic complexity and computational cost
both increase with the derivative order. In this work we exploit the
pseudospectral representation termed spectral element method, a spatial grid
representation presented for instance in Simoni  _et al._ (2017) and
references therein. In the present work we impose asymptotic $R$-matrix
boundary conditions, but one could equally well work with log-derivative or
scattering boundary conditions. In the discrete spatial grid representation,
the time-independent multichannel Schrödinger equation reads :

$[E{\mathbf{1}}-{\mathbf{H}}(B)]{\mathbf{\Psi}}(B)={\mathbf{c}}$ (25)

where the constant vector $\mathbf{c}$ on the rhs equal to zero everywhere but
at last grid point $r_{\rm max}$, where it equals the unit matrix of dimension
the total number of channels Simoni  _et al._ (2017). Elements of the
Hamiltonian matrix in the spectral element representation can be found in
Simoni  _et al._ (2017). By definition, the multichannel wavefunction
$\mathbf{\Psi}$ at $r_{\rm max}$ equals the $R$-matrix $\mathbf{R}$.

Once ${\mathbf{\Psi}}$ has been determined from equation (25), first
derivative of the matrix wavefunction is obtained from differentiating (25)
with respect to parameter $B$, and solving the resulting linear system:

${\mathbf{H}}(B)\frac{\partial\mathbf{\Psi}}{\partial
B}+\frac{\partial\mathbf{H}}{\partial B}{\mathbf{\Psi}}(B)=0.$ (26)

Again, the derivative of the $R$-matrix with respect to $B$ corresponds to the
computed ${\partial\mathbf{\Psi}}/{\partial B}$ at $r_{\rm max}$.

The $K$-matrix is ordinarily computed by imposing the asymptotic form of the
wavefunction ${\mathbf{\Psi}}\sim{\mathbf{f} }-{\mathbf{g}}{\mathbf{K}}$ in
terms of regular and irregular solutions $\mathbf{f}$ and $\mathbf{g}$ for the
free particle problem. Since in our algorithm $\mathbf{K}$ is not needed, we
rather determine directly its inverse $\mathbf{M}$ by matching to
${\mathbf{\Psi}}\sim{\mathbf{f} }{\mathbf{M}}-{\mathbf{g}}$. The linear
matching equations for $\mathbf{M}$ simply follow from equation (25) of Simoni
_et al._ (2017) and read :

$({\mathbf{f}}-{\mathbf{R}}{\mathbf{f}}^{\prime}){\mathbf{M}}=({\mathbf{g}}-{\mathbf{R}}{\mathbf{g}}^{\prime})$
(27)

the prime denotes the spatial derivative and all spatial arguments are tacitly
taken at last point of the discretization grid. The first derivative of
${\mathbf{M}}$ with respect to the field is solution to the linear system
obtained taking the derivative of (27)

$({\mathbf{f}}-{\mathbf{R}}{\mathbf{f}}^{\prime})\frac{\partial\mathbf{M}}{\partial
B}+\frac{\partial({\mathbf{f}}-{\mathbf{R}}{\mathbf{f}}^{\prime})}{\partial
B}{\mathbf{M}}=\frac{\partial\mathbf{g}}{\partial
B}-{\mathbf{R}}\frac{\partial{\mathbf{g}}^{\prime}}{\partial
B}-\frac{\partial\mathbf{R}}{\partial B}{\mathbf{g}}^{\prime}$ (28)

in which $\mathbf{M}$, $\mathbf{R}$, and $\dfrac{\partial\mathbf{R}}{\partial
B}$ are known. Derivatives of the reference functions with respect to the
field are transformed into derivatives with respect to energy using chain rule
and the known behavior of asymptotic channel energies with respect to $B$. The
latter are calculated using a stable formula for both open and closed channels
from Aquilanti  _et al._ (2005), equations (11) and (12).

Formal generalization to the $n$-th derivative of $\mathbf{\Psi}$ is
straightforward. One differentiates $n$ times the Schrödinger equation, and
solves for the highest order derivatives of the wavefunction

$\sum_{k}\binom{n}{k}\frac{\partial^{k}\mathbf{H}}{\partial
B^{k}}\frac{\partial^{n-k}\mathbf{\Psi}}{\partial B^{n-k}}=0$ (29)

knowing the derivatives of $\mathbf{\Psi}$ up to order $n-1$. Similarly, to
determine $\partial^{n}{\mathbf{M}}/\partial B^{n}$ one differentiates $n$
times the matching condition (27)

$\sum_{k}\binom{n}{k}\frac{\partial^{k}{(\mathbf{f}-{\mathbf{R}}{\mathbf{f}}^{\prime})}}{\partial
B^{k}}\frac{\partial^{n-k}\mathbf{M}}{\partial
B^{n-k}}=\frac{\partial^{n}\mathbf{g}}{\partial
B^{n}}-\sum_{k}\binom{n}{k}\frac{\partial^{k}{\mathbf{R}}}{\partial
B^{k}}\frac{\partial^{n-k}\mathbf{g}^{\prime}}{\partial B^{n-k}}$ (30)

and solves for the highest order derivative knowing the derivatives of
$\mathbf{M}$ up to order $n-1$ as well as the derivatives of $\mathbf{R}$ up
to order $n$. Higher order derivatives of reference functions with respect to
$B$ are obtained by successive differentiation of equations (11) and (12) in
Aquilanti  _et al._ (2005).

###  II.4 Observables : the energy-dependent scattering length

We now turn to the calculation of observable quantities. The scattering cross-
sections at arbitrary energy can be expressed exactly as a function of the so-
called energy-dependent scattering length Bolda  _et al._ (2003). Its
definition for a single channel in terms of the diagonal scattering matrix
element $S_{jj}$ :

$a(B)=q_{o}^{-(2\ell_{i}+1)}\frac{1-S_{jj}}{i(1+S_{jj})}$ (31)

where $q_{o}$ is the channel wavevector.

We consider here the general situation of a number $N_{o}$ of degenerate
channels occurring for example when several partial waves contribute to the
scattering process or in the presence of magnetic Zeeman sublevels in a
vanishing magnetic fields. Hence we define a square complex scattering length
matrix

${\mathbf{a}}(B)={{\mathbf{q}_{o}}^{-(\ell_{i}+1/2)}}\frac{{\mathbf{1}}_{oo}-{\mathbf{S}}_{oo}}{i\left({\mathbf{1}}_{oo}+{\mathbf{S}}_{oo}\right)}{{\mathbf{q}_{o}}^{-(\ell_{i}+1/2)}}.$
(32)

The off-diagonal elements of $\mathbf{a}$ are particularly relevant to
parametrize anisotropic pseudopotentials Derevianko (2003) or in the context
of spinor Bose-Einstein condensates Kawaguchi and Ueda (2012).

We now write the relevant square matrices
${\mathbf{X}}={\mathbf{K}},{\mathbf{A}},{\mathbf{1}}$ as well as the generally
non-square predissociation amplitudes $\mathbf{y}$ in blocked form for open
degenerate and non-degenerate channels :

${\mathbf{X}}=\left[\begin{array}[]{c|c}{\mathbf{X}}_{oo}&{\mathbf{X}}_{o{\bar{o}}}\\\
\hline\cr{\mathbf{X}}_{{\bar{o}}o}&{\mathbf{X}}_{{\bar{o}}{\bar{o}}}\end{array}\right]~{},~{}{\mathbf{\bar{y}}}=\left[\begin{array}[]{c}{\mathbf{\bar{y}}}_{o}\\\
\hline\cr{\mathbf{\bar{y}}}_{{\bar{o}}}\end{array}\right].$ (33)

With this decomposition and applying the projection technique of Dalitz (1961)
one can express the middle matrix ratio on the rhs of (32) as a function of
the blocks of the reactance matrix

$\frac{{\mathbf{1}}_{oo}-{\mathbf{S}}_{oo}}{i({\mathbf{1}}_{oo}+{\mathbf{S}}_{oo})}=-{\mathbf{K}}_{oo}-i
{\mathbf{K}}_{o{\bar{o}}}\left({\mathbf{1}}_{{\bar{o}}{\bar{o}}}-i
{\mathbf{K}}_{{\bar{o}}{\bar{o}}}\right)^{-1} {\mathbf{K}}_{{\bar{o}}o}$ (34)

Transforming to the normalized quantities according to (13)

${\mathbf{a}}(B)=-\bar{\mathbf{K}}_{oo}-i\bar{\mathbf{K}}_{o{\bar{o}}}\left({\mathbf{1}}_{{\bar{o}}{\bar{o}}}-i\bar{\mathbf{K}}_{{\bar{o}}{\bar{o}}}\right)^{-1}\bar{\mathbf{K}}_{{\bar{o}}o}$
(35)

and with the help of (14) we obtain :

$\displaystyle{\mathbf{a}}(B)$ $\displaystyle=$
$\displaystyle-{\mathbf{A}}_{oo}-\bar{\mathbf{y}}_{o}\left(B{\mathbf{1}}_{Q}-{\mathbf{b}}\right)^{-1}\bar{\mathbf{y}}^{t}_{o}$
(36) $\displaystyle-$ $\displaystyle
i\left({\mathbf{A}}_{o\bar{o}}+\bar{\mathbf{y}}_{o}\left(B{\mathbf{1}}_{Q}-{\mathbf{b}}\right)^{-1}\bar{\mathbf{y}}^{t}_{\bar{o}}\right)\left({\mathbf{1}}_{\bar{o}\bar{o}}-i{\mathbf{A}}_{\bar{o}\bar{o}}-i\bar{\mathbf{y}}_{\bar{o}}\left(B{\mathbf{1}}_{Q}-{\mathbf{b}}\right)^{-1}\bar{\mathbf{y}}^{t}_{\bar{o}}\right)^{-1}\left({\mathbf{A}}_{\bar{o}o}+\bar{\mathbf{y}}_{o}\left(B{\mathbf{1}}_{Q}-{\mathbf{b}}\right)^{-1}\bar{\mathbf{y}}^{t}_{\bar{o}}\right).$

Making use of the Woodbury matrix identity wik  one can further write :

$\displaystyle\left({\mathbf{1}}_{\bar{o}\bar{o}}-i{\mathbf{A}}_{\bar{o}\bar{o}}-i\bar{\mathbf{y}}_{\bar{o}}\left(B{\mathbf{1}}_{Q}-{\mathbf{b}}\right)^{-1}\bar{\mathbf{y}}^{t}_{\bar{o}}\right)^{-1}=\left({\mathbf{1}}_{\bar{o}\bar{o}}-i{\mathbf{A}}_{\bar{o}\bar{o}}\right)^{-1}$
(37) $\displaystyle+$ $\displaystyle
i\left({\mathbf{1}}_{\bar{o}\bar{o}}-i{\mathbf{A}}_{\bar{o}\bar{o}}\right)^{-1}\bar{\mathbf{y}}_{\bar{o}}\left(B{\mathbf{1}}_{Q}-{\mathbf{b}}-i\bar{\mathbf{y}}^{t}_{\bar{o}}\left({\mathbf{1}}_{\bar{o}\bar{o}}-{\mathbf{A}}_{\bar{o}\bar{o}}\right)^{-1}\bar{\mathbf{y}}_{\bar{o}}\right)^{-1}\bar{\mathbf{y}}^{t}_{\bar{o}}\left({\mathbf{1}}_{\bar{o}\bar{o}}-i{\mathbf{A}}_{\bar{o}\bar{o}}\right)^{-1}.$

Substituting (37) into (36), defining a new predissociation matrix
${\mathbf{u}}_{\bar{o}}=\left({\mathbf{1}}_{\bar{o}\bar{o}}-i{\mathbf{A}}_{\bar{o}\bar{o}}\right)^{-1}\bar{\mathbf{y}}_{\bar{o}}$
and the complex resonance position
${\mathbf{B}}_{c}={\mathbf{b}}+i\bar{\mathbf{y}}^{t}_{\bar{o}}\left({\mathbf{1}}_{\bar{o}\bar{o}}-{\mathbf{A}}_{\bar{o}\bar{o}}\right)^{-1}\bar{\mathbf{y}}_{\bar{o}}$,
equation (36) can be recast in the more compact form

$\displaystyle{\mathbf{a}}(B)$ $\displaystyle=$
$\displaystyle-{\mathbf{A}}_{oo}-\bar{\mathbf{y}}_{o}\left(B{\mathbf{1}}_{Q}-{\mathbf{b}}\right)^{-1}\bar{\mathbf{y}}^{t}_{o}$
$\displaystyle-$ $\displaystyle
i\left({\mathbf{A}}_{o\bar{o}}+\bar{\mathbf{y}}_{o}(B{\mathbf{1}}_{Q}-{\mathbf{b}})^{-1}\bar{\mathbf{y}}^{t}_{\bar{o}}\right)\left(({\mathbf{1}}_{\bar{o}\bar{o}}-i{\mathbf{A}}_{\bar{o}\bar{o}})^{-1}+i{\mathbf{u}}_{\bar{o}}(B{\mathbf{1}}_{Q}-{\mathbf{B}}_{c})^{-1}{\mathbf{u}}^{t}_{\bar{o}}\right)\left({\mathbf{A}}_{\bar{o}o}+\bar{\mathbf{y}}_{\bar{o}}(B{\mathbf{1}}_{Q}-{\mathbf{b}})^{-1}\bar{\mathbf{y}}^{t}_{o}\right).$

With some additional algebra and the definition
$\bar{\mathbf{u}}_{o}=i{\mathbf{A}}_{o{\bar{o}}}{\mathbf{u}}_{\bar{o}}+\bar{\mathbf{y}}_{o}$
it is possible to put (II.4) in a form that clearly exhibits the divergence of
the complex scattering length at the complex eigenvalues of the matrix
${\mathbf{B}}_{c}$ :

$\displaystyle{\mathbf{a}}(B)$ $\displaystyle=$
$\displaystyle-{\mathbf{A}}_{oo}-i{\mathbf{A}}_{o{\bar{o}}}\left({\mathbf{1}}_{{\bar{o}}{\bar{o}}}-i{\mathbf{A}}_{{\bar{o}}{\bar{o}}}\right)^{-1}{\mathbf{A}}_{{\bar{o}}o}-\bar{\mathbf{u}}_{o}(B{\mathbf{1}}_{Q}-{\mathbf{B}}_{c})^{-1}\bar{\mathbf{u}}_{o}^{t}.$
(39)

It should be noted that the matrix function ${\mathbf{B}}_{c}$ and therefore
its eigenvalues $B_{c,i}$ depend on the magnetic field $B$ through the
dependence of ${\mathbf{B}}_{c}$ on ${\mathbf{A}}_{{\bar{o}}{\bar{o}}}$, which
in general vary with $B$. In order to carry out the pole expansion one can
proceed in a similar way as it has been done for the $K$-matrix. First one
writes

$(B{\mathbf{1}}_{Q}-{\mathbf{B}}_{c})^{-1}=\dfrac{\text{adj}{~{}(B{\mathbf{1}}_{Q}-{\mathbf{B}}_{c})}}{\det{(B{\mathbf{1}}_{Q}-{\mathbf{B}}_{c})}}\equiv\dfrac{{\mathbf{C}}(B)}{\det{(B{\mathbf{1}}_{Q}-{\mathbf{B}}_{c})}}$
(40)

where the determinant can be expressed as the product of the eigenvalues of
its matrix argument

$\det{(B{\mathbf{1}}_{Q}-{\mathbf{B}}_{c})}={\prod}_{\alpha}(B-B_{c,\alpha}(B)).$
(41)

Hence, it is clear that poles occur at complex magnetic field locations
$B=b_{c,\alpha}$ such that $b_{c,\alpha}-B_{c,\alpha}(b_{c,\alpha})=0$. The
complex poles $b_{c,\alpha}$ in (39) can equivalently be defined as the
solutions of the nonlinear eigenvalue problem
${\mathbf{B}}_{c}(B)|\Omega_{\alpha}(B)\rangle=B|\Omega_{\alpha}(B)\rangle$.
In practice the dependence on $B$ is usually weak and such nonlinear problem
can be solved in a few iterations with excellent starting point the
eigenvalues and eigenvectors of the linear problem with constant
$\mathbf{A}=\mathbf{A}(b_{\alpha})$.

Using regularity of the adjugate and omitting for notational simplicity the
common energy dependence, the residue at the pole can be extracted as the
limit

$\displaystyle{\mathbf{p}}_{\alpha}$ $\displaystyle=$
$\displaystyle\mathop{\rm Res}(\mathbf{a},b_{c,\alpha})={\lim}_{B\to
b_{c,\alpha}}(B-b_{c,\alpha}){\mathbf{a}}(B)=-\dfrac{\bar{\mathbf{u}}_{o}(b_{c,\alpha})\text{adj}{\left(b_{c,\alpha}{\mathbf{1}}_{Q}-{\mathbf{B}}_{c}(b_{c,\alpha})\right)}\bar{\mathbf{u}}_{o}^{t}(b_{c,\alpha})}{(1-B_{c,\alpha}^{\prime}(b_{c,\alpha}))\prod_{\beta\neq\alpha}(b_{c,\alpha}-B_{c,\beta}(b_{c,\alpha}))}$
(42) $\displaystyle=$
$\displaystyle-\dfrac{\bar{\mathbf{u}}_{o}(b_{c,\alpha}){\mathbf{C}}(b_{c,\alpha})\bar{\mathbf{u}}_{o}^{t}(b_{c,\alpha})}{(1-B_{c,\alpha}^{\prime}(b_{c,\alpha}))D_{\alpha}(b_{c,\alpha})},$

where we have defined
$D_{\alpha}(B)=\prod_{\beta\neq\alpha}(B-B_{c,\beta}(B))$.

The local complex background at the pole $b_{c,\alpha}$ is defined by
subtracting the contribution of pole $\alpha$ from the full complex scattering
length :

${\mathbf{a}}_{{\rm bgl},\alpha}={\lim}_{B\to
b_{c,\alpha}}\left[{\mathbf{a}}(B)-\dfrac{{\mathbf{p}}_{\alpha}}{B-b_{c,\alpha}}\right]$
(43)

Substituting (42) for the pole strength, bringing to common denominator and
applying l’Hôpital’s rule twice, one obtains

$\displaystyle{\mathbf{a}}_{{\rm bgl},\alpha}$ $\displaystyle=$
$\displaystyle-{\mathbf{A}}_{oo}(b_{c,\alpha})-i{\mathbf{A}}_{o{\bar{o}}}(b_{c,\alpha})\left({\mathbf{1}}_{{\bar{o}}{\bar{o}}}-i{\mathbf{A}}_{{\bar{o}}{\bar{o}}}(b_{c,\alpha})\right)^{-1}{\mathbf{A}}_{{\bar{o}}o}(b_{c,\alpha})-\dfrac{B_{c,\alpha}^{\prime\prime}(b_{c,\alpha})\bar{\mathbf{u}}_{o}(b_{c,\alpha}){\mathbf{C}}(b_{c,\alpha})\bar{\mathbf{u}}_{o}^{t}(b_{c,\alpha})}{2(1-B_{c,\alpha}^{\prime}(b_{c,\alpha}))^{2}D_{\alpha}(b_{c,\alpha})}$
(44) $\displaystyle-$
$\displaystyle\dfrac{(\bar{\mathbf{u}}_{o}{\mathbf{C}}\bar{\mathbf{u}}_{o}^{t})^{\prime}(b_{c,\alpha})}{(1-B_{c,\alpha}^{\prime}(b_{c,\alpha}))D_{\alpha}(b_{c,\alpha})}+\dfrac{D_{\alpha}^{\prime}(b_{c,\alpha})\bar{\mathbf{u}}_{o}(b_{c,\alpha}){\mathbf{C}}(b_{c,\alpha})\bar{\mathbf{u}}_{o}^{t}(b_{c,\alpha})}{(1-B_{c,\alpha}^{\prime}(b_{c,\alpha}))D_{\alpha}^{2}(b_{c,\alpha})}.$

This quantity contains contributions from all poles but $\alpha$-th one. The
global complex background at the pole $b_{c,\alpha}$ is defined by subtracting
all other pole contributions :

${\mathbf{a}}_{{\rm bg},\alpha}={\mathbf{a}}_{{\rm
bgl},\alpha}-\sum_{\beta\neq\alpha}\dfrac{{\mathbf{p}}_{\beta}}{b_{c,\alpha}-b_{c,\beta}}$
(45)

resulting in

$\displaystyle{\mathbf{a}}_{{\rm bg},\alpha}$ $\displaystyle=$
$\displaystyle-{\mathbf{A}}_{oo}(b_{c,\alpha})-i{\mathbf{A}}_{o{\bar{o}}}(b_{c,\alpha})\left({\mathbf{1}}_{{\bar{o}}{\bar{o}}}-i{\mathbf{A}}_{{\bar{o}}{\bar{o}}}(b_{c,\alpha})\right)^{-1}{\mathbf{A}}_{{\bar{o}}o}(b_{c,\alpha})-\dfrac{B_{c,\alpha}^{\prime\prime}(b_{c,\alpha})\bar{\mathbf{u}}_{o}(b_{c,\alpha}){\mathbf{C}}(b_{c,\alpha})\bar{\mathbf{u}}_{o}^{t}(b_{c,\alpha})}{2(1-B_{c,\alpha}^{\prime}(b_{c,\alpha}))^{2}D_{\alpha}(b_{c,\alpha})}$
(46) $\displaystyle-$
$\displaystyle\dfrac{(\bar{\mathbf{u}}_{o}{\mathbf{C}}\bar{\mathbf{u}}_{o}^{t})^{\prime}(b_{c,\alpha})}{(1-B_{c,\alpha}^{\prime}(b_{c,\alpha}))D_{\alpha}(b_{c,\alpha})}+\dfrac{D_{\alpha}^{\prime}(b_{c,\alpha})\bar{\mathbf{u}}_{o}(b_{c,\alpha}){\mathbf{C}}(b_{c,\alpha})\bar{\mathbf{u}}_{o}^{t}(b_{c,\alpha})}{(1-B_{c,\alpha}^{\prime}(b_{c,\alpha}))D_{\alpha}^{2}(b_{c,\alpha})}-{\sum}_{\beta\neq\alpha}\dfrac{\mathbf{p}_{\beta}}{b_{c,\alpha}-b_{c,\beta}}$

Several quantities need to be computed in order to evaluate the rhs of this
equation. The derivatives $B_{c,\alpha}^{\prime}$ and
$B_{c,\alpha}^{\prime\prime}$ are obtained by first writing, for $B$ near
$b_{c,\alpha}$,
${\mathbf{B}}_{c}(B)\approx{\mathbf{B}}_{c}(b_{c,\alpha})+{\mathbf{B}}_{c}^{\prime}(b_{c,\alpha})\eta$
with $\eta=(B-b_{c,\alpha})$ and the derivative ${\mathbf{B}}_{c}^{\prime}$
promptly computed from the definition of ${\mathbf{B}}_{c}$ in terms of
${\mathbf{A}}_{{\bar{o}}{\bar{o}}}$ and of
$\mathbf{A}^{\prime}_{{\bar{o}}{\bar{o}}}$. Next, we apply perturbation theory
with small parameter $\eta$ to the eigenvalue problem
${\mathbf{B}}_{c}(B)|\Omega_{\alpha}(B)\rangle=B_{c,\alpha}(B)|\Omega_{\alpha}(B)\rangle$
as

$\left({\mathbf{B}}_{c}(b_{c,\alpha})+{\mathbf{B}}_{c}^{\prime}(b_{c,\alpha})\eta\right)(|\Omega_{\alpha}^{(0)}+|\Omega_{\alpha}^{(1)}\rangle\eta+\cdots)=(b_{c,\alpha}+b^{(1)}_{c,\alpha}\eta+b^{(2)}_{c,\alpha}\eta^{2}+\cdots)(|\Omega_{\alpha}^{(0)}\rangle+|\Omega^{(1)}_{\alpha}\rangle\eta+\cdots),$
(47)

where the superscript (k) denotes perturbation terms of order $k$ on
eigenvalues and eigenfunctions. Equating $\eta$ powers order by order, one can
compute terms up to perturbative order $k=2$ in the eigenvalues, and identify
$B_{c,\alpha}^{\prime}(b_{c,\alpha})=b^{(1)}_{c,\alpha}$ and
$2B_{c,\alpha}^{\prime\prime}(b_{c,\alpha})=b^{(2)}_{c,\alpha}$.

The only difference with the standard procedure is that matrix
${\mathbf{B}}_{c}$ is complex symmetric rather than self-adjoint such that a
generalization of perturbation theory to non-Hermitian operators must be
employed in order to calculate $b^{(1)}_{c,\alpha}$ and $b^{(2)}_{c,\alpha}$.
It should be stressed that the only approximation involved in this procedure
is the neglect of the second derivative in the Taylor expansion of
${\mathbf{B}}_{c}(B)$ near $b_{c,\alpha}$. This is justified since this term
is usually small and it could in principle be included if needed, even if at
the expenses of significant added computational complexity to obtain
$\mathbf{A}^{\prime\prime}_{{\bar{o}}{\bar{o}}}$. Further, derivatives of the
adjugate $\mathbf{C}$ are evaluated as explained below equation (21). The
terms $\bar{\mathbf{u}}_{o}^{\prime}$ are obtained based on the definition of
$\bar{\mathbf{u}}_{o}$ and from the knowledge of $\mathbf{A}$ and
${\mathbf{A}}^{\prime}$. Also note the useful formula

$\dfrac{D^{\prime}_{\alpha}}{D_{\alpha}}=\log^{\prime}(D_{\alpha})={\sum}_{\beta\neq\alpha}\log^{\prime}(B-B_{c,\beta})={\sum}_{\beta\neq\alpha}\dfrac{(1-B^{\prime}_{c,\beta})}{B-B_{c,\beta}}$
(48)

allowing the ratio in last but one term on the rhs of (46) to be easily
determined.

Once the pole strength ${\mathbf{p}}_{\alpha}$ known for each resonance, one
can cast the energy-dependent complex scattering length (39) in the intuitive
form

${\mathbf{a}}(B)={\mathbf{a}}_{\rm
bg}(B)+\sum_{\alpha}\dfrac{{\mathbf{p}}_{\alpha}}{B-b_{c,\alpha}}$ (49)

with

${\mathbf{a}}_{\rm
bg}(B)=-{\mathbf{A}}_{oo}-i{\mathbf{A}}_{o{\bar{o}}}\left({\mathbf{1}}_{{\bar{o}}{\bar{o}}}-i{\mathbf{A}}_{{\bar{o}}{\bar{o}}}\right)^{-1}{\mathbf{A}}_{{\bar{o}}o}-\bar{\mathbf{u}}_{o}(B{\mathbf{1}}_{Q}-{\mathbf{B}}_{c})^{-1}\bar{\mathbf{u}}_{o}^{t}-\sum_{\alpha}\dfrac{{\mathbf{p}}_{\alpha}}{B-b_{c,\alpha}}$
(50)

a quantity by construction regular at the poles. If $\mathbf{A}$ is strictly
constant, it is clear that last two terms exactly cancel and the background
reduces to

${\mathbf{a}}_{\rm
bg}(B)=-{\mathbf{A}}_{oo}-i{\mathbf{A}}_{o{\bar{o}}}\left({\mathbf{1}}_{{\bar{o}}{\bar{o}}}-i{\mathbf{A}}_{{\bar{o}}{\bar{o}}}\right)^{-1}{\mathbf{A}}_{{\bar{o}}o}.$
(51)

Generalizing the usual definition from the single-channel case Tiesinga  _et
al._ (1993), we define the partial magnetic widths associated to pole $\alpha$
via the equation

${\mathbf{\Delta}}_{\alpha}=-{\mathbf{a}}_{\rm
bg}^{-1/2}(b_{c,\alpha})~{}{\mathbf{p}}_{\alpha}~{}{\mathbf{a}}_{\rm
bg}^{-1/2}(b_{c,\alpha}).$ (52)

In the absence of inelastic channels, the matrix ${\Delta}_{\alpha}$ is real,
its imaginary part arises from inelastic processes.

As a final note, we have undertaken the approach to collect all degenerate
channels in the same open-channel manifold, which led us naturally to
introduce matrices of scattering lengths and of magnetic widths. This
description is most useful in the presence of degenerate partial waves or in
equilibrium situations where all degenerate internal states are equally
populated. In other experiments one is more interested in polarizing the atoms
in a single quantum state and in either monitoring the evolution of the
population or studying the collective behavior of the gaz in such specific
state. In this case, it will be more useful to identify the open channel
manifold with the unique state of interest, and to treat super-elastic
channels on the same footing as the inelastic ones. All derivation above
holds, the main difference being that the complex scattering length will not
be real even in the absence of inelastic channels, since the initial quantum
state will in general decay to degenerate channels.

Figure 1:  The real and the imaginary part of the complex scattering length
for collisions of two 39K atoms in hyperfine states $|1~{}1\rangle$ and
$|1~{}-1\rangle$. Full curves represent results of a close-coupled numerical
calculation, the dots the pole approximation of equation (49) evaluated on a
regular magnetic-field grid.

##  III Comparison with other approaches

The simplest case of an isolated elastic resonance in an elastic channel is
trivially recovered from equation (49) by dropping the $\alpha$ index, putting
$p=-a_{\rm bg}\Delta$ and realizing that the resonance center $b_{c}$ is a
real quantity.

For the case of several overlapping resonances in an elastic channel,
simplicity of equation (49) may seem at first sight surprising as compared to
equation (26) of Jachymski and Julienne Jachymski and Juilenne (2013). In
order to show that the present approach reduces to the quantum defect analysis
presented in the latter work, we give a brief alternative derivation of their
result starting from the present standard scattering theory. To this aim, let
us specialize (5) to the case of a single open channel and consider by the
sake of simplicity the low-energy limit :

${\mathbf{K}}=\cos({\xi_{\rm bg}})\left[\tan({\xi}_{\rm bg})+{\mathbf{K}}_{\rm
res}\right]\left[{\bm{1}}_{P}-\tan({\xi}_{\rm bg}){\mathbf{K}}_{\rm
res}\right]^{-1}\cos({\xi_{\rm bg}})^{-1}\underset{E\to 0}{\to}\tan({\xi}_{\rm
bg})+{\mathbf{K}}_{\rm res}$ (53)

In the approach of Jachymski and Julienne, one first factors the energy
dependence out of the predissociation amplitude ${\mathbf{y}}$ using the
energy-dependent quantum defect amplitude $C$ as
${\mathbf{y}}=C^{-1}\bar{\mathbf{y}}$. Then one recognizes that the shift
matrix in the denominator of ${\mathbf{K}}_{\rm res}$ contains an energy-
dependent part proportional to $\bar{\mathbf{y}}\bar{\mathbf{y}}^{t}$ via the
quantum defect function $\tan(\lambda)$ Bohn and Julienne (1999). This leads
to :

${\mathbf{K}}=\tan({\xi}_{\rm
bg})+C^{-2}\bar{\mathbf{y}}^{t}\left({\mathbf{D}-\tan(\lambda)\bar{\mathbf{y}}\bar{\mathbf{y}}^{t}}\right)^{-1}\bar{\mathbf{y}}\quad,\quad
E\to 0$ (54)

where $\mathbf{D}$ contains the Zeeman shift
$-{\bm{\mu}}\left(B{\mathbf{1}}_{Q}-{\mathbf{B}}_{r}\right)$ and all energy-
independent couplings between channels in $Q$ space. Making use of Woodbory’s
matrix identity to express the inverse, equation (54) is seen to be equivalent
to

${\mathbf{K}}=\tan({\xi}_{\rm
bg})+C^{-2}\frac{\bar{\mathbf{y}}^{t}{\mathbf{D}}^{-1}\bar{\mathbf{y}}}{1-\tan(\lambda)\bar{\mathbf{y}}^{t}{\mathbf{D}}^{-1}\bar{\mathbf{y}}}\quad,\quad
E\to 0$ (55)

Assuming absence of coupling between molecular states near threshold, matrix
$\mathbf{D}$ takes the simple diagonal form
${\mathbf{D}}=-{\bm{\mu}}\left(B{\mathbf{1}}_{Q}-{\mathbf{B}}_{r}\right)$.
Replacing in (55) leads after simple manipulations to equation (26) of
Jachymski and Juilenne (2013).

Note that in the present work one first diagonalizes the denominator of
${\mathbf{K}}_{\rm res}$ at each value of the magnetic field. Therefore, our
equation (49) can be considered in a sense the adiabatic version of equation
(26) of Jachymski and Julienne. Note however that equation (49) in the present
work has a somewhat broader domain of applicability, since it does not rely on
the assumption that the molecular states are uncoupled near threshold.

We now consider an inelastic isolated resonance. With reference to the
notation of equation (1) in van Kempen  _et al._ (2002) their lineshsape is
recovered by first writing our complex pole $p$ in the form $p=-\Delta a_{\rm
bg}$, then by expressing the complex width $\Delta$ in polar form as
$e^{2i\Phi_{R}}\Delta_{\rm el}$. Finally, one identifies $a_{\rm bg}$ with
$a_{\infty}$ and the complex magnetic field location $b_{c}$ in the
denominator of (49) with the term $B_{0}-i\dfrac{\Delta_{\rm inel}}{2}$ in van
Kempen  _et al._ (2002).

Barring from the choice of the notation, it is also easy to verify the
equivalence of the lineshapes in Hutson and Soldán (2007) and van Kempen  _et
al._ (2002) and thus in the present work.

##  IV Applications

###  IV.1 Overlapping inelastic $s$-wave resonances

As a first application, let us consider collisions of ultracold 39K atoms, a
system whose scattering properties have been determined very precisely first
based on photoassociation then on magnetic Feshbach resonance spectroscopy
D’Errico  _et al._ (2007). We consider here the case of $s$-wave collisions
for atoms in Zeeman sublevels correlating in zero magnetic field with
hyperfine quantum numbers $(F=1,~{}M=1)$ and $(F=1,~{}M=-1)$, respectively.
This combination is not stable under spin-exchange collisions, since it can
decay to a pair of $(F=1,~{}M=0)$ atoms.

Figure 2:  Trajectory in the complex plane of the complex scattering length
when the magnetic field spans the features of figure 1. The inset shows a
smaller loop described near the origin by $a(B)$ as the magnetic field varies
across the weaker of the two features. Note the largely different scales on
the axes of the main figure and the inset.

The complex scattering length shows two features, which have been reported
before Lysebo and Veseth (2008). Based on the much stronger variation of the
real part of $a$ near 475 G than near 500 G, the authors of Lysebo and Veseth
(2008) concluded that only one feature is of resonant origin. A closer
analysis of figure 1 shows however that the real part of $a$ does not diverge
near either feature, as expected on a general ground in the presence of
inelastic processes. While the feature at higher field is indeed significantly
more pronounced, one may still suspect the presence of two quasi-bound states
underlying the observed structure.

Figure 3:  The elements of the normalized matrix $\bar{\mathbf{K}}$, namely
$\bar{\mathbf{K}}_{11}$ (full line), $\bar{\mathbf{K}}_{12}$ (dashed line) and
$\bar{\mathbf{K}}_{22}$ (dash-dotted line) for for collisions of 39K atoms in
hyperfine states $|1~{}1\rangle$ and $|1~{}-1\rangle$. Vertical dotted lines
indicate the location of two neighboring poles. Indices $1$ and $2$ denote the
two-body collisions channels $|1~{}0\rangle\otimes|10\rangle$ and
$|1~{}1\rangle\otimes|1~{}-1\rangle$, respectively.

A complimentary picture is provided by the trajectory of $a(B)$ in the complex
scattering-length plane as a function of the real parameter $B$; See figure 2.
It it easy to verify that in the ideal case of an isolated resonance with
strictly constant background such trajectory is a circle Hutson and Soldán
(2007). One can actually observe in the figure one large circle associated
with the strong resonance, whereas closer inspection near the origin (see
inset) highlights the presence of a second deformed loop described as the
magnetic field crosses the weak feature, suggesting again the presence of a
resonance pair. One can also observe that the global trajectory does not
exactly close, since the values of the background scattering length away from
resonance are not identical on the left and right resonance sides.

Conclusive evidence about the nature of this feature is obtained by inspecting
the elements of the full $K$-matrix depicted in figure 3. In agreement with
equation (11), all elements of the reactance matrix diverge on resonance at
the same location. Two poles are clearly observed, confirming that the
relevant $Q$ subspace contains two bare bound states coupled to the incoming
channel. Following the numerical procedure outlined in section II, it is now
easy to extract well defined lineshape parameters without resorting to a fit.
Resulting numerical values for the complex background scattering length, field
location and magnetic width can be found in table 1.

It is also interesting to compare graphically the exact numerical result with
the pole approximation (49) represented by the dots in figure 1. Note from the
numerical values in the table that the background is nearly constant and can
be represented locally by an interpolating linear function. As expected, the
agreement with the non trivial behavior of the full numerical calculation is
excellent, confirming that all elements of the theory have been computed and
integrated properly.

Table 1:  Complex lineshape parameters (see text for definitions) appearing in
lineshape (49), with powers of tens indicated in brackets. Sample resonance
parameters for ultracold collisions of K atoms in hyperfine states and partial
wave specified in the first column are reported. For comparison, calculations
with and without the spin-spin interaction are reported for $p$-wave
scattering.
Channel | $a_{\rm bg}(a_{0}^{ \ell+1})$ |  $b_{c}$(G) | $p_{\alpha}(a_{0}^{\ell+1}{\rm G})$ |  $\Delta_{\alpha}$(G)
---|---|---|---|---
K($|1~{}1\rangle$)+K($|1~{}-1\rangle$) $s$-wave | $29.1+5.37[-2]i$ | $472.67+9.45i$ | $-639+26.9i$ | $21.9-0.965i$
| $29.3+4.97[-2]i$ | $501.34+1.70[-3]i$ | $-897+5.91i$ | $30.6-0.253i$
K($|1~{}1\rangle$)+K($|1~{}0\rangle$) $p$-wave (with spin-spin) | $-9.28[4]$ | $379.1$ | $-4.00[6]$ | $-43.1$
| $-9.60[4]$ | $400.7$ | $-1.27[6]$ | $-1.32$
| $-1.16[5]$ | $720.3$ | $-3.80[2]$ | $-3.26[-3]$
| $-1.17[5]$ | $757.8$ | $-2.38[3]$ | $-2.03[-2]$
| $-1.17[5]$ | $762.3$ | $-3.20[4]$ | $-0.272$
K($|1~{}1\rangle$)+K($|1~{}0\rangle$) $p$-wave (w/o spin-spin) | $2.81[5]$ | $379.2$ | $-4.01[6]$ | $14.3$
| $2.81[5]$ | $400.6$ | $-1.27[6]$ | $0.451$
| $2.87[5]$ | $758.0$ | $-2.38[3]$ | $8.30[-3]$
| $2.87[5]$ | $762.6$ | $-3.21[4]$ | $0.111$

###  IV.2 Energy-dependent $p$-wave scattering volume

As a further example, let us consider collisions in higher-order partial
waves. Specifically, we consider $p$-wave collisions between a pair of 39K
atoms prepared respectively in states $|1~{}1\rangle$ and $|1~{}0\rangle$ at a
finite collision energy of $E/k_{\rm B}=1\mu{\rm K}$. We consider the
projection $M=2$ of the total axial angular momentum, corresponding to a pair
of atoms with initial axial mechanical angular momentum $m=1$. In this
symmetry block, if $p$ waves only are included, scattering is purely elastic.
The energy-dependent scattering volume is defined according to the standard
definition (32).

Figure 4:  Energy dependent scattering volume for collisions of two 39K atoms
in hyperfine states $|11\rangle$ and $|10\rangle$ (full line). The dashed line
represents the background scattering volume (see text). The inset shows the
detail of a triplet of resonances above 700 G, of which two overlapping.
Figure 5:  Same as figure 4, but with the spin-spin interaction set
artificially to zero.

One can observe in figure 3 two series of overlapping resonances in the range
of magnetic fields below $1000$G, a doublet near 400 G and a triplet near 750
G. The largest resonance at 379 G is sufficiently broad to modify the local
background near 750 G. Visually, the background appears to vary in a non
trivial way with the magnetic field, in particular near the origin. The
present approach has the advantage that no modeling is required in order to
extract the resonance parameters. It is sufficient to locate precisely the
resonances, extract the corresponding pole residues, and finally the local and
global background at the resonance positions. The corresponding parameters for
the present case are reproduced in table 1.

If needed, the field-dependent background scattering length can be extracted
by subtracting the pole contribution from the full scattering volume, i.e.
last term on the rhs of (49) from the lhs. In spite of the fact that such
operation involves taking the difference of large numbers and as such is in
principle numerically flawed very close to the poles, we obtain a smooth
background that, as anticipated, varies relatively fast in a non trivial way
from negative to positive values.

Without performing a detailed analysis that is not the aim of the present
methodological work, it is interesting to stress the effect of the spin-spin
dipolar interaction. To this extent we set artificially to zero the dipolar
coupling constant and perform the same steps as above. In the absence of
anisotropic interactions, not only $M$ but also $M_{f}=M_{a}+M_{b}$ is a good
quantum number. This implies that resonances associated with molecular states
with $M_{f}$ different from the one of the incoming atoms are completely
uncoupled, i.e. have zero width. Indeed, one can remark in figure 4 the
absence of very narrow resonances observed in figure 3. The remaining
relatively broad resonances occur in blocks with given $M_{f}$ and are induced
by the interplay of spin-exchange and spin-spin interactions. In the present
case, spin-exchange is found to largely dominate over dipolar interactions, as
demonstrated by the fact that the pole strength $p$ in the presence and in the
absence of spin-spin couplings (see table) differ by no more less than $1\%$.
Similarly, the spin-spin interaction only weakly (below 300 mG) shifts the
resonance position. On the other side, the background scattering volume (and
a-fortiori the magnetic width), has a significantly stronger variation with
magnetic field in the presence of the dipolar interaction, highlighting the
strong influence of the spin-spin interaction on the incoming wavefunction.

##  V Conclusions

We have presented a formalism allowing resonance parameters to be extracted
from a numerical calculation in a well determined algorithmic way. In our
opinion the present approach fully solves a longstanding difficulty in the
analysis of resonance spectra involving multiple overlapping features or a non
trivial behavior of background scattering. Our results also provide a simple
and intuitive lineshape for inelastic overlapping resonances that has been
demonstrated starting from first principles. The machinery presented in this
work can be incorporated in existing computer packages, provided that the
underlying algorithm to solve the Schrödinger equation allows the wavefunction
derivatives to be easily evaluated. As a complimentary step to the present
paper that focuses on the lineshape parameters, it could be interesting to
analyze the wavefunction at the poles of the reactance matrix in order to gain
information on the nature of the resonant states.

###### Acknowledgements.

Useful exchanges with P. S. Julienne and K. Jachymski are gratefully
acknowledged.

##  Appendix A Proof of various equations

Proof of equation (9) goes as follows. Substituting (8) into (5), using
Woodbury matrix identity and the definition of $\hat{\bm{\mathcal{B}}}$, one
gets

$\displaystyle{\mathbf{K}}$ $\displaystyle=$ $\displaystyle\cos({\bm{\xi}_{\rm
bg}})\left[\tan({\bm{\xi}}_{\rm
bg})+\tilde{\mathbf{y}}\left(B{\bm{1}}_{Q}-{\bm{\mathcal{B}}}\right)^{-1}\tilde{\mathbf{y}}^{t}\right]\left[{\bm{1}}_{P}-\tan({\bm{\xi}}_{\rm
bg})\tilde{\mathbf{y}}\left(B{\bm{1}}_{Q}-{\bm{\mathcal{B}}}\right)^{-1}\tilde{\mathbf{y}}^{t}\right]^{-1}\cos({\bm{\xi}_{\rm
bg}})^{-1}$ (56) $\displaystyle=$ $\displaystyle\cos({\bm{\xi}_{\rm
bg}})\left[\tan({\bm{\xi}}_{\rm
bg})+\tilde{\mathbf{y}}\left(B{\bm{1}}_{Q}-{\bm{\mathcal{B}}}\right)^{-1}\tilde{\mathbf{y}}^{t}\right]\left[{\bm{1}}_{P}+\tan({\bm{\xi}}_{\rm
bg})\tilde{\mathbf{y}}\left(B{\bm{1}}_{Q}-{\bm{\mathcal{B}}}-\tilde{\mathbf{y}}^{t}\tan({\bm{\xi}}_{\rm
bg})\tilde{\mathbf{y}}\right)^{-1}\tilde{\mathbf{y}}^{t}\right]\cos({\bm{\xi}_{\rm
bg}})^{-1}$ $\displaystyle=$ $\displaystyle\cos({\bm{\xi}_{\rm
bg}})\left[\tan({\bm{\xi}}_{\rm
bg})+\tilde{\mathbf{y}}\left(B{\bm{1}}_{Q}-{\bm{\mathcal{B}}}\right)^{-1}\tilde{\mathbf{y}}^{t}\right]\left[{\bm{1}}_{P}+\tan({\bm{\xi}}_{\rm
bg})\tilde{\mathbf{y}}\left(B{\bm{1}}_{Q}-\hat{\bm{\mathcal{B}}}\right)^{-1}\tilde{\mathbf{y}}^{t}\right]\cos({\bm{\xi}_{\rm
bg}})^{-1}$

Carrying out the inner product, replacing the term
$\tilde{\mathbf{y}}^{t}\tan({\bm{\xi}}_{\rm bg})\tilde{\mathbf{y}}$ by
$(\hat{\bm{\mathcal{B}}}-{\bm{\mathcal{B}}})$ in the numerator and factoring
by grouping

$\displaystyle{\mathbf{K}}$ $\displaystyle=$ $\displaystyle\cos({\bm{\xi}_{\rm
bg}})\left\\{\tan({\bm{\xi}}_{\rm
bg})+\tilde{\mathbf{y}}\left[(B{\bm{1}}_{Q}-{\bm{\mathcal{B}}})^{-1}+(B{\bm{1}}_{Q}-{\bm{\mathcal{B}}})^{-1}(\hat{\bm{\mathcal{B}}}-{\bm{\mathcal{B}}})(B{\bm{1}}_{Q}-\hat{\bm{\mathcal{B}}})^{-1}\right]\tilde{\mathbf{y}}^{t}\right\\}\cos({\bm{\xi}_{\rm
bg}})^{-1}$ (57) $\displaystyle+$ $\displaystyle\cos({\bm{\xi}_{\rm
bg}})\tan^{2}({\bm{\xi}_{\rm
bg}})\tilde{\mathbf{y}}(B{\bm{1}}_{Q}-\hat{\bm{\mathcal{B}}})^{-1}\tilde{\mathbf{y}}^{t}.$

Finally noticing that

$\left[(B{\bm{1}}_{Q}-{\bm{\mathcal{B}}})^{-1}+(B{\bm{1}}_{Q}-{\bm{\mathcal{B}}})^{-1}(\hat{\bm{\mathcal{B}}}-{\bm{\mathcal{B}}})(B{\bm{1}}_{Q}-\hat{\bm{\mathcal{B}}})^{-1}\right]=(B{\bm{1}}_{Q}-{\bm{\mathcal{B}}})^{-1}\left[B{\bm{1}}_{Q}-\hat{\bm{\mathcal{B}}}+\hat{\bm{\mathcal{B}}}-{\bm{\mathcal{B}}}\right](B{\bm{1}}_{Q}-\hat{\bm{\mathcal{B}}})^{-1}=B{\bm{1}}_{Q}-\hat{\bm{\mathcal{B}}},$
(58)

substituting in (57) and using simple trigonometric identities, equation (9)
easily follows.

To prove (39), we start by expanding (II.4) as

$\displaystyle{\mathbf{a}}(B)$ $\displaystyle=$
$\displaystyle-{\mathbf{A}}_{oo}-i{\mathbf{A}}_{o\bar{o}}({\mathbf{1}}_{\bar{o}\bar{o}}-i{\mathbf{A}}_{\bar{o}\bar{o}})^{-1}{\mathbf{A}}_{\bar{o}o}+{\mathbf{A}}_{o\bar{o}}{\mathbf{u}}_{\bar{o}}(B{\mathbf{1}}_{Q}-{\mathbf{B}}_{c})^{-1}{\mathbf{u}}^{t}_{\bar{o}}{\mathbf{A}}_{{\bar{o}}o}-\bar{\mathbf{y}}_{o}(B{\mathbf{1}}_{Q}-{\mathbf{b}})^{-1}\bar{\mathbf{y}}^{t}_{\bar{o}}$
(59) $\displaystyle+$
$\displaystyle\bar{\mathbf{y}}_{o}(B{\mathbf{1}}_{Q}-{\mathbf{b}})^{-1}\bar{\mathbf{y}}_{\bar{o}}^{t}{\mathbf{u}}_{\bar{o}}(B{\mathbf{1}}_{Q}-{\mathbf{B}}_{c})^{-1}{\mathbf{u}}_{\bar{o}}^{t}\bar{\mathbf{y}}_{\bar{o}}(B{\mathbf{1}}_{Q}-{\mathbf{b}})^{-1}\bar{\mathbf{y}}_{o}^{t}-i\bar{\mathbf{y}}_{o}(B{\mathbf{1}}_{Q}-{\mathbf{b}})^{-1}\bar{\mathbf{y}}_{\bar{o}}^{t}({\mathbf{1}}_{\bar{o}\bar{o}}-i{\mathbf{A}}_{\bar{o}\bar{o}})^{-1}\bar{\mathbf{y}}_{\bar{o}}(B{\mathbf{1}}_{Q}-{\mathbf{b}})^{-1}\bar{\mathbf{y}}_{o}^{t}$
$\displaystyle-$ $\displaystyle
i{\mathbf{A}}_{o\bar{o}}({\mathbf{1}}_{\bar{o}\bar{o}}-i{\mathbf{A}}_{\bar{o}\bar{o}})^{-1}\bar{\mathbf{y}}_{\bar{o}}(B{\mathbf{1}}_{Q}-{\mathbf{b}})^{-1}\bar{\mathbf{y}}_{o}^{t}+{\mathbf{A}}_{o\bar{o}}{\mathbf{u}}_{\bar{o}}(B{\mathbf{1}}_{Q}-{\mathbf{B}}_{c})^{-1}{\mathbf{u}}^{t}_{\bar{o}}\bar{\mathbf{y}}_{\bar{o}}(B{\mathbf{1}}_{Q}-{\mathbf{b}})^{-1}\bar{\mathbf{y}}^{t}_{o}$
$\displaystyle-$ $\displaystyle
i\bar{\mathbf{y}}_{o}(B{\mathbf{1}}_{Q}-{\mathbf{b}})^{-1}\bar{\mathbf{y}}_{\bar{o}}^{t}({\mathbf{1}}_{\bar{o}\bar{o}}-i{\mathbf{A}}_{\bar{o}\bar{o}})^{-1}{\mathbf{A}}_{\bar{o}o}+{\mathbf{A}}_{o\bar{o}}{\mathbf{u}}_{\bar{o}}(B{\mathbf{1}}_{Q}-{\mathbf{b}})^{-1}{\mathbf{u}}_{\bar{o}}^{t}\bar{\mathbf{y}}_{\bar{o}}(B{\mathbf{1}}_{Q}-{\mathbf{B}}_{c})^{-1}\bar{\mathbf{y}}_{o}^{t}.$

Let us now transform second line. Using definitions of ${\mathbf{B}}_{c}$ and
${\mathbf{u}}_{\bar{o}}$, one can replace the terms
$\bar{\mathbf{y}}_{\bar{o}}^{t}{\mathbf{u}}_{\bar{o}}$,
${\mathbf{u}}_{\bar{o}}^{t}\bar{\mathbf{y}}_{\bar{o}}$, and
$\bar{\mathbf{y}}_{\bar{o}}^{t}({\mathbf{1}}_{\bar{o}\bar{o}}-i{\mathbf{A}}_{\bar{o}\bar{o}})^{-1}\bar{\mathbf{y}}_{\bar{o}}$
by $({\mathbf{B}}_{c}-{\mathbf{b}})/i$. Regrouping and simplifying, second
line transforms into :

$\displaystyle-$
$\displaystyle\bar{\mathbf{y}}_{o}(B{\mathbf{1}}_{Q}-{\mathbf{b}})^{-1}\left[({\mathbf{B}}_{c}-{\mathbf{b}})(B{\mathbf{1}}_{Q}-{\mathbf{B}}_{c})^{-1}({\mathbf{B}}_{c}-{\mathbf{b}})+({\mathbf{B}}_{c}-{\mathbf{b}})\right](B{\mathbf{1}}_{Q}-{\mathbf{b}})^{-1}\bar{\mathbf{y}}^{t}_{o}$
$\displaystyle=$ $\displaystyle-$
$\displaystyle\bar{\mathbf{y}}_{o}(B{\mathbf{1}}_{Q}-{\mathbf{b}})^{-1}(B{\mathbf{1}}_{Q}-{\mathbf{B}}_{c})^{-1}(B{\mathbf{1}}_{Q}-{\mathbf{b}})^{-1}\bar{\mathbf{y}}^{t}_{o}.$
(60)

Adding last term in first line of (59) to result (A) yields after simple
matrix algebra similar to the manipulations in (58)

$\displaystyle-$
$\displaystyle\bar{\mathbf{y}}_{o}(B{\mathbf{1}}_{Q}-{\mathbf{b}})^{-1}(B{\mathbf{1}}_{Q}-{\mathbf{B}}_{c})^{-1}(B{\mathbf{1}}_{Q}-{\mathbf{b}})^{-1}\bar{\mathbf{y}}^{t}_{o}-\bar{\mathbf{y}}_{o}\left(B{\mathbf{1}}_{Q}-{\mathbf{b}}\right)^{-1}\bar{\mathbf{y}}^{t}_{o}$
$\displaystyle=$ $\displaystyle-$
$\displaystyle\bar{\mathbf{y}}_{o}(B{\mathbf{1}}_{Q}-{\mathbf{B}}_{c})^{-1}\bar{\mathbf{y}}^{t}_{o}$
(61)

Along the same lines, one can simplify third and last line in (59). For
instance, for the third line one has

$\displaystyle-$ $\displaystyle
i{\mathbf{A}}_{o\bar{o}}({\mathbf{1}}_{\bar{o}\bar{o}}-i{\mathbf{A}}_{\bar{o}\bar{o}})^{-1}\bar{\mathbf{y}}_{\bar{o}}(B{\mathbf{1}}_{Q}-{\mathbf{b}})^{-1}\bar{\mathbf{y}}_{o}^{t}+{\mathbf{A}}_{o\bar{o}}{\mathbf{u}}_{\bar{o}}(B{\mathbf{1}}_{Q}-{\mathbf{B}}_{c})^{-1}{\mathbf{u}}^{t}_{\bar{o}}\bar{\mathbf{y}}_{\bar{o}}(B{\mathbf{1}}_{Q}-{\mathbf{b}})^{-1}\bar{\mathbf{y}}^{t}_{o}$
$\displaystyle=$ $\displaystyle-$ $\displaystyle
i{\mathbf{A}}_{o\bar{o}}{\mathbf{u}}_{\bar{o}}\left[(B{\mathbf{1}}_{Q}-{\mathbf{b}})^{-1}+(B{\mathbf{1}}_{Q}-{\mathbf{B}}_{c})^{-1}({\mathbf{B}}_{c}-{\mathbf{b}})(B{\mathbf{1}}_{Q}-{\mathbf{b}})^{-1}\right]\bar{\mathbf{y}}^{t}_{o}$
$\displaystyle=$ $\displaystyle-$ $\displaystyle
i{\mathbf{A}}_{o\bar{o}}{\mathbf{u}}_{\bar{o}}(B{\mathbf{1}}_{Q}-{\mathbf{B}}_{c})^{-1}\bar{\mathbf{y}}^{t}_{o}$
(62)

where first equality follows upon using the definition of
${\mathbf{u}}_{\bar{o}}$ and using again
$\bar{\mathbf{y}}_{\bar{o}}^{t}{\mathbf{u}}_{\bar{o}}=-i({\mathbf{B}}_{c}-{\mathbf{b}})$,
and collecting common terms. Second equality results from simple matrix
algebra. Finally, last line in (59) can be expressed as

$\displaystyle-$ $\displaystyle
i\bar{\mathbf{y}}_{o}(B{\mathbf{1}}_{Q}-{\mathbf{b}})^{-1}\bar{\mathbf{y}}_{\bar{o}}^{t}({\mathbf{1}}_{\bar{o}\bar{o}}-i{\mathbf{A}}_{\bar{o}\bar{o}})^{-1}{\mathbf{A}}_{\bar{o}o}+{\mathbf{A}}_{o\bar{o}}{\mathbf{u}}_{\bar{o}}(B{\mathbf{1}}_{Q}-{\mathbf{b}})^{-1}{\mathbf{u}}_{\bar{o}}^{t}\bar{\mathbf{y}}_{\bar{o}}(B{\mathbf{1}}_{Q}-{\mathbf{B}}_{c})^{-1}\bar{\mathbf{y}}_{o}^{t}$
$\displaystyle=$ $\displaystyle-$ $\displaystyle
i\bar{\mathbf{y}}_{o}(B{\mathbf{1}}_{Q}-{\mathbf{B}}_{c})^{-1}{\mathbf{u}}_{\bar{o}}^{t}{\mathbf{A}}_{{\bar{o}}o}$
(63)

Substituting equations (A), (A) and (A) into (59), collecting terms with
common factor $(B{\mathbf{1}}_{Q}-{\mathbf{B}}_{c})^{-1}$ and defining,
according to the main text,
$\bar{\mathbf{u}}_{o}=i{\mathbf{A}}_{o{\bar{o}}}{\mathbf{u}}_{\bar{o}}+\bar{\mathbf{y}}_{o}$
, (39) promptly follows.

## References

## References

  * Tiesinga  _et al._ (1993) E. Tiesinga, B. J. Verhaar, and H. T. C. Stoof, Phys. Rev. A 47, 4114 (1993).
  * Bloch  _et al._ (2008) I. Bloch, J. Dalibard, and W. Zwerger, Rev. Mod. Phys. 80, 885 (2008).
  * Chin  _et al._ (2010) C. Chin, R. Grimm, P. Julienne, and E. Tiesinga, Rev. Mod. Phys. 82, 1225 (2010).
  * Roy  _et al._ (2013) S. Roy, M. Landini, A. Trenkwalder, G. Semeghini, G. Spagnolli, A. Simoni, M. Fattori, M. Inguscio, and G. Modugno, Phys. Rev. Lett. 111, 053202 (2013).
  * Jachymski and Juilenne (2013) K. Jachymski and P. S. Juilenne, Phys. Rev. A 88, 052701 (2013).
  * van Kempen  _et al._ (2002) E. G. M. van Kempen, S. J. J. M. F. Kokkelmans, D. J. Heinzen, and B. J. Verhaar, Phys. Rev. Lett. 88, 093201 (2002).
  * Hutson and Soldán (2007) J. M. Hutson and P. Soldán, Int. Rev. Phys. Chem. 26, 1 (2007).
  * Naik  _et al._ (2011) D. Naik, A. Trenkwalder, C. Kohstall, F. M. Spiegelhalder, M. Zaccanti, G. Hendl, F. Schreck, R. Grimm, T. M. Hanna, and P. S. Julienne, Eur. Phys. J. D 65, 55 (2011).
  * Ciurylo  _et al._ (2005) R. Ciurylo, E. Tiesinga, and P. S. Julienne, Phys. Rev. A 71, 030701(R) (2005).
  * Frye and Hutson (2017) M. D. Frye and J. M. Hutson, Phys. Rev. A 96, 042705 (2017).
  * Frye  _et al._ (2019) M. D. Frye, B. C. Yang, and J. M. Hutson, Phys. Rev. A 100, 022702 (2019).
  * Feshbach (1967) H. Feshbach, Ann. Phys. 43, 410 (1967).
  * Taylor (2006) J. R. Taylor,  _Scattering Theory: The Quantum Theory of Nonrelativistic Collisions_ (Dover, 2006).
  * Feshbach (1958) H. Feshbach, Ann. Phys. 5, 357 (1958).
  * Domcke (1983) W. Domcke, Phys. Rev. A 28, 2777 (1983).
  * Cohen-Tannoudji  _et al._ (1992) C. Cohen-Tannoudji, J. Dupont-Roc, and G. Grynberg,  _Atom-Photon Interactions: Basic Processes and Applications_ (Wiley, 1992).
  * Simoni  _et al._ (2002) A. Simoni, P. S. Julienne, E. Tiesinga, and C. J. Williams, Phys. Rev. A 66, 063406 (2002).
  * Wilkinson (1965) J. H. Wilkinson,  _The Algebraic Eigenvalue Problem_ (Oxford University Press, New York, 1965).
  * Dalitz (1961) R. H. Dalitz, Rev. Mod. Phys. 33, 471 (1961).
  * Press  _et al._ (1996) W. H. Press, S. A. Teukolsky, W. T. Vettering, and B. P. Flannery,  _Numerical recipes (2nd edition)_ (Cambride University Press, 1996).
  * Simoni  _et al._ (2017) A. Simoni, A. Viel, and J.-M. Launay, J. Chem. Phys. 146, 244106 (2017).
  * Aquilanti  _et al._ (2005) V. Aquilanti, S. Cavalli, D. De Fazio, A. Simoni, and T. V. Tscherbul, J. Chem. Phys. 123, 054314 (2005).
  * Bolda  _et al._ (2003) E. L. Bolda, E. Tiesinga, and P. S. Julienne, Phys. Rev. A 68, 032702 (2003).
  * Derevianko (2003) A. Derevianko, Phys. Rev. A 67, 033607 (2003).
  * Kawaguchi and Ueda (2012) Y. Kawaguchi and M. Ueda, Phys. Rep. 520, 253 (2012).
  * (26) “Woodbury matrix identity — Wikipedia, the free encyclopedia,” https://en.wikipedia.org/wiki/Woodbury_matrix_identity, [Online; accessed 27-July-2021].
  * Bohn and Julienne (1999) J. L. Bohn and P. S. Julienne, Phys. Rev. A 60, 414 (1999).
  * D’Errico  _et al._ (2007) C. D’Errico, M. Zaccanti, M. Fattori, G. Roati, M. Inguscio, G. Modugno, and A. Simoni, New Jour. Phys. 9, 223 (2007).
  * Lysebo and Veseth (2008) M. Lysebo and L. Veseth, Phys. Rev. A 77, 032721 (2008).