Phase sensitivity of spatially broadband high-gain SU(1,1) interferometers

Nonlinear interferometers are promising tools for quantum metrology, as they are characterized by an improved phase sensitivity scaling compared to linear interferometers operating with classical light. However, the multimodeness of the light generated in these interferometers results in the destruction of their phase sensitivity, requiring advanced interferometric configurations for multimode light. Moreover, in contrast to the single-mode case, time-ordering effects play an important role for the high-gain regime in the multimode scenario and must be taken into account for a correct estimation of the phase sensitivity. In this work, we present a theoretical description of spatially multimode SU(1,1) interferometers operating at low and high parametric gains. Our approach is based on a step-by-step solution of a system of integro-differential equations for each nonlinear interaction region. We focus on interferometers with diffraction compensation, where focusing elements such as a parabolic mirror are used to compensate for the divergence of the light. We investigate plane-wave and Gaussian pumping and show that for any parametric gain, there exists a region of phases for which the phase sensitivity surpasses the standard shot-noise scaling and discuss the regimes where it approaches the Heisenberg scale. Finally, we arrive at insightful analytical expressions for the phase sensitivity that are valid for both low and high parametric gain and demonstrate how it depends on the number of spatial modes of the system.


I. INTRODUCTION
In recent years, nonlinear SU(1, 1) interferometers have become an important subject in quantum metrology [1][2][3][4][5][6] as they can beat the shot noise limit (SNL), which defines the best phase sensitivity that can be achieved using classical light [1,2].Compared to linear interferometers, such as the traditional Mach-Zehnder interferometer operating with coherent light and bounded by SNL in estimating phase sensitivity, they employ nonlinear processes, namely, parametric down-conversion (PDC) or four-wave mixing (FWM) to create quantum states of light that allow for achieving phase sensitivities below the SNL [1].
A general sketch of the SU(1, 1) interferometer with two nonlinear crystals where the PDC process takes place is presented in Fig. 1.Drawing an analogy with a linear interferometer, one can determine the SNL for a nonlinear SU(1, 1) interferometer using the integral intensity of the light passing through the phase object placed inside the interferometer [2,3] (the light generated by the first crystal): where is the total number of photons generated by the first crystal and N (1) (q) is the corresponding photon number distribution over the transverse momentum q.Throughout this paper, we will use the superscript (1) to refer to quantities related to the first crystal and, analogously, (2)  for the second crystal.
Using nonlinear interferometers, it is possible to surpass the shot-noise scaling defined by Eq. (1.1) and achieve the Heisenberg limit or Heisenberg scaling given by [1,4,7,8] ∆φ which beats the shot-noise scaling for large intensities of light N (1) tot [8].Most of the early studies focuses on the theoretical description of single-mode and two-mode interferometers [1,3,4].In the single-mode regime, the phase sensitivity of nonlinear interferometers can be improved by increasing the parametric gain of the nonlinear processes and by unbalancing their gains [1,5].However, in general, PDC and FWM couple many plane-wave modes and result in the generation of multimode light [1,9,10], which requires a proper engineering of nonlinear interferometers based on such multimode sources.Indeed, in the multimode case, to achieve a phase sensitivity below SNL, an appropriate dispersion compensation in the frequency domain [6,10] or diffraction compensation in the spatial domain [3] must be performed.At the same time, in the multimode scenario, increasing the parametric gain brings time-ordering effects into play [11,12] requiring a proper theoretical description of nonlinear interferometers at high gain.
In this work, we theoretically study the phase sensitivity of spatially multimode SU(1, 1) interferometers at low and high parametric gains.Our approach is based on the solution of the system of integro-differential equations for the plane-wave operators.We investigate various interferometric configurations and pump widths in order to show how the multimodeness of the light affects the phase sensitivity.
For the multimode scenario discussed in this paper, we consider the phase sensitivity ∆φ based on the output integral intensity of the interferometer.The phase sensitivity is usually defined via the error propagation relation as [1,10,13]: where Ntot , analogously to Eq. (1.2), is the output integral light intensity of the interferometer and the covariance is given by cov(q, q ′ ) = N (q) N (q ′ ) − N (q) N (q ′ ) . (1.5) In the following discussion, we will consider the phase sensitivity normalized with respect to the shot noise limit: This paper is organized as follows.Sec.II presents our theoretical approach for the description of high-gain SU(1, 1) interferometers based on the integro-differential equations and explains the concept of compensated SU(1, 1) interferometers.Throughout this paper, we will assume that the signal and idler photons are distinguishable in some degree of freedom.However, as is shown in Appendix A, this assumption does not affect the presented results.Additional details regarding the approach for the numerical solution of the integrodifferential equations are provided in Appendix B. In Sec.III, we first derive analytical solutions of the integrodifferential equations in the limit of a plane-wave pump.Using these equations, the behavior of the optimal phase sensitivity is analyzed for both compensated and noncompensated SU (1,1) interferometers.This discussion is then extended in Sec.IV to a finite-width Gaussian pump.Here, we additionally obtain intriguing connections between the transfer functions describing the PDC process in each of both crystals and the Schmidt modes of such squeezers (see also Appendix E).In Sec.V, we then continue with a comparison of the optimal phase sensitivities of plane-wave and finite-width pumping with respect to the parametric gain, and investigate the width of the phase range for which the supersensitivity is achieved.Finally, we draw our conclusion in Sec.VI.

II. HIGH-GAIN SU(1,1) INTERFEROMETER
In the Heisenberg representation, for both low and high parametric gain, the PDC process can be described by solving a set of coupled integro-differential equations for the signal/idler plane-wave operators âs /â i [11].In the case of a Gaussian pump E p (x, z, t) = E 0 e − x 2 2σ 2 e i(kpz−ωpt) , these equations are given by [11]: where L is the integration variable, that is, the Cartesian coordinate axis in the collinear (longitudinal) direction parallel to the pump radiation.Γ is the theoretical gain parameter proportional to the field amplitudes, pump width σ and the nonlinear susceptibility; q s/i are the transverse wavevector s of the signal/idler photons and h(q s , q i , L) is a function describing the phase matching of the PDC process.Note that the full width at half maximum (FWHM) of the intensity distribution of the Gaussian pump is given by 2 √ ln 2σ.
The plane-wave operators obey the bosonic commutation relations: Note that the last commutation relation implies that the signal and idler photons are distinguishable with respect to some degree of freedom (frequency, polarization).More precisely, we will assume that the frequencies of the signal and idler photons are not identical but sufficiently close to each other so that their refractive indices are coincide.
With this assumption, the solution to this set of integro-differential equations has the form [11]: âs (q s , L, ω s ) = âs (q s ) + dq ′ s η(q s , q ′ s , L) âs (q ′ s ) where âs (q s ) and âi (q i ) are the initial signal and idler plane-wave operators, while η and β are the complexvalued transfer functions depending on the transverse wavevectors and the crystal length.
Note that this form of the solution with only two transfer functions η and β requires that the phase matching function is symmetric with respect to the wave vectors, namely, h(q s , q i , L) = h(q i , q s , L), which is fulfilled for the type-I PDC process considered in this paper since, as we mentioned above, the refractive indices of the signal and idler photons are identical.The general solution to the system of integro-differential equations (2.1) requires four transfer functions, so that the functions appearing in Eq. (2.3b) are no longer the complex conjugates of the functions appearing in Eq. (2.3a) [12,14].
Plugging the solution in the form of Eqs.(2.3a) and (2.3b) into the integro-differential equations (2.1) yields two equivalent 1 sets of coupled integro-differential equations for η, β and their complex conjugates, one of which reads × h(q s , q i , L) η * (q i , q ′ i , L) , where η * (q i , q ′ i , L) This system of equations is easier to solve numerically than the one given by Eqs.(2.1a) and (2.1b) because it is no longer operator-valued.Note that this system is also similar to the ones already derived in Refs.[12,14] for the frequency domain.Using Eq. (2.3a), the mean photon number distribution (intensity distribution) of the signal photon can be expressed in terms of the transfer function β [11]: 1 This is due to the fact that the input and the output operators are only connected via two distinct functions and their complex conjugates.

2
+ δ(q s − q ′ s ) Ns (q s ) . (2.7) Note that this expression does not contain a signal-idler cross-correlation term, since we assumed that the signal and idler photons are distinguishable, see Eq. (2.2c).The second term, proportional to the intensity distribution, can be identified as the shot noise term and follows from the Dirac-delta commutation relation of the operators due to the light energy quantization, while the argument of the modulus squared corresponds to the field amplitude (first-order) correlation function G (1) (q s , q ′ s ) [15].In the completely degenerate case, where the photons of the signal and idler beam are indistinguishable, Γ in the system of the integro-differential equations (2.1) must be replaced with 2Γ due to the appearance of an additional term in the equations during the evaluation of the commutators in the Heisenberg equations (see the derivation in Ref. [11]).However, as will be seen later, this replacement does not affect the experimentally relevant parametric gain G defined in Secs.III and IV.Hence, the transfer functions η and β will remain unchanged, given the same value of G.
Following from that, as shown in Appendix A, the integral covariance is increased by a factor of 2 due to the appearance of the cross-correlation term.Additionally, due to the photon indistinguishability, the total intensity spectrum should be considered instead of the signal beam intensity, which leads to an increase in the integrated intensity by a factor of 2: N (q) = Ns (q) + Ni (q) = 2 Ns (q) .Ultimately however, the normalized phase sensitivity f will remain unchanged because these additional factors of 2 cancel each other out, see Eqs. (1.4) and (1.6).
To describe the entire SU(1, 1) interferometer, we solve the systems of integro-differential equations separately for the first and for the second crystals with functions h (1) (q s , q i , L) and h (2) (q s , q i , L), respectively, taking into account that the output operators of the first crystal are the input operators for the second crystal.In both cases the system is solved for the region L ∈ [0, L 1 ], where L 1 is the (single) crystal length.Further details are given in Appendix B.

A. Non-compensated interferometer
An SU(1, 1) interferometer in its conventional form (without any focusing optical elements) is presented in Fig. 2(a).It consists of two PDC sections, for example nonlinear crystals, separated by some spatial region in which the measured object is placed and the pump signal and idler radiation acquire some phase φ.This phase can be simply induced by an air gap between the two crystals.However, more generally, electro-optical modu-lation or any kind of linear material between the crystals can be used to induce such a phase shift.
Moreover, in a general description, the phase φ(q s , q i ) might depend on the transverse wave vectors q s and q i [for example, for large distances between the crystals in Fig. 2(a)].However, for simplicity, we assume that the distance between two crystals is small enough, so that we can restrict ourselves to phase shifts given by an effective constant phase φ.This allows us to provide a direct comparison with the compensated interferometer and to obtain analytical expressions for the phase sensitivity in the finite-width pump case, see Sec.IV B and Appendix B.
In the multimode case, due to the divergence of the light, the radiation generated before (in the first crystal) and after (in the second crystal) the phase shift does not completely overlap, even if the distance between the crystals approaches zero, resulting in imperfect interference, which can be regarded as an internal loss mechanism.In the following, we will refer to such an interferometer as a non-compensated SU(1, 1) interferometer.
For the non-compensated interferometer, the functions h that determine the dynamics of the field operators of the first and the second crystal are then given by h (1) (q s , q i , L) = e i∆k(qs,qi)L , (2.8a) respectively, where ∆k(q s , q i is the wavevector mismatch inside the PDC section and we suppose that the distance between the crystals is d = 0.

B. Compensated interferometer
To compensate for this divergence, various focusing optical elements can be used, such as a spherical mirror as shown in Fig. 2(b) or the 4f-optical system of lenses [16].Indeed, such focusing optical elements change the wavefront of light and add a quadratic phase depending on the transverse coordinate x, namely, e −i x 2 2F , where F is the focal length [16].Therefore, if the crystal is placed at the 2F position with respect to the spherical mirror (or two crystals placed at the F positions with respect to the 4f-optical system of lenses), the quadratic phase introduced by the mirror compensates for the quadratic phase of light, which leads to a change in the wavefront from convex to concave.In this case, taking into account the phase compensation and the wavefront modification, the h-function for the second crystal is given by h (2) (q s , q i , L) = e −i∆k(qs,qi)[L−L1] e iφ . (2.9) We define this kind of SU(1, 1) interferometer with divergence compensation as a compensated SU(1, 1) interferometer.

III. PLANE-WAVE PUMP
We start our analysis by considering the plane-wave pump case, where analytical expressions for the output plane-wave operators can be obtained.Formally, the transition to the plane-wave pump case can be performed by taking the limit σ → ∞ in the system of integrodifferential equations (2.1): The resulting delta-function eliminates the integrals in the integro-differential equations (2.1) and allows for an analytical treatment of the problem.More rigorously, the following set of coupled differential equations can be obtained by repeating the derivation presented in Ref. [11] for a plane-wave pump: where h(q s , L) def.
= h(q s , −q s , L).Note that in this case, each signal transverse momentum q s is connected with only one fixed idler transverse momentum, q i = −q s .The commutation relations for the plane-wave operators are given by Eqs.(2.2a)-(2.2c).
Similarly to the finite-width case presented in Sec.II, the solutions to the system of differential equations (3.2) can be written in the form [17][18][19][20]: âs (q s , L, ω s ) = ηpw (q s , L) âs (q s ) + β pw (q s , L) â † i (−q s ) , Again, this form of the solution requires h(q s , −q s ) = h(−q s , q s ), or, more specifically, h(q s ) = h(−q s ) and leads to two equivalent sets of differential equations for the β pw and η * pw functions, one of which reads: dη * pw (−q s , L) dL = Γ 0 h * (q s , L) β pw (q s , L) .
(3.4b) Using Eq. (3.3a), the distribution of the mean number of signal photons is given by Ns (q s ) = |β pw (q s , L)| 2 δ(0) . (3.5) The divergent factor δ(0) appears due to the Dirac-delta commutation relations for the plane-wave operators and In the non-compensated case, the radiation generated in the first crystal diverges, which can be considered as an internal loss mechanism since the spatial overlap of the radiation generated in two crystals is reduced, as indicated by the black arrows.Importantly, note that even if the distance between the crystals vanishes, such an overlap is still not perfect (the radiation also diverges inside the crystals).
In the compensated case, the spherical mirror reflects the radiation back into the crystal and compensates for the divergence.
the infinite transverse size of the system [20].Below, to avoid this divergence, we instead consider the density of quantities, namely, the intensity and covariance per transverse length of the system, which are well-defined even in the limit of L x → ∞.
For example, the density of the intensity spectrum of the signal radiation is given by Similarly, the covariance density of the signal radiation can be expressed via the intensity as (see also Ref. [21]): where A more detailed discussion of this divergence problem and the derivation of Eqs.(3.6), (3.7a) and (3.7b) is shown in Appendix D. Solving the coupled differential equations (3.4) for a single crystal with Eq. (2.8a) and the initial conditions ηpw (q s ) = 1 and β pw (q s ) = 0 yields the transfer functions connecting the input and output operators [21]: η( 1) pw (q s , L) = cosh where and ∆k(q s ) def.
Then, the density of the intensity distribution of the signal beam after the first crystal is given by: Note that although g(q s ) can take both real and purely imaginary values, the intensity distribution is always real and positive due to the sinh-term.Therefore, it is not necessary to write | • | 2 on the right hand side of Eq. (3.10).Note that Eqs.(3.8)-(3.10)derived above coincide with the results already found in a similar fashion in Refs.[17][18][19]21].
To analyze the effect of the focusing element (spherical mirror) on the output spectra and the phase sensitivity, we first start by extending the plane-wave-pump analytical treatment to the SU(1, 1) interferometer without any compensation elements and then compare the noncompensated scheme with its compensated counterpart.
In this work, we consider an SU(1, 1) interferometer consisting of BBO crystals of length L 1 = 2 mm, pumped by a laser with a wavelength of 354.6 nm.To obtain a connection between the theoretical gain parameter Γ 0 and the experimental gain, the collinear output inten-sity N (1) s (0) of a single crystal is fitted by the function y(Γ 0 ) = B sinh 2 (AΓ 0 ).The experimental gain G is then given by G = AΓ 0 [11].For plane-wave pumping and perfect phase matching in the collinear direction [that is, ∆k(0) = 0], it is immediately clear from Eq. (3.10) that A = L 1 and B = 1.
If instead of Γ 0 , the prefactor in the system of coupled differential equations (3.4) is cΓ 0 def.

= Γ ′
0 , where c > 0, the fitting constant would still be given by A = L 1 for a fit with respect to Γ ′ 0 and the experimental gain would be defined as G = AΓ ′ 0 .The fitting constant B would be given by c 2 .Therefore, the definition of the experimental gain is independent on the scaling factor c. For the degenerate case mentioned in Sec.II, the scaling factor is c = 2. Therefore, all results presented in this section would be the same for the degenerate case where the signal and idler photons are indistinguishable.

A. Non-compensated scheme
In order to find a solution in the case of the noncompensated SU(1, 1) interferometer, we analytically solve the system of differential equations (3.2) using Eqs.(2.8a) and (2.8b), and obtain the output plane-wave operators at the output of the SU(1, 1) interferometer.In this case, the form of the solution solution (3.3) has the following transfer functions: These functions allow us to calculate the density of the output intensity distribution of the signal beam: (3.12) The covariance density for the signal beam can be found by plugging Eq. (3.12) into Eq.(3.7b) with the use of Eq. (3.7a).
Figs. 3(a)-3(c) present the intensity distributions for different gains and phases between the crystals according to Eq. (3.12).They are plotted over the external angle θ s , which, in the case of small angles, is connected to the transverse wavevector via θ s ≈ q s /k vac s , where k vac s is the wavevector of the signal photons in vacuum.Note that throughout this work, we will use the external angles instead of the transverse wave vectors for all relevant plots.
One can observe that the intensity profiles broaden as the gain increases.The phase between the crystals strongly modifies the intensity profiles, leading to destructive interference for a certain range of angles.This is due to the fact that the radiation generated in each crystal has a quadratic phase with respect to the angle, see Eqs. (3.11a) and (3.11b).Even if the distance between the crystals is zero, such a quadratic phase is present and leads to a modification of the intensity distribution when various additional constant phases are applied.The dis-tributions of the covariance have a similar behavior to the intensity profiles and are shown in Appendix F 1. The normalized phase sensitivity is presented in Fig. 4(a).It can be seen that the sensitivity is destroyed as the gain increases.This is directly related to the fact that the quadratic phase of the radiation entering the second crystal is not compensated, that effectively acts as internal losses and destroys the phase sensitivity.With increasing gain, such losses become more pronounced.A similar behavior of the phase sensitivity of the frequency multimode SU(1, 1) interferometer was observed in Refs.[6,10].

B. Compensated scheme
In the case of the compensated SU(1, 1) interferometer, we analytically solve the system of differential equations (3.2) with the use of Eq. (2.8a) and Eq.(2.9), and find the following transfer functions which connect the output and the input plane-wave operators of the entire interferometer: One can observe that in this case, both the functions β pw (q s ) and ηpw (q s ) have a simple dependence on the phase φ, which allows the output intensity to be written in a much more compact and simple form compared to Eq. (3.12): where ξ pw (q s ) = N (1)  s (q s ) 1 + N (1) s (q s ) . (3.14b) Note that the entire intensity spectrum now scales as ∼ cos 2 (φ/2) and is therefore identically zero for φ = π, which implies perfect destructive interference for all angles of emission.
The normalized intensity distribution for different gains is shown in Fig. 3(d).For the chosen interval of gains, the width of the distribution remains almost unchanged.However, since the distribution of a single crystal broadens with increasing gain, then, according to Eqs. (3.14a) and (3.14b), the intensity distribution of the entire interferometer will also be broadened for higher values of gain.
Substituting the expression for the intensity distribution given by Eqs. where is the density of the integral intensity of the SU(1, 1) interferometer and Note that in the plane-wave pump case, the phase sensitivity defined via the integral quantities [Eq.(1.4)] can be obtained by dividing the phase sensitivity density [Eq.(3.15a)] by √ L x , which means that the integral phase sensitivity is length-dependent.However, the SNL-normalized phase sensitivity f , given by Eq. (1.6), no longer depends on the considered transverse size of the system and coincides with the normalized phase sensitivity density.
The SNL-normalized phase sensitivity f is presented in Fig. 4(b).In this case, due to the full compensation, all spatial modes are deamplified simultaneously, therefore, the phase sensitivity behaves similarly to the singlemode case [1] and beats the shot noise limit.Moreover, from Fig. 4(b) it becomes apparent that the width ∆ of the phase range for which the SNL is overcome (f < 1) gradually narrows as the parametric gain increases.Furthermore, this region is defined for all gains, which means that for all gain values (even for low gains), there exists some range of phases where the phase sensitivity beats the SNL, which is similar to early studies of the singlemode interferometer [7].
From Eqs. (3.15a)-(3.15c)it can be easily seen that the best phase sensitivity (minimal values of ∆φ and f ) is achieved for φ = π, that is, at the dark fringe of the interferometer.Similarly, the optimal working point for more elaborate models (including losses) was found to be near the dark fringe [1,4,22].The value of the minimized SNL-normalized phase sensitivity is then given by: where N s,2 = dq s N (1) s (q s ) 2 . (3.16c) Notably, the expressions derived above only depend on the intensity spectra of the first crystal.This is due to the fact that both crystals have the same parameters (length, gain etc.).At the same time, the perfect compensation induces a symmetry, so that the parameters of the first crystal fully describe the system.
The behavior of the optimal phase sensitivity as a function of the parametric gain G is discussed in Sec.V where a comparison to the finite-width pump case and the Heisenberg scaling of the phase sensitivity is also drawn.

IV. FINITE-WIDTH GAUSSIAN PUMP
In this section, we extend our studies by considering a pump beam with a Gaussian profile with a FWHM of the intensity distribution of 50 µm.All other parameters are the same as in the plane-wave case, see Sec.III.For a finite-width pump, the output signal and idler planewave operators can only be found by numerical integration of the system of integro-differential equations (2.4), as described in Appendix B. From there, the intensity distributions, covariances and the normalized phase sensitivities are calculated numerically.
Analogously to the plane-wave case, the collinear intensity2 N (1) s (0) dq from a single crystal must be fitted in order to obtain a connection between the theoretical gain parameter Γ and the experimental gain G.However, in the multimode regime with finite-width pumping, this procedure is not straightforward as for the planewave pump case, especially if a large range of gain values should be considered.To resolve this, instead of a single fit, many fits are performed, capturing different ranges of the parametric gain.More details regarding this procedure and exact values of A for the results presented in the following sections are given in Appendix C. Similarly to the plane-wave pump (see Sec. III), the definition of the experimental gain is independent of any additional scaling factors c in the system of integrodifferential equations (2.1), since the theoretical gain constant can be redefined to absorb such factors.The remaining integro-differential equations are identical to the ones without the scaling factor c, so making a fit with respect to the newly defined theoretical gain constant Γ ′ def.
= cΓ yields the same A. The definition of G is then given by with respect to the scaled gain constant G = AΓ ′ , and therefore remains unchanged.

A. Non-compensated scheme
In the case of the non-compensated scheme, the profiles of the intensity distribution for different phases and gains are presented in Figs.5(a)-5(c).One can observe that the Gaussian pump profile modifies the intensity at the zero phase, bringing it closer to the Gaussian shape.Overall, for different phases, the graphs show a similar behavior to the plane-wave pump case.However, since the number of modes is smaller in the case of a finite-width pump compared to the case of a planewave pump, the intensity distribution in the case of finite width pumping broadens more slowly.The plots of the covariance are presented and discussed in Appendix F 2. Fig. 6(a) presents profiles of the SNL-normalized phase sensitivity for different gains where no supersensitivity regions are observed.

B. Compensated scheme
For the compensated configuration, similarly to the plane-wave pump, it is possible to obtain simplified expressions for the intensity, the covariance and the phase sensitivity.To do this, relationships between the transfer functions of the first and the second crystal are required.In the case of the perfect compensation, these relations are derived in Appendix E and are given by: η(2) φ (q, q ′ ) = η(1) (q ′ , q) * , (4.1a) φ (q, q ′ ) = e iφ β (1) (q ′ , q) .(4.1b) To be precise, in the above expression, the transfer functions with the index (1) (index (2) ) connect the input and the output plane-wave operators of the first (second) crystal according to Eqs. (2.3a) and (2.3b).The index φ indicates the dependence on the interferometer phase for the second crystal 3 .
In Appendix B, we derive expressions for the composite transfer functions of the entire interferometer (η (SU) and β (SU) ) in terms of the transfer functions of both sep-arate crystals, see Eqs. (B.13a) and (B.13b), which can be used to obtain a simplified expression for the intensity spectrum of the entire interferometer via Eq.(2.6): where ξ(q s , q ′ ) = dq β (1) (q, q s ) η(1) (q, q ′ ) * . (4.2b) Importantly, note that the function ξ(q s , q ′ ) no longer depends on φ.This implies that, as in the plane-wave As mentioned in the text, the minimal value of f is always reached at φ = π.Furthermore, f diverges for φ → 2πk, k ∈ Z, see Eq. (4.7).In both cases, the pump is given by a Gaussian pump with an intensity-FWHM of 50 µm.
pump case, all spatial modes are amplified or deamplified simultaneously, and the intensity distribution is identically zero for φ = π.Furthermore, the expression for the output intensity depends only on the transfer functions of the first crystal.This is due to the fact that the full compensation induces a symmetry in the system, see Appendix E.
Plots of the intensity distributions for several parametric gain values can be found in Fig. 5(d).For the chosen gain parameters, the width of the intensity distribution remains almost unchanged, however, as the gain increases, its shape becomes more Gaussian.This is associated with a decrease in the number of modes with increasing gain.
Further simplifications are possible by considering the Schmidt decomposition of the transfer functions.As it was shown in Refs.[12,14], there exists a joint Schmidt decomposition (Bloch-Messiah reduction) for the transfer functions β and η, which reads where √ Λ n and Λn are the singular values of β and η, respectively, u n (q) are the functions associated with the output Schmidt operators (we denote these functions as the output Schmidt modes) and ψ n (q ′ ) are associated with the input Schmidt operators (the input Schmidt modes) of the considered system.See also Appendix A 1 for further details.Eqs.(4.3a) and (4.3b) extend the results for the Schmidt decomposition already found for the frequency domain in Refs.[12,14] to the spatial domain.Note that the eigenvalues Λ n and Λn are connected via Λn = 1 + Λ n [12].
Applying Eqs.(4.3a) and (4.3b) to Eqs. (4.2a) and (4.2b), the total intensity of the SU(1, 1) interferometer can be written as where Similarly, the covariance can be written as cov(q s , q ′ s ) = 16 cos 4 φ 2 dq |ξ(q s , q) ξ * (q ′ s , q)| 2 + δ(q s − q ′ s ) Ns (q s ) .The integral covariance can be written in the following form: where Finally, combining the formulas for the integral intensity and covariance presented in this section, it is possible to derive an analytic expression for the phase sensitivity in the case of perfect compensation and finite-width pumping: Note that the expression above has a similar form to the phase sensitivity in terms of the density of quantities presented in Eqs.(3.15a)-(3.15c)for the plane-wave case.Profiles of the SNL-normalized phase sensitivity f for different gains are shown in Fig. 6(b).One can observe that as the gain increases, the phase supersensitivity range width ∆ (the range of phases for which the phase sensitivity beats the shot noise limit) gradually decreases.However, for all presented gains, there exists some region where the phase sensitivity beats the shot noise level.From Eq. (4.7) it can be seen that, similarly to the plane-wave pump case, the best (optimal) phase sensitivity is achieved for the optimal point φ = π, that is, at the dark fringe of the interferometer (see also Refs.[1,4,22]).For this phase, the SNL-normalized phase sensitivity is given by A more detailed analysis of the behavior of f min and ∆ as a function of G is given in the next Sec.V, where we also draw a comparison to the plane-wave pump case and the Heisenberg scaling of the phase sensitivity.

V. COMPARISON OF THE OPTIMAL PHASE SENSITIVITIES
In this section, to compare the plane-wave and finitewidth pumping regimes, we only consider the compensated setup, because the non-compensated setup shows no supersensitivity for any phase, as seen in the previous sections.
As it was discussed in Sec.III, in the plane-wave pump case, the integral quantities depend on the transverse length L x of the system.Therefore, to ensure a proper comparison of the two cases, we assume that the planewave pump system has the required transverse length L x , so that its integral number of photons after the first crystal is equal to the integral number of photons in the case of finite-width pumping for the same value of the Optimal phase sensitivity f min FIG.7. Optimal phase sensitivities fmin of the compensated SU(1, 1) interferometer for the plane-wave pump and the finite-width Gaussian pump with a FWHM of the intensity distribution of 50 µm versus the parametric gain G. fH corresponds to the expected asymptotic Heisenberg scaling normalized with respect to the SNL, see Eq. (5.2).Lx (right-hand vertical axis) is the transverse length of the system pumped by the plane-wave pump to reach the same output intensity after the first crystal as in the case of the system pumped by the finite-width pump.Note that due to the assumed equality of intensities, the Heisenberg scalings coincide in both cases.See also Appendix D for more details on Lx. parametric gain.The dependence of L x on the gain is presented in Fig. 7.For more details see Appendix D.
The assumed equality of the light intensities after the first crystal leads to the same Heisenberg limit for both cases which, in accordance with Eq. (1.3), we define as (5.1) Thus, the SNL-normalized Heisenberg scaling is given by To compare the optimal phase sensitivity in the planewave and the finite-width cases, we show f min for both cases over a range of parametric gains in Fig. 7.One can observe that, as the parametric gain increases, the optimal phase sensitivity approaches the Heisenberg scaling for both pumps.However, the phase sensitivity for the finite width pump slightly surpasses the phase sensitivity for the plane-wave pump.The reason for this is the different number of modes in two cases, namely, the number of modes in the finite-width case is smaller compared to the plane-wave case.
To see this more clearly, we define the Schmidt number, which is a measure of the effective number of modes, for the first crystal in the finite-width case as [15,23,24]: where are the normalized eigenvalues obtained from the decomposition of the transfer function β (1) of the first crystal, see Eq. (4.3a).We can then rewrite Eq. (4.4b) as s,tot k .This allows us to express the optimal phase sensitivity via the Schmidt number as s,tot K (1)   . (5.5) From this expression it becomes clear that the phase sensitivity is optimized for high intensities (high parametric gains) and small effective mode numbers.In the extreme case of a single mode (K = 1) with an eigenvalue Λ (1) 0 , the optimal phase sensitivity is given by which approaches the Heisenberg limit f H for Λ (1) 0 ≫ 1.In the plane-wave case, due to the strong correlations between the signal and idler photons, it is no longer possible to define Schmidt modes as in the finite-width pump case.However, by analogy, we find from Eqs. (3.16a)-(3.16c): where s (q s ) 2 (5.7b) can be defined as the Schmidt number in the plane-wave pump case.Alternatively, we can write it as with the following functional defined for some function u: This provides us with an expression for K pw which has already been found and discussed in Ref. [25] for the lowgain regime.
Plots of the Schmidt numbers for plane-wave and finite-width pumping are shown in Fig. 8.As expected, the Schmidt number gradually decreases as the parametric gain increases [11].For G > 1, the Schmidt number for the finite-width pump is lower than for the planewave pump.Therefore, assuming equal intensities after the first crystal in both cases, the phase sensitivity for finite-width pumping surpasses the phase sensitivity for plane-wave pumping as it is shown in Fig. 7.
The phase supersensitivity range width ∆ as the function of parametric gain G is shown in Fig. 9 for both finite-width pumping and plane-wave pumping.Clearly, at low gains the width for both cases coincides.As the gain increases, plane-wave pumping leads to a slightly narrower phase sensitivity region.However, in both cases, ∆ drops quickly and decreases by about one order of magnitude with each increase of the gain G by two.

VI. CONCLUSION
In this work, we present a theoretical description of the high-gain multimode SU(1, 1) interferometer, which, however, is also valid for the low-gain regime.Our approach is based on the solution of the system of integro-

Supersensitivity region width ∆ [rad]
plane-wave 50 µm FIG. 9.The phase supersensitivity range width ∆ (length of the phase interval where f < 1) versus the parametric gain for the plane-wave and finite-width Gaussian pump.
differential equations for the plane-wave operators and allows us to introduce a set of gain-dependent Schmidt modes.We investigate the spatial properties and phase sensitivity of multimode nonlinear interferometers in different configurations (with and without diffraction compensation) for various pump widths.
In the case of a plane-wave pump, the system of integro-differential equations is reduced to the set of ordinary differential equations with an analytical solution.At finite-width pumping, the system of integrodifferential equations is evaluated numerically to find the corresponding transfer functions.In addition, the use of the derived relations between the transfer functions of individual crystals and the entire interferometer, allows us to obtain compact expressions for the phase sensitivity for finite-width pumping.Finally, the minimized phase sensitivity for the plane-wave and the finite-width pump is compared for different gains.
We demonstrate that for the non-compensated interferometers, the phase sensitivity worsens as the parametric gain increases and never surpasses the shot noise scaling, which is due to incomplete overlapping of the radiation of the first and the second crystals.However, for the compensated interferometers, both finite-width and plane-wave pumping result in a phase sensitivity that always exceeds SNL, so that the optimal phase sensitivity approaches the Heisenberg scaling for large parametric gains.However, the phase width of the supersensitivity region becomes increasingly narrow as the gain increases.To counteract this narrowing, as well as to optimize the phase sensitivity in the case of optical and detection losses, unbalanced SU(1, 1) interferometers [1] can be considered, which may be the subject of future research.tal â(1,out) to the vacuum plane-wave operators: where for simplicity, we have dropped the length L and the frequency dependence ω from the arguments of the operators.The index (1) indicates that these functions are the solutions for the first crystal.The initial conditions for the functions read Next, we can set up a similar set of equations for the second crystal: i which now connects the input operators of the second (not the vacuum operators) to its output operators.We have added the index φ to emphasize that the solution of this system depends on the interferometer phase φ via the function h (2) describing the second crystal, see Eqs. (2.9) and (2.8b).
Plugging this back into the system of integrodifferential equations (2.1) yields, again, the same system of integro-differential equations (2.4) with η and β replaced by η(2) φ and β φ , respectively.This system of the integro-differential equations is solved for L ∈ [0, L 1 ], where L 1 is the crystal length, and has the same initial conditions as for the first crystal: The phase matching function h (2) which describes the second crystal can be factorized as h (2) (q s , q i , L) = p (2) (q s , q i , L) e iφ , (B.5) where p (2) is a complex-valued function.
Plugging this form of h (2) into the system of integrodifferential equations (2.4) and defining (q, q ′ ) def.

= η(2)
φ (q, q ′ ) , (B.7) one can obtain the following system of integro-differential equations for the second crystal: where the initial conditions for the transfer functions are given by This system, importantly, results in a solution for η(2) and β (2) that is independent of φ.

(B.14b)
Hence, it is possible to solve the set of the integrodifferential equations independently for both crystals.Using p (2) instead of h (2) for the second crystal makes it possible to obtain a phase independent solution for the second crystal.Then, the phase independent solutions of both crystals can be combined for an arbitrary value of φ according to Eqs. (B.14a) and (B.14b) without having to solve the integro-differential equations for each value of phase.
Appendix C: Behavior of the fitting constant A with increasing gain As it was shown in Sec.III, for the plane-wave pump, the collinear intensity behaves as ∝ sinh 2 (G), where G = ΓL 1 is the parametric gain and L 1 is the crystal length.Similarly, the sinh 2 -behavior is expected for the collinear emission in the case of pumping with an arbitrary width [11].However, in general, for the multimode light, this tendency is not preserved for all gain intervals, namely, the fitting parameter is no longer a constant over the full gain range.
In order to obtain the correct gain-dependence, we therefore perform a series of fits of collinear intensities, one for each small interval of theoretical gain Γ under consideration.From this procedure we effectively obtain the fitting parameter A Γ for different values of Γ.Then, the experimental gain is given by G = A Γ Γ.We have added a subscript Γ to show that the fit parameter A now depends on Γ.Thus, the fitting formula now reads y(Γ) ∝ sinh 2 (A Γ Γ).
The values of A Γ corresponding to the experimental gains used for the plots shown in Sec.IV are given in Table I.Furthermore, we show the full behavior of A over an entire range of gain values in Fig. 10.It can be seen that in the limit of high gain values (when the radiation tends to a single mode), the fitting parameter approaches the constant value.The uncertainties indicated for A Γ and G are the uncertainties introduced by the fitting procedure.Note that since the fit parameter A Γ has an uncertainty assigned, G can also be provided only within some uncertainty.Note that the density of points is higher for G < 1.9 that results in decrease of the error below G = 1.9.See also Table I for the selected values.

Appendix D: Divergence in the plane-wave pump case
After quantizing the electromagnetic field in a cavity of length L x in terms of the discrete plane-wave operators âq obeying the commutation relations [13] [â q , â † q ′ ] = δ q,q ′ , (D.1) it is common to take the continuous limit L x → ∞, which allows the transition from the discretely spaced wave vectors with spacing ∆q = 2π Lx to a continuous spectrum (∆q → 0) [26].
Here, we have assumed a fixed frequency for the cavity modes of interest.This is why we only consider the commutator in the transverse variables, since the zcomponents of the wave vectors are strictly fixed by given q and the frequency.
As was seen in Sec.III, taking the continuous limit for the plane-wave pump approximation results in a divergence of the intensity spectrum, namely, Ns (q s ) = 0| â † s (q s ) âs (q s ) |0 ∝ δ(0).Physically, if the quantization size of the system becomes unbounded while the system is pumped by a plane wave covering all of space, an infinite amount of PDC photons is created along the transverse dimension.
To avoid this divergence and make the description more mathematically rigorous, it is possible to repeat the derivation for the integro-differential equations presented in Ref. [11] with a plane-wave pump for the discrete operators.In that case, the integrals over wave vectors have to be replaced with the sums over wave vectors with the correspondence dq ↔ ∆q.The system which generates the PDC radiation then consists of the quantization box with the periodic boundary conditions and is again described by the system of differential equations (3.4) with the form of the solution (3.3).
In this case, evaluating the expectation value in Eq. (3.5) yields

Ns (q
Clearly, the quantity defined by describes the intensity per length Lx 2π of the system.The important observation is now that N s (q s ) remains unchanged and finite when taking the continuous limit Similarly, the covariance evaluated from its definition [Eq.(1.5)] is given by cov(q s , q ′ s ) = After defining cov pw (q s , q ′ s ) def.
= ∆q cov(q s , q ′ s ) , (D.6) the continuous limit can be taken (L x → ∞) and results in the well-behaved quantity cov pw (q s , q ′ s ) = N s (q s ) [1 + N s (q s )] δ(q s − q ′ s ) , (D.7) which corresponds to the covariance per the transverse length of the system.
In order to make the results for the plane-wave and the finite-width pump comparable, we consider such a transverse size L x of the system in the case of planewave pumping, which gives the same integral intensity after the first crystal as in the finite-width pump case.Note that L x is different to the quantization length L x .Then, where Ns,tot is the integral intensity in the case of the finite-width pump and is the integral photon density (integral intensity per transverse length), analogously to Eq. (1.2).
Note that this step is not necessary for the SNLnormalized phase sensitivity f , since according to Eqs. (1.4) and (1.6), it does not depend on the factor L x .Indeed, obtaining the integral number of photons and covariance by multiplying their densities [Eqs.(D.4) and (D.7)] by L x and substituting them into Eqs.(1.4) and (1.6) shows that the additional factors L x are cancelled out.

Appendix E: Analytical implications of the compensation
More explicitly, the solutions of the integro-differential equations [Eqs.(2.3a) and (2.3b)] can be written as âs (q s , L 1 ) = dq ′ s η(q s , q ′ s ; L 1 , L 0 ) âs (q ′ s , L 0 ) which corresponds to a more general case where the transfer functions η(q, q ′ ; L 1 , L 0 ) and β(q, q ′ ; L 1 , L 0 ) connect the plane-wave operators at L 0 to the ones at L 1 .
The inverse transform of the solutions (E.1) is given by [12,27,28]: which can be immediately verified by plugging the form of the solutions (E.2) into the solutions (E.1) and using Eqs.(A.3a) and (A.3b).

(E.3b)
These have also been found in Ref. [12] for a more general PDC process in the frequency domain.
For a single crystal ranging from the coordinates L 0 = 0 (beginning of the crystal) and L 1 (end of the crystal), the system of the integro-differential equations for the inverse transform [Eqs.(E.2a) and (E.2b)] is given by: Here we have assumed a more general case than in the main text, where r(q s , q i , L) contains the pump term and might additionally depend on the integration coordinate L. Note that the integration interval for the solution of the system starts at the upper bound L 1 of the PDC section.The initial value conditions state that the planewave operators are known at L 1 .After performing the integration, the obtained functions β (1) (q ′ i , q s ; L 1 , L 0 ) and η(1) (q ′ i , q i ; L 1 , L 0 ) connect the plane-wave operators at L 1 to the ones at the beginning of the crystal at L 0 .
Next, we consider the second crystal for the compensated configuration with the same length as the first crystal ranging from L 1 to L 2 = 2L 1 .As shown in Appendix B, it is not necessary to solve the system of integro-differential equations with the phase term e iφ included in r (2) (q s , q i , L).Instead, one can first integrate the equations without the phase and then multiply β (2)  with the phase afterwards to obtain the solution for the second crystal including the phase.Furthermore, note that we will assume that r (2) (q s , q i , L) = r(q s , q i , L).As will be seen below, this is a necessary requirement to obtain compensation.The system of the integro-differential equations and the initial value conditions for the phase independent transfer functions of the second crystal read: dβ (2) (q s , q ′ i ; L, L 1 ) dL = Γ dq i r(q s , q i , L) e −i∆k(L−2L1) η(2) (q i , q ′ i ; L, L 1 ) * d η(2) (q i , q ′ i ; L, L 1 ) * dL = Γ dq s r * (q s , q i , L) e i∆k(L−2L1) β (2) (q s , q ′ i ; L, L 1 ) η(2) (q i , q ′ i ; L 1 , L 1 ) * = δ(q i − q ′ i ) β (2) (q s , q ′ i ; L 1 , L 1 ) = 0. (E.5) Assuming that r(q s , q i , 2L 1 − L) = r(q s , q i , L) , (E. 6) which means that r is mirror-symmetric around L 1 , it is possible to make the substitution L = 2L 1 − L so that after redefining η(2) (q i , q ′ i ; 2L 1 − L, L 1 ) def.
From that, we can identify the following properties relating the transfer functions of the first and the second crystal: η(2) (q, q ′ ; 2L 1 , L 1 ) * = η(1) (q ′ , q; L 1 , 0) , (E.9a) β (2) (q, q ′ ; 2L 1 , L 1 ) = β (1) (q ′ , q; L 1 , 0) .(E.9b) Note that these solutions still do not include the interferometer phase.When including it, the second relation instead reads β (2) (q, q ′ ; 2L 1 , L 1 ) = e iφ β (1) (q ′ , q; L 1 , 0) .In this work, we have assumed a perfect compensation of the quadratic phase of the light, and therefore, we have only considered an interferometer phase which has no dependence on the transverse signal and idler wave vectors q s and q i .Generally however, when the compensation is not performed, that is, when the phase is q-dependent, Eqs.(E.10), (E.9a) and (E.9b) do not hold.Nevertheless, in systems with a complicated q-dependent phase, as long as it is possible to perform the compensation using a set of lenses, mirrors or spatial light modulators (SLMs), the necessary symmetry is induced and, therefore, Eqs.(E.10), (E.9a) and (E.9b) hold.In that case however, proper attention must be given to the losses introduced by these additional optical components since they reduce the phase sensitivity [1].Figs.11(a)-11(f) present the covariance profiles for different gains in the plane-wave pump case.The profiles are described by the function C(q s ) defined in Eq. (3.7b).As one can see, the plane-wave covariance cov pw (q s , q ′ s ) [Eq. (3.7a)] vanishes unless q s = q ′ s , therefore, the function C contains all the information about the covariance.This is a direct consequence of the plane-wave pumping which results in perfect correlations between the planewave modes.As was already mentioned in the main text (Sec.III), the covariance profiles calculated for the noncompensated setup [Figs.11(

Finite-width pump
Contrary to the plane-wave pump case, the covariance profiles for a finite-width pump are in general only fully described by considering the entire range for both arguments of the covariance function cov(q s , q ′ s ).The plots for the non-compensated setup are shown in Fig. 12(a) and for the compensated setup in Fig. 12(b).Since we have to show these plots as colormesh-plots, we restrict ourselves to only showing two values for the parametric gain for each setup.Nevertheless, important features of the covariance functions become apparent.
First, it should be noted that for both setups and for the small gain G = 1.25, a sharp line can be observed along the diagonal.This feature is the shot-noise term δ(q s − q ′ s ) Ns (q s ) appearing in Eq. (2.7).Numerically, the Dirac-delta is expressed via δ(x − y) ≡ δ ix,iy /∆x, where i x and i y are the indices on the numerical grid corresponding to the points x and y, respectively, and ∆x is the step size of the grid.Therefore, in the discretized case, we observe a line of single-pixel width and finite value along the diagonal.As the parametric gain increases, this term becomes negligible [15], which results in the ellipse-shaped covariance for the larger gain value G = 3.75.However, due to the perfect destructive interference in the compensated setup, the shot noise term becomes visible again as φ → π.This is due to the fact that the term resulting in the ellipse-shaped contribution scales as ∼ cos 4 (φ/2) with the interferometer phase φ, while the shot noise term scales as ∼ cos 2 (φ/2), see Eq. (4.5) and Eqs.(4.2a) and (4.2b), respectively.
For the non-compensated setup, the covariance-ellipse narrows along the diagonal as the phase increases.Then, similar to the intensity profiles shown in Fig. 5, a splitting of the covariance is observed as φ → π which is caused by the imperfect destructive interference.For the compensated setup this is not the case.Instead, the widths of the ellipse along the diagonal and anti-diagonal stay unchanged.However, as mentioned above, the shot noise term becomes visible for phases close to π.Clearly, the ellipse shape is determined by the integral term in Eq. (4.5) and is therefore independent of φ.

FIG. 1 .
FIG.1.Simplified sketch of an SU(1, 1) interferometer consisting of two crystals.The PDC radiation (signal, idler) generated in the first crystal acquires some phase and is amplified/deamplified in the second crystal depending on the relative phase difference φ = φp − φs − φi of the pump, signal and idler radiation.

FIG. 2 .
FIG. 2. (a) Non-compensated and (b) compensated configurations of the SU(1, 1) interferometer.In the non-compensated case, the radiation generated in the first crystal diverges, which can be considered as an internal loss mechanism since the spatial overlap of the radiation generated in two crystals is reduced, as indicated by the black arrows.Importantly, note that even if the distance between the crystals vanishes, such an overlap is still not perfect (the radiation also diverges inside the crystals).In the compensated case, the spherical mirror reflects the radiation back into the crystal and compensates for the divergence.

FIG. 3 .
FIG. 3. Calculated intensity distributions for different phases between the crystals and gains in the case of plane-wave pumping.(a)-(c) Non-compensated interferometer and (d) compensated interferometer.In the compensated case, the intensity profile is only scaled as ∼ cos 2 (φ/2), see Eqs. (3.14a) and (3.14b).

) 3 FIG. 4 .
FIG. 4. (a) Calculated normalized phase sensitivity f for different gains for the non-compensated setup.For each gain, the crosses mark the points where f is minimized.The horizontal thin gray line indicates the shot noise level (f = 1).The divergences occur at points where the derivative of the integral intensity (density) vanishes: dNs,tot dφ = 0 [see Eqs.(1.4) and (D.9)].(b) Calculated normalized phase sensitivity f for different gains for the compensated setup.Note that the entire plotted region is below the shot noise level.The inset shows the full phase interval [0, 2π].As mentioned in the text, the minimal value of f is always reached at φ = π.Furthermore, f diverges for φ → 2πk, k ∈ Z, see Eqs. (3.15a)-(3.15c).In both cases, the pump is given by a plane wave.
(3.14a) and (3.14b) into Eqs.(3.7a) and (3.7b), one can calculate the covariance in the compensated case.The plots of the covariances for different interferometric phases and gains are presented in Appendix F 1 and show a similar behaviour to the intensity profiles.Finally, by replacing the intensity and the covariance in Eq. (1.4) with their densities given in Eqs.(3.14a) and (3.14b) and Eqs.(3.7a) and (3.7b), respectively, it is possible to obtain an analytical expression for the phase sensitivity density:

FIG. 5 .
FIG.5.Calculated intensity distributions for a finite-width Gaussian pump with FWHM of 50 µm for the (a)-(c) noncompensated interferometer and for the (d) compensated interferometer for different phases.Same as for the plane-wave pump case, varying the phase, the intensity profiles of the compensated SU(1, 1) interferometer are only scaled as ∼ cos 2 (φ/2), see Eqs. (4.2a) and (4.2b), therefore, only one plot for all φ = π is presented.

3 FIG. 6 .
FIG. 6.(a) Calculated normalized phase sensitivity f for different gains for the non-compensated setup.For each gain, the crosses mark the points where f is minimized.The horizontal thin gray line indicates the shot noise level (f = 1).The divergences occur at points where d Ns,tot dφ = 0 [see Eq. (1.4)].(b) Calculated normalized phase sensitivity f for different gains for the compensated setup.Note that the entire plotted region is below the shot noise level.The inset shows the full phase interval [0, 2π].As mentioned in the text, the minimal value of f is always reached at φ = π.Furthermore, f diverges for φ → 2πk, k ∈ Z, see Eq. (4.7).In both cases, the pump is given by a Gaussian pump with an intensity-FWHM of 50 µm.

(4. 5 )
Plots of the covariance for different gains and phases are presented and discussed in Appendix F 2.

FIG. 8 .
FIG. 8.The Schmidt number for the plane-wave pump case [from Eqs.(5.7a) and (5.7b)] and finite-width pump case [from Eqs.(5.3a) and (5.3b)] for several parametric gains.Additionally, on the right-hand vertical axis, the plane-wave Schmidt number per length Lx (K (1) pw /Lx) is presented.The required values of Lx are shown in Fig. 7.
between the output plane-wave operators of the first crystal and the input plane-wave operators of the second crystal to connect Eqs.(B.1a) and (B.1b) with Eqs.(B.3a) and (B.3b), one can obtain the following connection relations between the transfer functions of the entire interferometer and each crystal:

FIG. 10 .
FIG.10.Values of the fitting parameter A versus the parametric gain.The error bars correspond to the standard error of the fitting procedure.The bars are scaled up by a factor of 200 along G and by a factor of 20 along A for better visibility.Note that the density of points is higher for G < 1.9 that results in decrease of the error below G = 1.9.See also TableIfor the selected values.
Appendix F: Covariances for a plane-wave and finite-width pump 1. Plane-wave pump Figs. 11(a)-11(f) present the covariance profiles for different gains in the plane-wave pump case.The profiles are described by the function C(q s ) defined in Eq. (3.7b).As one can see, the plane-wave covariance cov pw (q s , q ′ s ) [Eq. (3.7a)] vanishes unless q s = q ′ s , therefore, the function C contains all the information about the covariance.This is a direct consequence of the plane-wave pumping which results in perfect correlations between the planewave modes.As was already mentioned in the main text (Sec.III), the covariance profiles calculated for the noncompensated setup [Figs.11(a)-11(c)],and for the compensated setup [Figs.11(d)-11(f)], have shapes similar to the intensity profiles and broaden as the gain increases.Note that the figures are plotted in the external angles.

TABLE I .
Fitting parameter A for the selected gains G.