Squeezing oscillations in a multimode bosonic Josephson junction

,

Understanding the role of quantum fluctuations and entanglement in interacting many-body systems is of critical importance for the development of quantum technologies, such as quantum metrology [1,2] and quantum simulation [3].Especially the ability to prepare entanglement [4] in quantum many-body systems is pivotal.In this context, ultracold atoms have proven to be a versatile platform [1].
Internal spin degrees of freedom (DoFs) offer a high level of control and using a single spatial mode spin squeezed states could be generated experimentally [5][6][7][8].Creating entangled quantum states in motional DoFs, for example in tunnel-coupled BECs in double wells (DWs) [9,10] is, on the contrary, less explored; the impact of spin squeezing and entanglement has only been studied with respect to the global observables.
For spatially extended systems, where more than one longitudinal mode is occupied, the interplay between transverse (spin) and longitudinal DoFs leads to new physical phenomena.For coupling internal states interesting dynamical phenomena as well as a high degree of control have been shown [11,12]; for tunnel-coupled systems, the influence of the multimode situation on Josephson oscillations has been studied [13].However, so far the quantum regime has not been accessible.
Interestingly, in the one-dimensional (1D) regime these systems are excellent quantum simulators for the sine-Gordon field theory [14,15], which was mainly explored in tunnel-coupled systems [16,17].It is worth pointing out that BECs in DWs at equilibrium can be experimentally achieved by direct cooling in DWs [16].But fluctuations introduced by the finite temperature make it challenging to prepare the system in the quantum correlated regime [18].However, it has been shown that using a splitting routine from a single condensate, quantum correlated states can be prepared in a multimode scenario [10].In this work, we combine the ability to prepare quantum correlated states with spatially resolved measurements of the sine-Gordon fields entering a new regime.
We realise a multimode bosonic Josephson junction (BJJ) [19,20] by splitting 1D quasi-BECs into a DW arXiv:2304.02790v2[cond-mat.quant-gas]8 Feb 2024 (see Fig. 1 upper panel).Instead of studying Josephson oscillations [13,[21][22][23][24][25], we prepare the system at the stable fixed point of the classical phase space and access its quantum properties.Despite the stationary expectation values of both relevant observables, relative number, and relative phase, we observe oscillatory dynamics of the quantum fluctuations which result from the evolution of an initial quantum state different from the ground state at zero temperature.We employ these non-equilibrium dynamics of the fluctuations in the BJJ to foster strong spin squeezing in decoupled condensates.Furthermore, we study its impact on spatial phase correlations in the multimode system.
To probe the multimode properties of the prepared states we use a spatially resolved detection of the relative phase (see Fig. 1 lower panel).The observations in this observable, enable us to demonstrate how the improved number squeezing enhances spatial correlations between two decoupled 1D condensates.We are able to connect quantum properties resulting from the preparation process to enhanced phase correlations in a regime where a sine-Gordon model is realized.The enhanced squeezing can in principle be directly related to lower effective temperatures in the prethermalised regime [26,27] and our work thus shows a pathway to prepare sine-Gordon field simulators in a regime dominated by quantum fluctuations.

REALISATION AND READOUT OF A MULTIMODE BJJ
We realise a multimode BJJ with 1D tunnel-coupled quasi-BECs of 87 Rb atoms trapped magnetically under an Atom Chip [28].The DW potential is generated with the Radio-Frequency (RF) dressing [29] technique.The dressing is calibrated such that the precise transformation from a single well to a DW is achieved by ramping up the amplitude A of the RF field whose value is normalised to the maximal current.We display in Fig. 1 the transverse trap configuration and typical experimental readout with our single-atom-sensitive fluorescence imaging system [30] after Time-of-Flight of t F = 43.4ms.The aspect ratio between the radial and axial trap frequency is ∼ 100.
Each of the quasi-condensate in a DW can be described in the density-phase representation by Ψ j = ρ j (z) × exp(−iϕ j (z)), where j ∈ [L, R] labels the left (L) and the right (R) condensate.The local density ρ j and local phase ϕ j are spatially varying and fluctuating fields.This is illustrated by the wiggly brown lines in Fig. 1.We are interested in the local fields in the relative DoF, namely the local relative density ρ − (z) = ρ R (z) − ρ L (z) and local relative phase ϕ(z) = ϕ R (z) − ϕ L (z).For finitely tunnel-coupled condensates, the spatially dependent relative phase field ϕ(z) can be described by the sine-Gordon field theory [14].
Our imaging allows us to measure the relative phase with spatial resolution (see Fig. 1 lower panel).From this, we can access the two-point phase correlations to study the sine-Gordon physics in the multimode regime, as will be discussed later.Valuable insights into the dynamical quantum properties of BJJs can already be obtained from a single-mode description of the BJJ [31].Therefore, we will first investigate the properties of the two global observables of the 1D BJJ, the relative atom number, N − = z ρ − (z) and the relative phase, Φ = z ϕ(z), which correspond to the spatial zeroth mode.The macroscopic dynamics of the global observables N − and Φ in a BJJ can be described by the two-mode Bose-Hubbard (BH) model where N = N L + N R is the total atom number, n = N − /N is the relative imbalance, U is the interaction strength and J is the single particle tunnel coupling strength.The Josephson regime of the BJJ is indicated by 1/N ≪ U/2J ≪ N .From Eq. ( 1), we can obtain the Josephson oscillation frequency [19,20] of the mean values ⟨Φ⟩ and ⟨N − ⟩ also known as the plasma frequency.Here cos Φ 0 indicates the initial phase coherence factor.From Eq.(2), we see that f p depends explicitly on the single particle tunnel coupling strength J.
For the multimode many-body dynamics dominating in decoupled DWs, the two-mode BH model does definitely not give a full description.However, we will see in the following that it nicely captures the dynamics of the squeezing in the spatial zeroth mode.Exact modelling of the splitting process together with the ensuing dynamics are challenging when taking into account the true multimode quantum dynamics; thus, quantum simulations can provide new insights.

PREPARING QUANTUM-CORRELATED STATES
To prepare strongly correlated BECs, we split a single BEC into two by transforming a single well into a DW (see Fig. 2a).BEC splitting as a pathway to generate quantum correlated states was investigated in many earlier works [9,32,33]; here, we introduce and implement a new scheme which is based on understanding the nonequilibrium dynamics after a non-adiabatic splitting into a tunnel-coupled DW.In the following, we introduce the relevant quantities and summarize our findings which lay the basis for observing the non-equilibrium dynamics of the prepared squeezed state.N and ξ 2 Φ , right after a linear ramp-up at a constant ramp speed to various tunnelcoupled DWs, denoted with the rescaled RF dressing amplitude A. In the upper x-axis, the parameter 1/N ≪ U/2J ≪ N indicates that we are in the Josephson regime of BJJ.The coloured bands mark the expected ground state squeezing factors of BJJ at zero temperature (Eq.( 6)) with atom number N ∈ [2, 5] × 10 3 .b. Schematics of state evolution in the classical phase space.Left column illustrates the expected ground state fluctuation of BJJ according to Eq. ( 7), where orange (green) arrows represent projection noise in the relative number (phase) quadrature.The right column represents the out-of-equilibrium quantum state evolution.The triangle markers indicate single classical realizations orbiting along the equipotential lines at the plasma frequency fp.In the Josephson regime, the distribution of out-of-equilibrium states rotates and deforms which results in squeezing oscillations.A full period of the quantum state evolution around the stable fixed point (n = N−/N = 0 and Φ = 0) conforms to a half period of evolution of single realisations (mean values) in the phase space.
To quantify the fluctuations of the observables, we define the squeezing factors where ∆ 2 N − and ∆ 2 Φ represent the statistical variance, evaluated as in Ext.Data Fig. 1.We use the quadrature projection noise of spin coherent states in the denominators.Hence ξ 2 N = ξ 2 Φ = 1 represents the standard quantum limit.Furthermore, spin-squeezed states [34,35] representing a class of entangled states are characterised by spin squeezing factor ξ 2 s = ξ 2 N /⟨cos Φ⟩ 2 < 1, where ⟨cos Φ⟩ is the phase coherence factor.
The ground state of the many-body Hamiltonian [36] exhibits a growing degree of number squeezing in less coupled DWs owing to repulsive interatomic interactions.In Fig. 2a, we show the experimentally inferred squeezing factors in N − and Φ as a result of BEC splitting with finite duration, denoted with rescaled RF dressing amplitude A. We conduct ∼ 200 repetitions for each measurement to ensure reliable statistics.At the beginning of the Josephson regime, U/2J ∼ 1/N , the measured relative number fluctuations of split BECs approach the expected ground state number squeezing ξ 2 N,0 (see Ap-pendix Eq. ( 7)).
During further ramp-up of the RF amplitude A, the adiabaticity condition can no longer be satisfied.This breakdown is shown explicitly in Fig. 2a, where ξ 2 N become increasingly larger than ξ 2 N,0 in less tunnel-coupled DWs.This initialised out-of-equilibrium state with BEC splitting is expected to evolve dynamically in the Josephson regime (see Fig. 2b).As we will discuss in the next sections, this evolution of the quantum state in phase spaces with number-squeezed ground states is not just a rotation, but also additional shearing.
The fluctuations in the relative phase are always much above the expected ground state phase squeezing factor ξ 2 Φ,0 as shown in Fig. 2b.This is due to the interatomic interaction induced phase diffusion [37] during the ramp and can be prevented by increasing the splitting speed of linear ramps.The trade-off of faster splitting is poorer number squeezing right after the splitting.Thus alternative splitting routines beyond linear single ramps are sought after to enhance spin squeezing in decoupled traps.

SQUEEZING OSCILLATIONS IN CONJUGATE QUADRATURES
We prepare two BECs in a strongly coupled DW (A = 0.5) by linearly ramping up from a single well; this prepares the system at the stable fixed point of the classical phase space, i. e. ⟨N − ⟩ = 0 and ⟨Φ⟩ = 0 (see Ext. Data Fig. 3) with phase space fluctuations different from the ground state.By tracking the evolution of this out-of-equilibrium quantum state in the strongly coupled double well, we observe the dynamics of the quantum fluctuations of the conjugate observables (see Fig. 3).
split hold time FIG. 3. Squeezing oscillations in a coupled double well.a. Dynamics of the fluctuations of relative phase (green) and relative number (orange markers) for variable hold time in the coupled trap after the splitting.We observe oscillatory dynamics with comparable frequency and relative phase shift of π.The oscillatory dynamics is not a simple linear rotation of the state in phase space but results from the interplay of tunnel-coupling and onsite interaction and leads to a rotation plus deformation of the state (see inset).From a sinusoidal fit (solid lines) we extract the frequencies which are roughly twice the plasma frequency as expected from the rotation and shearing of the distribution in the phase space (see inset and The squeezing factors in both quadratures undergo oscillatory dynamics.Strikingly, the number quadrature stays always squeezed while the phase never gets squeezed below the SQL, i.e. the oscillatory dynamics is not a simple linear rotation of the state in phase space.This results from the interplay of tunnel coupling and onsite interaction and leads to a rotation and deformation of the state (see Fig. 2b and Fig. 3 inset).We can understand the frequency of the oscillation with the intuitive picture that a π rotation of a single realisation corresponds to a 2π rotation of quantum state distribution in phase space (see Fig. 2b and Appendix).We thereby deduce that the squeezing oscillations are twice as fast as the Josephson oscillations of the mean (with plasma frequency f p ); this is in accordance with a semiclassical analysis for Raman coupled BECs [38].
We fit the observed squeezing factors in Fig. 3 with a sine function to determine the squeezing oscillation frequency and obtain f ξ = 567 (29) Hz in relative number N − quadrature with total atom number N = 4154 (35) and f ξ = 649 (33) Hz in relative phase Φ quadrature with N = 4302 (45) (see Appendix).The measured squeezing oscillations in both quadratures as expected from the nonlinear Josephson dynamics match with twice the experimentally measured plasma frequency f p and have a π phase shift with respect to each other.Combining the complementary measurements in Fig. 3, we infer the mean quantum state fluctuations in the phase space to be ξ 2 N • ξ 2 Φ ≈ 3.5 with ξ 2 N = 0.44(2) and ξ 2 Φ = 8.0 (7).

TWO-STEP SEQUENCE FOR OPTIMISED SPIN SQUEEZING
So far mostly simple linear ramps have been used to generate squeezing.In single-mode situations theoretical studies showed that non-linear ramps generated by optimal control can lead to better squeezing [39,40]; for our multimode BJJ the development of a many-body simulation, required for performing open loop optimisation on spin squeezing, is not possible, as the splitting process is hard to calculate.Thus we develop an experimentally tractable two-step approach (see Fig. 4) based on the observed squeezing oscillations in the coupled DW; in single-mode spinor BECs a similar approach has been used to achieve spin-squeezed ground states [41,42].
Here, we employ it for optimising the spin squeezing in a decoupled DW (J = 0) in a multimode situation (see Fig. 4).This is of particular interest for their application as a sensitivity-enhanced matter-wave interferometer [1,2,10].
We show in Fig. 4 that the generated squeezing oscillation in the coupled DW is successfully preserved after the ramp to the decoupled DW with two-step splitting; this means we are able to manipulate the quantum properties in the decoupled trap by a holding time in the coupled trap.Furthermore, the number squeezing is even further enhanced during the second ramp.With this approach, we obtain a phase coherence factor of ⟨cos Φ⟩ = 0.86 +0.01 −0.02 .This leads to a spin squeezing factor [5] of ξ 2 s = −9.2+1.9 −3.0 dB (with no detection noise correction, ξ 2 s = −4.0(1.1)dB, see Appendix), which witnesses many-body entanglement [35].
In contrast, a single linear ramp with a similar duration yields lower spin squeezing of ξ 2 s ≈ −1.5 dB, demonstrat-split split hold time FIG. 4. Two-step sequence for optimised spin squeezing For optimising the spin squeezing in the decoupled trap (A = 0.65), we perform a two-step sequence (black line in left schematic) instead of a simple linear ramp (blue line).By varying the hold time in the coupled trap (A = 0.5) and performing the final measurement after adding a second linear ramp to the decoupled trap, we observe an oscillatory behaviour.This oscillation is similar to Fig. 3, but overall better number squeezing.For the optimal spin squeezed point (red marker) we find −9.2 dB spin squeezing instead of −1.5 dB for a single linear ramp with a similar total duration as the two-step sequence.
ing a significant gain with two-step splitting.A simple way to understand this is that two-step splitting enables us to achieve optimal number squeezing at an earlier stage of the splitting procedure and, as a result, suppressing the phase diffusion more efficiently.

TUNING SQUEEZING OSCILLATIONS
To achieve controlled preparation of strongly correlated states with two-step splitting, it is crucial to understand the control parameters and tunability of the squeezing oscillation frequency.We identify two strategies: the first one involves tuning the parameters in a static BJJ, and the second one involves introducing time dependence on the control parameters through transversal motions.
First, we explore experimentally the frequency scaling in elongated BJJ based on Eq.( 2).In Fig. 5a, we show the measured number squeezing oscillations in different DWs with decreasing single particle tunnel coupling J (bottom to top).We observe squeezing oscillations with frequencies spanning more than one order of magnitude.In Fig. 5b, we plot collectively the extracted frequencies and as a comparison the expectation of 2f p (grey band) estimated directly from Eq. (2).Here, J for each A is inferred from simulated DW potential and experimental atom number N ∈ [2, 5] • 10 3 .We find good agreement between the experimentally observed fluctuation dynamics on the global observable N − (zeroth mode) of the multimode BJJ and the simple calculations based on the two-mode BH model.This confirms first of all that the observed squeezing oscillations indeed originate from stationary nonlinear Josephson dynamics.
On top of which, we show in Fig. 5c  of the squeezing oscillations in Fig. 5a match quantitatively with the expected evolution of the quantum state during the linear ramp.The experimental results in Fig. 5a are obtained after single linear ramps with the same splitting speed.We can estimate the additionally gained phase during the ramp compared to the first experimental data point at A = 0.48 (set as t = 0) as the sum of constantly evolving plasma frequencies f p (t i ) (Eq. ( 2)) over a ramp time ∆T , namely ∆T • ∆T t=0 2πf p (t).Despite the quantum state evolution in the phase space (Fig. 2b) not being a simple rotation, the differences between the initial phases in different DWs built up during the ramp can be surprisingly well predicted with this simple estimation.
The other DoF for tuning squeezing oscillation is the total atom number N .In Fig. 5c we show experimentally measured squeezing oscillation frequency f ξ as a function of N (for fixed A = 0.5) and compare them directly with experimentally measured plasma frequency 2f p and a solid line derived from Eq. ( 2).In this strongly coupled DW, we can tune squeezing oscillation frequency from 300 Hz to 800 Hz by only adjusting the atom number.
To expand the tuning capabilities, we incorporate an active modulation on the tunneling mechanism.By performing a splitting quench, we excite out-of-phase transverse sloshing between the two BECs at the transverse trap frequency f x = 1418(10) Hz.We depict in Fig. 6a this induced motion with the inferred distance d.This motional excitation drives the effective tunnel coupling J periodically at the trap frequency f x .We show in Fig. 6b how this periodic drive enforces ξ 2 N to oscillate at frequencies comparable with f x .Furthermore, we can reproduce this oscillation frequency with a two-step splitting quench to the decoupled DW (see Fig. 6c).
By utilizing this method, we boost the preparation of spin-squeezed states with two-step from two perspectives: reduced ramp times and faster squeezing dynamics.Our observations potentially facilitate investigations into captivating phenomena like parametric resonance [43] and Floquet engineering [44].

SQUEEZING-PROTECTED SPATIAL CORRELATIONS
With the gained insight on how to optimise spin squeezing, we now investigate how the spin squeezing influences the spatial correlations in our multimode system which is a quantum simulator for the sine-Gordon field theory.In a decoupled DW, number squeezing prolongs the global phase coherence by reducing the relative phase diffusion [10,33].Indeed, we find a good qualitative agreement between experimentally extracted global phase diffusion rates and the levels of global number squeezing (see Ext. Data Fig. 6).
The spatial resolution of our imaging system grants direct access to the local relative phase (see Appendix).This allows us to explore how the observed number squeezing impacts the local dephasing [26,45].To study a.
FIG. 7. Influence of number squeezing on multimode phase correlation function.a. correlation function PCF ⟨cos θ(z, z ′ )⟩ between two spatial positions z and z ′ along the condensates.Upper and lower panels show PCF of a state with 0.31(3) (orange) and ξ 2 N = 0.44(4) (green) and at time t after splitting to DW at A = 0.6.b.Spatially averaged PCF, ⟨cos θ(z)⟩, with z = |z − z ′ | from a. visualises how number squeezing suppresses decay of PCF.c. Evolution of ∆θ(z, −z) for z = 8µm (big circle) and z = 24 µm (small circle).We infer the dephasing rate ∂t∆θ(z, −z) with a linear fit and plot the extracted rates in d (see Ext. Data Fig. 7 for all distances z).The shaded regions represent the extracted global phase diffusion rate.We observe a spatial dependence of ∂t∆θ(z, −z) along the condensates, which hints at local number squeezing originating from the multimode dynamics of quasicondensates.
In Fig. 7a, we compare the two-point PCF in the decoupled DW at two time instances, t = 0 ms and t = 4 ms, after initial preparation with ξ 2 N = 0.31(3) (upper panels) and 0.44(4) (lower panels).Despite higher PCF in the beginning, the decay is faster with weaker number squeezing (lower panels).For better visualisation, we plot in Fig. 7b  For prethermalised states [27], it was found that the local number squeezing is directly linked to the effective temperature T − of the implemented sine-Gordon model according to the relation T − ∝ ξ 2 N (see [46] and Appendix); thus for decoupled DWs enhanced squeezing leads to a larger phase coherence length λ T ∝ 1/T − .
In multimode systems, different spatial modes can in principle feature different levels of squeezing; this would result in mode-dependent effective temperatures as for example in a generalized Gibbs ensemble [47].To examine the spatial dependence of squeezing in our experiment, we track the time evolution of the fluctuations ∆θ(z, −z) of the relative phase field between two symmetric points.We extract linear rates ∂ t ∆θ(z, −z) and show an exemplary at two distances z = 8µm and 24µm in Fig. 7c.We observe a spatial dependence of dephasing rates ∂ t ∆θ(z, −z) (see Fig. 7d).As expected, the experimental set with better global number squeezing yields an overall lower dephasing rate ∂ t ∆θ(z, −z) and additionally we observe slower dephasing at small distances.

CONCLUSION AND OUTLOOK
We have observed oscillatory dynamics of quantum fluctuations on conjugate variables of a multimode BJJ and based on this observation, developed a more efficient approach for achieving enhanced spin-squeezed states.We envision a more efficient preparation of spin-squeezed states with the help of optimal control algorithms optimizing the classical external dynamics after rapid twostep sequences.In consideration of our 1D multimode system, we have demonstrated the influence of number squeezing on prohibiting local dephasing in decoupled DWs.In the future, the ability to track the squeezing dynamics provides a new way for optimizing the preparation of strongly correlated Sine-Gordon field simulators with lower effective temperatures; by experimentally approaching a regime that is dominated by quantum fluctuations, measurements of the entanglement entropy in quantum fields will become possible [48].

Semiclassical simulation on squeezing oscillations
To visualize the squeezing oscillation in the coupled DW (see Fig. 3), we set up a semi-classical simulation.In Ext.Data.Fig. 4, we plot the equipotential lines in the classical phase space given by the two-mode Bose-Hubbard model where Λ = U N/2J signifies the interplay between the inter-atomic interaction energy and tunnel coupling energy.With suitable parameter values for DW A = 0.5: single particle tunnel coupling energy J = 41 Hz, interaction energy U = 0.33 Hz and total atom number N = 3500.Here J is estimated using the energy difference between the two lowest single-particle eigenstates in the simulated DW potential, J = (E 1 − E 0 )/2, and experimentally measured atom number N .Around the stable fixed point and given small fluctuations in the lower energy states, BH Hamiltonian can be further linearized and expressed in harmonic approxima- tion as where ∆ 2 GS Φ = √ 1 + Λ/N and ∆ 2 GS N − = N/ √ 1 + Λ are ground state fluctuations.With these parameters, we can estimate the expected ground state squeezing factors in Eq. ( 6) to be shown as shaded bands in Fig. 1b.
We sample 1000 realisations from two normal distributions (one for each observable)) with variances larger than the ground state fluctuations ∆ 2 GS (deduced from Eq. ( 6)) and propagate them with equations of motion deduced from Eq. 5 in the classical limit.We show in Ext.Data Fig. 4 the simulation result of quantum state propagation in a time span of T = 1/f p , where T corresponds to a period of Josephson oscillation of the expectation value of the obsevables (star marker).To make this more explicit, we plot the evolution of N − and Φ of a single realisation and the projected fluctuations as squeezing factor ξ 2 in each quadrature in Ext.Data Fig. 5.As one can see, the projection noise in each observable oscillates at twice the frequency as the expectation values, namely f ξ = 2f p (see also [38]).

Impact of number squeezing on global phase diffusion
In decoupled trap (J ≈ 0), phase diffusion after symmetric splitting [49,50] can be expressed as where is the phase diffusion rate and ∆ 2 Φ 0 is the initial variance of Φ right after splitting and µ(N ) is the chemical potential of BEC with atom number N .Eq. ( 8) implies slower phase diffusion with stronger number squeezing.
We investigate experimentally the influence of number squeezing on global phase diffusion rate by splitting into effectively decoupled DW, A = 0.6.For different split speed κ = δA/δt we measure the phase spread ∆Φ(t) and deduce the phase diffusion rate from a linear fit.The extracted rates ∂ t ∆Φ match the trend of the measured ξ N .(see Ext. Data.Fig. 6).In previous works studying the phenomenon of prethermalisation [26,27], the effective temperature T − between the two condensates was connected to the relative density fluctuations ∆ 2 ρ − = ⟨ρ 2 − ⟩ by

Impact of number squeezing on spatial dephasing
where ρ 0 is the peak atomic density in each condensate after splitting and ξ 2 ρ is the local number squeezing factor.

FIG. 1 .
FIG. 1. System and readout with spatial resolution Two 1D quasi-BECs with locally fluctuating quantum phases trapped magnetically below an Atom Chip with RF dressing technique.The local relative phase ϕ(z) and relative atom density ρ−(z) between the two BECs are detected with fluorescence imaging after a long Time-of-Flight.The spatial resolution allows us to probe the spatial phase correlations along the condensates.We integrate over the longitudinal direction to obtain the global observables, namely the relative phase Φ and relative atom number N−.

FIG. 2 .
FIG. 2. Preparation and dynamics of spin-squeezed quasi-BECs in the Josephson regime a. Measured squeezing factors in N− and Φ quadrature, denoted as ξ 2N and ξ 2 Φ , right after a linear ramp-up at a constant ramp speed to various tunnelcoupled DWs, denoted with the rescaled RF dressing amplitude A. In the upper x-axis, the parameter 1/N ≪ U/2J ≪ N indicates that we are in the Josephson regime of BJJ.The coloured bands mark the expected ground state squeezing factors of BJJ at zero temperature (Eq.(6)) with atom number N ∈ [2, 5] × 10 3 .b. Schematics of state evolution in the classical phase space.Left column illustrates the expected ground state fluctuation of BJJ according to Eq. (7), where orange (green) arrows represent projection noise in the relative number (phase) quadrature.The right column represents the out-of-equilibrium quantum state evolution.The triangle markers indicate single classical realizations orbiting along the equipotential lines at the plasma frequency fp.In the Josephson regime, the distribution of out-of-equilibrium states rotates and deforms which results in squeezing oscillations.A full period of the quantum state evolution around the stable fixed point (n = N−/N = 0 and Φ = 0) conforms to a half period of evolution of single realisations (mean values) in the phase space.
FIG.3.Squeezing oscillations in a coupled double well.a. Dynamics of the fluctuations of relative phase (green) and relative number (orange markers) for variable hold time in the coupled trap after the splitting.We observe oscillatory dynamics with comparable frequency and relative phase shift of π.The oscillatory dynamics is not a simple linear rotation of the state in phase space but results from the interplay of tunnel-coupling and onsite interaction and leads to a rotation plus deformation of the state (see inset).From a sinusoidal fit (solid lines) we extract the frequencies which are roughly twice the plasma frequency as expected from the rotation and shearing of the distribution in the phase space (see inset and Fig 2b).Bands indicate 68% prediction confidence interval of the fits and errorbars represent one s.e. m.

FIG. 5 .
FIG. 5. Tuning squeezing oscillation frequency a. Observed number squeezing factor ξ 2 N evolution in coupled DWs after a constant splitting speed with increasing A (bottom to top).The solid line is a fitted sine function with 68% simultaneous prediction bounds.b.The extracted squeezing frequency f ξ (diamond) from a. together with the calculated prediction 2fp (grey shade) with N ∈ [2, 5]•10 3 .c. Fitted initial phases of squeezing oscillations in a. compared to the expected accumulated phases (grey band) of the quantum state during the linear ramp with constant ramp speed to different DWs at A with atom number N ∈ [2, 5] • 10 3 as in b.Here the band uses the fitted phase of squeezing oscillation at A = 0.48 (in a) as the initial value.d.Dependence of f ξ (diamond) on total atom number N in coupled DW at A = 0.5, and for comparison experimentally measured plasma frequency 2fp (circle) and solid line marks the inferred 2fp from Eq. (2).

FIG. 6 .
FIG. 6. Driving number squeezing oscillation.a. Intercondensate distance d in 0.5 trap resulting from splitting quench with κ = 0.085 ms −1 into A = 0.5 trap.d oscillates at the transverse trap frequency fx.With this effective periodic modulation of tunnel coupling, we observe in b. that ξ 2 N (orange circle) is driven at the trap frequency fx > fp.c.Transferred squeezing oscillation after a two-step quench.τ indicates hold time in coupled DW.
the averaged PCF, ⟨cos θ(z)⟩, where z = |z − z ′ |, in the central region z = [−36, 36] µm.It is evident that enhanced global number squeezing slows down the decay of PCF over large spatial separations z.

Ext. Data Fig. 3 .
Evolution of expectation value of relative imbalance ⟨n⟩ = ⟨N−/N ⟩ in coupled DW.Grey markers, ⟨np⟩, show Josephson oscillation with an imprinted nonzero initial imbalance.Orange markers, ⟨n ξ ⟩, indicate that the expectation values are at equilibrium after symmetrical splitting.Ext.Data Fig. 4. Propagation of imprinted initial fluctuations in the two-mode BH model.The star marker signifies the evolution of a single realisation, representing the mean field value.A π rotation of a single realisation corresponds to a 2π rotation of the phase space fluctuations.Here T is one period of Josephson oscillation.

. 5 .
Squeezing oscillation frequency Evolution of a single realisation in Ext.Data Fig. 4 in the Φ (in a) and N− quadrature (in b).Evolution of quantum fluctuation in the conjugate quadratures (in c) and its fitted frequency corresponding to twice the plasma frequency fp.

Ext. Data Fig. 6 .
Global number squeezing suppresses global phase diffusion in decoupled trap.Left panel: phase diffusion in A = 0.6 trap for different splitting speeds κ.Right panel: The linearly fitted phase diffusion rates ∂t∆Φ agree qualitatively with the measured ξN with different κ.

Ext. Data Fig. 7 .
Extraction of local dephasing rates in decoupled trap in relation to global number squeezing.Evolution of spatial relative phase fluctuations ∆θ(z, −z) between two symmetrically located points around the longitudinal center of the BEC.Upper panel: ∆θ(z, −z) evolution with global number squeezing ξ 2 N = 0.31(3) and a linear fit on the four first time instances.Lower panel: ∆θ(z, −z) evolution with global number squeezing ξ 2 N = 0.44(4) and a linear fit on the three first time instances.