# Saturating interaction in coherently coupled two-component Bose-Einstein condensates

R Eid, S Tiengo, M Lévy, T Bourdel

<figure><img src="image_1.png" /></figure>

## To cite this version:

R Eid, S Tiengo, M Lévy, T Bourdel. Saturating interaction in coherently coupled two-component Bose-Einstein condensates. 2025. hal-05190202v1

HAL Id: hal-05190202

https://hal.science/hal-05190202v1

Preprint submitted on 28 Jul 2025 (v1), last revised 28 Nov 2025 (v2)

HAL is a multi-disciplinary open access archive for the deposit and dissemination of scientific research documents, whether they are published or not. The documents may come from teaching and research institutions in France or abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est destinée au dépôt et à la diffusion de documents scientifiques de niveau recherche, publiés ou non, émanant des établissements d'enseignement et de recherche français ou étrangers, des laboratoires publics ou privés.

<figure><img src="image_2.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
Saturating interaction in coherently coupled two-component Bose-Einstein condensates

R. Eid, S. Tiengo, M. Lévy, and T. Bourdel*

Université Paris-Saclay, Institut d'Optique Graduate School,
CNRS, Laboratoire Charles Fabry, 91127 Palaiseau, France

(Dated: July 28, 2025)

Rabi-coupled spinor Bose-Einstein condensates, with competing intra- and interspecies interac-
tions, enable independent control of two- and three-body interactions. We show that coupling can
also drive the system into a strongly nonlinear regime of saturating interaction. More precisely,
the equation of state interpolates between low- and high-density regimes described by two different
two-body scattering lengths. Interestingly, the transition can be determined by the strength of the
coupling. We experimentally demonstrate this saturation phenomenon by measurements of the in-
teraction energy of a Bose-Einstein condensate as a function of the detuning and of the strength of
the Rabi coupling in spin mixtures of potassium 39.

# I. INTRODUCTION

For small perturbations, a physical system usually re-
sponds according to linear response theory. As the per-
turbation strength increases, the response may become
non-linear giving rise to a wide variety of physical phe-
nomena (solitons, bistability, chaotic behavior, etc...).
[1]. Nonlinear response is often modeled by the addi-
tion of other response terms in powers of the pertur-
bation. Such a Taylor expansion approach explains the
higher harmonic generation in a system driven at a single
frequency, as is commonly observed in nonlinear optics
[2]. The Taylor expansion of the nonlinear response is
not always valid for large perturbations, as many physi-
cal systems exhibit saturation of their response at large
driving amplitudes. Typical examples are the magnetic
or electric responses of materials that saturate when the
maximum polarization is reached. In optics, saturation
of the non-linear index of refraction, of absorption, or of
the gain in amplifiers is also commonly observed [3].

In the physics of dilute atomic Bose-Einstein condensates, the non-linearity primary comes from the interaction between atoms [4] which are characterized by the scattering length $a$. In the case of a single-component Bose-Einstein condensate in the weak interaction regime $na^3 \ll 1$, where $n$ is the atomic density, the energy density is dominated by the mean-field interaction scaling as $n^2$. Beyond this paradigm, quantum fluctuations are responsible for the Lee-Huang Yang term scaling non-analytically as $n^{5/2}$ in the energy density [5]. In bosonic mixtures, this term is responsible for the stabilization of quantum droplets [6, 7]. Higher-order interaction terms, for example, corresponding to a three-body interaction, are usually negligible in the dilute regime $na^3 \ll 1$ [8–10].

However, this paradigm is no longer valid in the pres-
ence of coupling between two spin states, i.e. when the
condensate is in a dressed state that is a coherent super-
position of two spin states. A condensate in a dressed

state may experience both reduced usual two-body in-
teractions [11] and large three-body interactions that are
due to possible virtual excitations during collisions [12].
Such three-body processes are increasingly relevant as the
ratio $\gamma$ of the two-body interaction strength to the Rabi
coupling strength increases. Their effect on the equation
of state, that is, the addition of a $n^3$ term, was experi-
mentally measured through modifications in the conden-
sate dynamics [13].

In this paper, we consider the situation where $\gamma \gtrsim 1$ and demonstrate a saturation of the interaction energy as a function of $\gamma$, an effect that can be captured in a mean field framework. In this regime of strong nonlinear behavior, the equation of state can no longer be described as a Taylor series expansion as a function of the density, with terms scaling as $n^p$ describing p-body interactions. In contrast, one needs to compute the exact solution.

The paper is organized as follows. We first describe the
coupled system Hamiltonian and provide a theoretical ex-
planation for the saturation of the interaction energy in a
homogeneous system (with a fixed density $n$). Second, we
present our experiment using $^{39}$K atoms. We explain the
preparation of the condensates in the dressed state and
the measurement of the interaction energy through one-
dimensional time-of-flight expansions. As the conden-
sate lies in the one-dimensional (1D) to three-dimensional
(3D) crossover regime, we finally compare our results
with Gross-Pitaevskii simulations that include the two
spin states, the Rabi coupling, as well as the radial trap.
A good match, with no fit parameter, confirms our inter-
pretation of the data.

II. THEORETICAL DESCRIPTION

We consider a Bose gas in a volume $V$ consisting of
$N$ atoms of mass $m$ with two coupled internal states,
$\sigma = \uparrow, \downarrow$. The three relevant interaction constants $g_{\sigma\sigma'} =$
$4\pi\hbar^2 a_{\sigma\sigma'}/m$ are related to the scattering lengths $a_{\sigma\sigma'}$
and $\hbar$ the reduced Planck constant. In an homogeneous
system with density $n = N/V$, the mean field energy for

<footnote>* thomas.bourdel@institutoptique.fr</footnote>
a condensate in the spinor state ($\phi_↑, \phi_↓$) reads

$$ \frac{E_{MF}}{V} = -\frac{\hbar\Omega}{2}(\phi_{\uparrow}^{*}\phi_{\downarrow} + \phi_{\downarrow}^{*}\phi_{\uparrow}) + \frac{\hbar\delta}{2}(|\phi_{\uparrow}|^2 - |\phi_{\downarrow}|^2) + \sum_{\sigma\sigma'} \frac{g_{\sigma\sigma'}}{2} |\phi_{\sigma}|^2 |\phi_{\sigma'}|^2, \quad (1) $$

where $\delta$ is the coupling detuning and $\Omega$ the Rabi coupling strength. The ground state is found upon minimization of the energy with respect to the internal state. The first term fixes the relative phase of the spinor that we can thus write ($\phi_↑, \phi_↓$) = $\sqrt{n}(\sin(\theta/2), \cos(\theta/2))$ with $\theta \in [0, \pi]$. The energy is then

$$ \frac{E_{MF}}{N} = -\frac{\hbar\Omega}{2}\sin\theta - \frac{\hbar\delta}{2}\cos\theta + \frac{gn}{2} - \frac{\Delta gn}{4}\cos\theta - \frac{\bar{g}n}{2}\sin^2\theta, \quad (2) $$

where $g = (g_{↑↑} + g_{↓↓})/2$, $\Delta g = g_{↑↑} - g_{↓↓}$, and $\bar{g} = (g_{↑↑} + g_{↓↓} - 2g_{↑↓})/4$.

The minimization of Eq. 2 as a function of $\theta$ leads to the following implicit equation for the spin composition

$$ \cot\theta = \frac{\delta}{\Omega} + \alpha - 2\gamma \cos\theta, \quad (3) $$

where $\alpha = \frac{\Delta gn}{2\hbar\Omega}$ and $\gamma = \frac{\bar{g}n}{\hbar\Omega}$. The spin magnetization and, equivalently, the mixing angle $\theta$ thus depend not only on $\delta/\Omega$ but also on the density $n$. In the low-density regime, i.e. for $\alpha, \gamma \ll 1$, the density dependence of the spin composition can be calculated in perturbation. This effect was experimentally observed [14]. At the same order in perturbation, we find two- and three-body corrections on the interaction energy, generalizing the results from [13] to the asymmetric case

$$ \frac{E_{MF}}{N} \approx \epsilon_{-} + g_{2}\frac{n}{2} + g_{3}\frac{n^{2}}{3} \quad (4) $$
with $\epsilon_{-} = -\hbar\sqrt{\Omega^{2} + \delta^{2}/2}$,
$$ g_{2} = g - \frac{\bar{g}}{1 + \delta^{2}/\Omega^{2}} - \frac{\Delta g}{2} \frac{\delta/\Omega}{\sqrt{1 + \delta^{2}/\Omega^{2}}}, $$
$$ g_{3} = -\frac{3\bar{g}^{2}}{\hbar\Omega} \frac{\delta^{2}/\Omega^{2}}{(1 + \delta^{2}/\Omega^{2})^{5/2}} + \frac{3\bar{g}\Delta g\delta/\Omega}{2\hbar\Omega(1 + \delta^{2}/\Omega^{2})^{2}} - \frac{3\Delta g^{2}}{16\hbar\Omega(1 + \delta^{2}/\Omega^{2})^{3/2}}. $$

$\epsilon_{-}$ is the one-body Rabi energy of the lower-dressed state $|-\rangle$. The two-body coupling constant $g_2(\delta/\Omega)$ corresponds to the interaction between two atoms in $|-\rangle$. It interpolates between $g_{↑↑}$ and $g_{↓↓}$, reaching intermediately a minimum $g_{min} = (g_{↑↑}g_{↓↓} - g_{↑↓}^2)/4\bar{g}$ in the experimentally relevant case $g_{↑↓} < g_{↑↑}, g_{↓↓}$. The variation of this two-body term was experimentally measured in [11]. The three-body coupling constant $g_3(\delta/\Omega, \Omega)$ scales as $1/\Omega$ and was also experimentally demonstrated [13].

We now consider the non-pertubative case where the interactions are not small as compared to the Rabi coupling strength $\Omega$ corresponding to $\alpha, \gamma \gtrsim 1$. In this

<figure><img src="image_3.png" /><figcaption>FIG. 1. Interaction energy per particle $E_{int}/N = (E_{MF}/N - \epsilon_{-})$ as a function of $\gamma$. The specific parameters are $\delta/\Omega = 1$, $a_{↑↑} = 86.8 a_0$, $a_{↓↓} = 33.3 a_0$, and $a_{↑↓} = -53.2 a_0$, where $a_0$ is the atomic Bohr radius. The solid blue curve is the exact mean-field calculation, the dash-dot red line is the linear two-body results, and the dash green curve is the quadratic two-plus three-body approximation (Eq. 4).</figcaption></figure>

regime, the energy is evaluated by numerically solving Eq. 3, and inserting the result in Eq. 2. In Fig. 1, we plot the mean field interaction energy per particle $E_{int}/N = (E_{MF}(n)/N - \epsilon_{-})$ as a function of the parameter $\gamma$ for a detuning $\delta/\Omega = 1$ and for scattering lengths that are relevant to our experiment with $^{39}$K. We observe that the interaction energy first grows linearly with the density because of a positive coupling constant in the dressed state. Then it bends down according to the attractive three-body interaction term $g_3$. When $\gamma$ increases above 1, the energy per particle tends to saturate. A Taylor expansion as a function of density is no longer valid for large $\gamma$. In fact, in the limit $\gamma \gg 1$, the system minimizes its energy by choosing the condensate spin state corresponding to a minimal interaction $\cos\theta = \alpha/2\gamma = \Delta g/4\bar{g}$ and a pure two-body interaction is recovered with the coupling constant $g_{min}$. In the case $g_{↑↓} = -\sqrt{g_{↑↑}g_{↓↓}}$, $g_{min} = 0$ and the condensate in the dressed state has an interaction energy that saturates to a constant value at large $\gamma$.

The previous results can be experimentally investigated by modifying $\Omega$ in order to change $\gamma$ over a wide range. In Fig. 2, we thus plot the interaction energy as a function of $\delta/\Omega$ for different values of $\gamma$ assuming a homogeneous density $n$. For small values of $\gamma$, the energy is only weakly modified by nonlinear effects, whereas for $\gamma \gtrsim 1$, the effect is large. The saturation of the interaction energy at large $\gamma$ manifests itself as a broadening of the drop in interaction energy as a function of $\delta/\Omega$. In particular, for $\gamma \gg 1$, the width of the dip increases linearly with $\gamma$. Note that in this regime, the Taylor expansion up to three-body interaction is very wrong. Already for $\gamma = 1$, this approximation (blue dashed curve in Fig. 2) leads to negative mean-field energies, although
<figure><img src="image_4.png" /><figcaption>FIG. 2. Interaction energy per particle as a function of the detuning $\delta/\Omega$ for different values of $\gamma$. The scattering lengths are the same as in Fig. 1. The interaction energy is normalized such that it gives the scattering length in unit of $a_0$ in the regime where the interactions are two-body ($g_0 = 4\pi\hbar^2 a_0/m$). The bare scattering lengths $a_{\uparrow\uparrow}$ and $a_{\downarrow\downarrow}$ are thus found for large values of $|\delta/\Omega|$. The solid lines correspond to the exact mean-field calculation for different values of $\gamma$. The blue dashed-line corresponds to the two- plus three-body approximation given by Eq. 4 for $\gamma = 1$.</figcaption></figure>

the interaction energy is always positive for the chosen scattering lengths.

### III. EXPERIMENTAL RESULTS AND ADJUSTMENT

The experiment starts with a $^{39}$K Bose-Einstein condensate with $N = 1.2 \times 10^5$ atoms in the $|F = 1, m_F = -1\rangle = |\downarrow\rangle$ spin state in a crossed optical dipole trap with oscillation frequencies $169 \times 169 \times 26$ Hz. The coupling to the state $|F = 1, m_F = 0\rangle = |\uparrow\rangle$ is performed through a radio frequency (RF) field. The magnetic field $B = 56.85$ G is chosen so that the scattering lengths are $a_{\uparrow\uparrow} = 84.3$, $a_{\downarrow\downarrow} = 33.3$, $a_{\uparrow\downarrow} = -53.2$, where $a_0$ is the atomic Bohr radius [15]. In short, the measurement idea consists in preparing the dressed state $|\text{--}\rangle$ at a specific $\delta/\Omega$ and then monitoring the longitudinal expansion of the cloud once the weak longitudinal confinement is removed. At long times of flight $t$, the expansion is ballistic and we can simply deduce the expansion energy from the relation $E_{\text{exp}}/N = m\sigma^2/2t^2$, where $\sigma$ is the measured longitudinal rms size of the cloud. For our condition, this expansion energy predominantly reflects the initial interaction energy as the initial kinetic energy is negligible (see below).

The experiment is repeated with different radio frequency powers or equivalently different Rabi coupling frequencies $\Omega/2\pi$ between 0.95 and 30 kHz. These values of $\Omega$ are chosen in order to explore both the regimes of weak ($\gamma \ll 1$) and strong ($\gamma \gtrsim 1$) nonlinear inter-

action. Measurements at 0.95 kHz Rabi frequency are made possible by precise control of detuning $\delta$ down to $\sim 100$ Hz. This level of control requires a ppm magnetic field stability and is obtained thanks to a feedforward compensation of both the current in the bias coils and the external magnetic field [16].

The dressed state $|\text{--}\rangle$ is prepared with an adiabatic passage from the state $|\downarrow\rangle$. The RF detuning is linearly swept from $9\Omega$ to the desired final value $\delta$. The total sweep time $\tau=9$ ms is chosen in order to be adiabatic with respect to the spin degree of freedom even at the lowest $\Omega$ values. In addition, this time is sufficiently large as compared to $1/\omega_{\perp}$ such that the radial wave-function can adiabatically adapt itself to the interaction change during the sweep [17].

In the longitudinal direction, the cloud dynamics is much slower, and the interaction change induced by the RF sweep can almost be considered as a sudden quench. However, this approximation is not strictly true. If the RF sweep is done with (resp. without) the longitudinal trap on, the atomic cloud starts to slightly shrink (resp. expands). We reduce this residual longitudinal dynamics by adjusting the time at which the longitudinal trap is turned off during the RF sweep. In principle, the most appropriate time depends on the final value of $\delta/\Omega$ and also on $\Omega$ (through the nonlinear effects). However, we can consider that the additional kinetic energy becomes relevant only when the interaction energy is small, that is, around $\delta/\Omega \approx 0.24$, where $a = a_{\text{min}} \sim 0$. To minimize the longitudinal dynamics during the RF sweep, we thus adjust the longitudinal trap extinction time for each value of $\delta$, ensuring that the induced kinetic energy is minimized when $\delta/\Omega \approx 0.24$. In practice, this is done by searching for the minimal longitudinal expansion at a long time of flight of 300 ms. We find that the proper trap extinction time varies from 0.5 ms before the end of the sweep for large $\Omega$ to 5 ms before for $\Omega = 0.95$ kHz. The precision of this adjustment allows us to reduce the initial kinetic energy per particle to below 6 Hz, a negligible quantity compared to the measured interaction energies.

In Fig. 3, we plot the longitudinal rms size of the condensate (obtained from a Gaussian fit) after a time of flight $t = 62.6$ ms (during which the radial confinement is kept on) as a function of $\delta/\Omega$ and for different values of $\Omega$. Each data set shows a minimum around the value $\delta/\Omega \approx 0.24$ where the two-body coupling constant $g_2$ in the $|\text{--}\rangle$ dressed state is minimum. For large and positive $\delta/\Omega$, the expansion corresponds to the one of the condensate in $|\downarrow\rangle$ with the scattering length $a_{\downarrow\downarrow}$. For large and negative $\delta/\Omega$, the expansion is even faster, as the scattering length $a_{\uparrow\uparrow}$ is larger. A striking feature is that the width of the expansion drop is significantly larger as $\Omega$ decreases. This behavior qualitatively follows the theoretical expectation from Fig. 2 although the initial density is not homogeneous. Even without quantitative analysis, we thus see that in our experiment the strongly nonlinear regime corresponding to $\gamma \gtrsim 1$ is reached. This experimental finding is in line with the estimation $\gamma = 1.6$.
<figure><img src="image_5.png" /><figcaption>FIG. 3. Size after 1D expansion of condensates in the $|−⟩$ dressed state as a function of the detuning $\delta/\Omega$ and for different values of $\Omega$ (see legend). The points correspond to the experimental measurements and some typical error bars are shown. In the main graph, the solid lines correspond to the mean-field numerical simulations (see text). At $\delta/\Omega \approx 0.24$, the two-body scattering length is already minimal and all theoretical mean-field curves are superimposed. The inset shows an enlargement of the experimental data at low $\delta//\Omega$. The points are linked by dashed lines for clarity. The theoretical lines are omitted in the inset because the theoretical model does not capture the vertical offsets between the curves (see text).</figcaption></figure>

for the initial condensate peak density and $\Omega = 0.95$ kHz.

Let us turn to a more quantitative analysis with a numerical simulation of the experimental situation. We have to take into account that the nature of the confinement changes depending on the parameters. Around the expansion minimum, the interaction energy is small and the condensate is in a quasi-1D regime (the radial wavefunction is close to the trap radial ground state). In contrast, for large $|\delta/\Omega|$ and for the initial condensate, the interactions are stronger and the condensate is effectively 3D. In order to account for these two behaviors as well as the dimensional crossover regime, we numerically calculate the radial profile of the spinor condensate wavefunction when the interactions are modified during the sweep. This is done in spinor Gross-Pitaevskii simulations. In this framework, the previously discussed mean-field effects giving rise to the saturating interactions in the equation of states are included.

In more detail, we assume an adiabatic following of the radial wavefunction, and thus we search for the radial ground state as a function of the local longitudinal 1D densities. This is done through 2D imaginary-time Gross-Pitaevskii simulations of the spinor wavefunction using the GPELab toolbox [18]. We then calculate the total energy, which includes potential, kinetic, Rabi coupling, and interaction contributions. Finally, we integrate this energy along the longitudinal trap direction, assuming that the density profile is the one of the initial condensate, i.e. a 3D Thomas-Fermi profile.

The total energy $E_{\text{tot}}$ transforms into kinetic energy during expansion. After a significant expansion, the expected size is thus found from the following relation $m\sigma^2/2t^2 = E_{\text{tot}}/N - \hbar\omega_{\perp}/2 - \epsilon_0$, where the last two terms come from the fact that the condensate adiabatically evolves to the noninteracting spinor ground state (with energy $\epsilon_0$) and to the radial ground state during expansion (with energy $\hbar\omega_{\perp}/2$).

The results of the above simulations match quite well the experimental findings with no fit parameter (Fig. 3). The width of the minimum expansion drop as a function of $\Omega$ is nicely reproduced, showing that the saturation of the mean-field interaction energy is captured in our mean-field framework. The slightly smaller size at large and negative $\delta/\Omega$ can be attributed to the ~20% losses that are experimentally observed and are due to three-body recombination processes. A smaller atom number reduces the interaction energy as compared to the simulation where no losses are included.

Looking carefully at the experimental results, one can notice that the observed minimum sizes (at $\delta/\Omega \approx 0.24$) are higher when $\Omega$ increases (see inset in Fig. 3). This behavior is not reproduced in the above mean-field approach, where, at this point, the interaction corresponds to a pure two-body interaction with a coupling constant $g_{\text{min}}$, which is independent of $\Omega$. We attribute this behavior to beyond mean field effects due to virtual excitations to the $|+⟩$ states, leading to a slight increase of the minimal interaction energy when $\Omega$ increases [12, 19]. In the context of this paper, we do not focus on this residual energy shift, as its magnitude is much lower than the dominant mean-field effects studied above. Its quantitative measurement is demanding and would require an improved experimental control on the initial kinetic energy of the cloud.

## IV. CONCLUSION

In this paper, we have shown that the density dependence of the equation of state of a coherently driven two-component condensate can be tuned. For some parameters, the interaction term evolves from a two-body repulsive interaction at low densities, through a regime that can be described with additional three-body interactions at intermediate densities, to yet another regime of two-body interaction with a reduced two-body constant at large densities. This second two-body interaction can be positive or negative, depending on the initial spin-dependent scattering lengths. If it is close to zero (as in our work), then it leads to an unusual equation of state with an energy saturation as a function of the density.

Our findings open the way to the study of novel nonlinear effects in the physics of condensates. Already, at equilibrium, we anticipate density profiles that differ from the usual Thomas-Fermi profiles. The condensate excitation dynamics (vortices, solitons, low energy excitation modes) will also be modified. For example, in our con-
figuration, the speed of sound is expected to decrease at high densities, in contrast to the usual behavior for a Bose-Einstein condensate with repulsive interactions.

ACKNOWLEDGMENTS

This research was supported by CNRS, Ministère de l'Enseignement Supérieur et de la Recherche, Quantum

Paris-Saclay, Région Ile-de-France in the framework of Domaine d'Intérêt Majeur Quantip (H3B), France 2030 in the framework of PEPR Quantique (ANR-23-PETQ-0001), the Simons Foundation (Award No. 563916, localization of waves), France and Chicago collaborating in the sciences.

[1] S. H. Strogatz, *Nonlinear Dynamics and Chaos: With Applications to Physics, Biology, Chemistry, and Engineering*, 3rd ed. (CRC Press, Boca Raton, FL, 2024) p. 616.

[2] Y.-R. Shen, *The Principles of Nonlinear Optics* (Wiley-Interscience, New York, 1984).

[3] R. W. Boyd, *Nonlinear Optics*, 4th ed. (Academic Press, 2020).

[4] P. Pitaevskii, Lev and S. Stringari, *Bose-Einstein Condensation* (Clarendon Press, Oxford Univ. Press, Oxford, UK, 2003) p. 382.

[5] T. D. Lee, K. Huang, and C. N. Yang, Eigenvalues and eigenfunctions of a bose system of hard spheres and its low-temperature properties, Physical Review **106**, 1135 (1957).

[6] D. S. Petrov, Quantum mechanical stabilization of a collapsing Bose-Bose mixture, Physical Review Letters **115**, 155302 (2015).

[7] C. R. Cabrera, L. Tanzi, J. Sanz, B. Naylor, P. Thomas, P. Cheiney, and L. Tarruell, Quantum liquid droplets in a mixture of Bose-Einstein condensates, Science **359**, 301 (2018).

[8] E. Braaten and A. Nieto, Renormalization effects in a dilute Bose gas, Phys. Rev. B **55**, 8090 (1997).

[9] T. Köhler, Three-body problem in a dilute Bose-Einstein condensate, Phys. Rev. Lett. **89**, 210404 (2002).

[10] P. M. A. Mestrom, V. E. Colussi, T. Secker, and S. J. J. M. F. Kokkelmans, Scattering hypervolume for ultracold bosons from weak to strong interactions, Phys. Rev. A **100**, 050702 (2019).

[11] J. Sanz, A. Frölian, C. S. Chisholm, C. R. Cabrera, and L. Tarruell, Interaction control and bright solitons in

coherently-coupled Bose-Einstein condensates, Physical Review Letters **128**, 013201 (2022).

[12] S. Tiengo, R. Eid, and T. Bourdel, Three-body interactions in Rabi-coupled bose gases: A perturbative approach, Phys. Rev. A **111**, 053319 (2025).

[13] A. Hammond, L. Lavoine, and T. Bourdel, Tunable three-body interactions in driven two-component Bose-Einstein condensates, Phys. Rev. Lett. **128**, 083401 (2022).

[14] A. Farolfi, A. Zenesini, R. Cominotti, D. Trypogeorgos, A. Recati, G. Lamporesi, and G. Ferrari, Manipulation of an elongated internal Josephson junction of bosonic atoms, Phys. Rev. A **104**, 023326 (2021).

[15] E. Tiemann, P. Gersema, K. K. Voges, T. Hartmann, A. Zenesini, and S. Ospelkaus, Beyond Born-Oppenheimer approximation in ultracold atomic collisions, Phys. Rev. Res. **2**, 013366 (2020).

[16] S. Tiengo, R. Eid, M. Apfel, G. Brulin, and T. Bourdel, A simple magnetic field stabilization technique for atomic Bose-Einstein condensate experiments, Review of Scientific Instruments **96**, 063201 (2025).

[17] This is not strictly the case in particular for negative $\delta$ where residual radial oscillations have been observed after the sweep. However, we estimate that these radial oscillations do not significantly affect our longitudinal expansion results.

[18] X. Antoine and R. Duboscq, Gpelab, a Matlab toolbox to solve Gross-Pitaevskii equations I: computation of stationary solutions, Computer Physics Communications **185**, 2969 (2014).

[19] L. Lavoine, A. Hammond, A. Recati, D. S. Petrov, and T. Bourdel, Beyond-mean-field effects in Rabi-coupled two-component Bose-Einstein condensate, Phys. Rev. Lett. **127**, 203402 (2021).