Extracting the field theory description of a quantum many-body system from experimental data
Torsten V. Zache Thomas Schweigler Sebastian Erne Jörg Schmiedmayer Jürgen Berges



Institut für Theoretische Physik, Ruprecht-Karls-Universität Heidelberg, Philosophenweg 16, 69120 Heidelberg, Germany

Vienna Center for Quantum Science and Technology, Atominstitut, TU Wien, Stadionallee 2, 1020 Vienna, Austria

Vienna Center for Quantum Science and Technology, Atominstitut, TU Wien, Stadionallee 2, 1020 Vienna, Austria
School of Mathematical Sciences, Centre for the Mathematics and Theoretical Physics of Quantum Non-Equilibrium Systems, University of Nottingham, University Park, Nottingham NG7 2RD, UK

Vienna Center for Quantum Science and Technology, Atominstitut, TU Wien, Stadionallee 2, 1020 Vienna, Austria

Institut für Theoretische Physik, Ruprecht-Karls-Universität Heidelberg, Philosophenweg 16, 69120 Heidelberg, Germany





Quantum field theory is a powerful tool to describe the relevant physics governing complex quantum many-body systems. Here we develop a general pathway to extract the irreducible building blocks of quantum field theoretical descriptions and its parameters purely from experimental data. This is accomplished by extracting the one-particle irreducible (1PI) correlation functions from which one can construct all physical observables. To match the capabilities of experimental techniques, our approach employs a formulation of quantum field theory based on equal-time correlation functions only. We illustrate the theoretical foundations of our procedure by applying it to the sine-Gordon model in thermal equilibrium, and then demonstrate explicitly how to extract these quantities from an experiment where we quantum simulate the sine-Gordon model by two tunnel-coupled superfluids. We extract all 1PI correlation functions up to the 1PI four-point function (interaction vertex) and their variation with momentum, encoding the `running' of the couplings. The measured 1PI correlation functions are compared to the theoretical estimates, verifying our procedure. 
Our work opens new ways of addressing complex many-body questions emerging in a large variety of settings from fundamental science to practical quantum technology.


Extracting the field theory description of a quantum many-body system from experimental data
    Jürgen Berges
    December 30, 2023
============================================================================================




Quantum Field Theory (QFT) has a wide range of very successful applications from early-universe cosmology and high-energy physics to condensed matter physics. A central aspect of QFT is that it describes the many-body limit of complex interacting quantum systems, which is also relevant for quantum technology if devices become large. 
Present large-scale analog quantum simulators using ultra-cold atoms explore the many-body limit described by QFT, e.g. <cit.>. 
Therefore, they may also be used to solve outstanding theoretical problems of QFT that are beyond classical computational techniques.
 
One of the big experimental challenges is probing the complex many-body states. One strategy is to detect every constituent (atom, superconducting qubit, quantum dot …) and its state. Such detections constitute a projective measurement of the many-body wave function in the constituent basis. For large systems such a measurement contains way too much information to be ever analysed fully. This is reflected by the exponential complexity of ‘tomography’ that prevents a complete characterization of the many-body quantum states <cit.>.
 
By contrast, there are important simplifications occurring in the many-body limit described by QFT. In QFT, often only a small subset of the microscopic details of the underlying theory is relevant for the computation of measurable physical properties. This effective loss of details has its mathematical foundation in the renormalization program of QFT <cit.>. As a result, for a quantum simulation of such a theory, many of the detailed properties of the microscopic quantum device have no effect on the simulation outcome for quantities of interest <cit.>.

This raises the important question of how to extract from experimental data the relevant information content of QFT. It is well known that an efficient description of QFT can be based on one-particle irreducible (1PI) correlation functions, called irreducible or proper vertices <cit.>. They represent the irreducible building blocks from which all physical observables may be constructed. This can be, e.g. the effective Hamiltonian determining the macroscopic dynamics, a possible spectrum of quasi-particles and their effective interaction strength. In a general setting, these vertices are functions of space and time or momentum and frequency, encoding the ‘running' of couplings prominently discussed in high-energy physics in the framework of the Standard Model of particle physics.
 
In principle, the irreducible vertices can be extracted from higher-order correlation functions <cit.>. The standard procedure employs correlation functions involving large time differences. While this is very suitable for high-energy collider experiments, where an analysis is based on the concept of asymptotic states in the infinite past and future, this is not adequate for many realizations of strongly interacting many-body systems where the notion of an initial state `long before' and a final state `long after' the collision is not physical. Moreover, often these systems are studied at a given snapshot in time, without any direct reference to states in the asymptotic past or future. This is especially true for cold-atom experiments where one takes pictures, for example measuring every atom either after time of flight <cit.> or in-situ <cit.>.
 
In this paper, we develop a pathway to extract the irreducible vertices of a quantum many-body system from experimental measurements. Our approach employs a formulation of QFT based on equal-time correlation functions only <cit.>. Equal-time correlation functions can be extracted from snapshot measurements <cit.> and, therefore, match well with experimental capabilities. We lay out the theoretical foundations of this approach, and illustrate the derivations using the sine-Gordon model. The irreducible vertices at equal times are estimated for this model both analytically and using numerical simulations. In particular, we show how to recover from the vertices the effective Hamiltonian underlying the dynamics. These theoretical results provide the basis for the benchmark verification of the QFT description extracted from experimental measurement. In the experiment, the sine-Gordon model is quantum simulated with two tunnel-coupled superfluids in thermal equilibrium <cit.>. We show how to extract the irreducible vertices from the experimental setup and compare the measurements to the theoretical estimates. The agreement of the experimental results with the theoretical expectations within errors provide a proof-of-principle verification of the approach. This represents an important step towards quantum simulator applications that are beyond reach of classical computational techniques. A first example of such an application is the recent experimental extraction of the irreducible two- and four-vertices for a strongly correlated spin-1 Bose condensate far from equilibrium <cit.>, where no theoretical solution is available and which has been performed in parallel to this work.

The paper is organised as follows. 
We start in section <ref> with a self-contained description of an equal-time formulation of quantum field theory and  equal-time correlation functions as they arise naturally in experiments. In particular, we show how the one-particle irreducible (1PI) vertices, which constitute the fundamental building blocks of the QFT description of the many-body system, can be extracted from the measured equal-time correlation functions. 
In section <ref>, we illustrate these theoretical foundations in the framework of the quantum sine-Gordon (SG) model <cit.>  and calculate the 1PI correlation functions and the effective action in the classical field theory limit in thermal equilibrium and compare it to numerical simulations. 
As a proof-of-principle, we show in section <ref> an application to an experiment with two tunnel-coupled superfluids, which realises the SG model  <cit.>. 
We conclude our work in section <ref>.
An extensive appendix contains detailed calculations.




§ EXTRACTING THE IRREDUCIBLE VERTICES  FROM EQUAL-TIME CORRELATIONS 


In the standard formulation of quantum field theory one starts from a typical scattering experiment which gives access to the transition amplitude between an initial state at times long before the collision and its final state at much later times. These transition amplitudes determine the S-matrix elements, which can be expressed in terms of time-ordered correlation functions of the underlying quantum field theory <cit.>. Knowledge of all time-ordered correlation functions is then equivalent to solving the quantum theory <cit.>.

However, time-ordered correlation functions and the description by an S-matrix formulation are  conceptually less suitable in the analysis of strongly correlated complex quantum systems, which are often studied at a given snapshot in time. Such measurements at a given instant of time lead to the notion of equal-time correlation functions. In quantum field theory, these can be represented by expectation values of Weyl ordered products of field operators <cit.>. Knowledge of all equal-time correlation functions at a given time t contains all information about the many-body system at this instant of time. For example, the factorisation properties of higher-order correlation functions directly reveal if the system is free (factorising) or interacting (non factorising) <cit.>.  To extract the interaction constants of the underlying (effective) Hamiltonian one has to extract the so-called one-particle irreducible (1PI) correlation  functions <cit.> which represent the full non-perturbative interaction vertices of the quantum system.

While there  are standard textbook concepts to extract the 1PI correlation functions from time-ordered correlation functions, the possibility to extract them from equal-time correlation functions is much less explored. Here we illustrate how to extract them from the equal-time correlations and thereby show how to determine the effective Hamiltonian from experiment at a snapshot in time.

We start in section <ref> with an introduction to quantum field theory in an equal-time formalism. At the example of a scalar field theory, we show the relation to Wigner's phase-space formalism commonly used e.g. in quantum optics. We further summarise how to extract connected correlation functions (<ref>) and one-particle irreducible vertices (<ref>) by introducing suitable generating functionals. Finally, we approximately calculate the 1PI effective action in thermal equilibrium in section <ref>, which provides a direct connection to the parameters of the microscopic Hamiltonian, and give a recipe on how to proceed (<ref>).




 §.§ Equal-time formulation of quantum field theory


The use of equal-time correlations is motivated by the progress of cold atomic setups which nowadays allow to extract highly resolved images at a given instant in time. It has long been known that QFT can be set up by only employing such equal-time information, without relying on multi-time correlations  <cit.>. This formulation has, however, never been widely used.  Theoretical progress in solving the equal-time formalism is  hampered by the lack of appropriate approximation schemes. Nevertheless, an equal-time formulation is perfectly suited to extract the irreducible vertices from experimental data representing a snapshot of the system at a fixed time.

Setting up an equal-time formulation relies on measurements of conjugate elementary operators that are non-commuting. As a consequence, one has to choose an ordering prescription <cit.>. Since correlations in cold atom systems are straightforwardly obtained by multiplying and averaging single shot results, the obtained correlations correspond to a fully symmetrised (so-called Weyl) ordering of the quantum operators. Moreover, as we show below, this choice of ordering leads to a definition of 1PI correlators that is directly related to Hamiltonian parameters. For the rest of this paper, we thus focus on Weyl-ordered correlation functions. A short discussion of other ordering prescriptions is given in the appendix <ref>.

For simplicity, we start with a real scalar field theory with Schrödinger field operators Φ̂(𝐱) and Π̂(𝐱) that fulfill the canonical commuation relation

    [Φ̂(𝐱), Π̂(𝐲)] = i ħδ(𝐱- 𝐲)   .

A general quantum state at time t is in the Schrödinger picture described by the density operator ρ̂_t.
Equivalently, knowing all correlations characterises the state ρ̂_t (see appendix <ref> for more details). Formally, all correlations can be conveniently summarised in the generating functional

    Z_t[J] =  Tr[ ρ̂_t exp(J^φ_𝐱Φ̂_𝐱 + J^π_𝐱Π̂_𝐱)]   .

Here we have introduced a notation where repeated indices are integrated over, e.g. J^φ_𝐱Φ̂_𝐱 = ∫ d^dx   J^φ(𝐱) Φ̂(𝐱). The J's are so-called source fields, i.e. auxiliary variables that encode the dependence of ρ̂_t on Φ̂ and Π̂, as indicated by φ and  π.
In the definition of Z_t, we have implemented the choice of ordering by treating the conjugate fields Φ̂ and Π̂ symmetrically. The resulting correlation functions, which are obtained taking functional derivatives are Weyl-ordered. For example at second order, we have

    G^(2)_𝐱,𝐲(t) = [ ⟨φ_𝐱φ_𝐲⟩_W_t ⟨φ_𝐱π_𝐲⟩_W_t; ⟨π_𝐱φ_𝐲⟩_W_t ⟨π_𝐱π_𝐲⟩_W_t ] ,

which consists of the three independent correlators

	
    ⟨φ_𝐱φ_𝐲⟩_W_t   = .δ^2 Z_t[J]/δ J^φ_𝐱δ J^φ_𝐲|_J=0 = 1/2Tr[ ρ̂_t (Φ̂_𝐱Φ̂_𝐲 + Φ̂_𝐲Φ̂_𝐱)  ]  ,
    ⟨φ_𝐱π_𝐲⟩_W_t   = .δ^2 Z_t[J]/δ J^φ_𝐱δ J^π_𝐲|_J=0 = 1/2Tr[ ρ̂_t (Φ̂_𝐱Π̂_𝐲 + Π̂_𝐲Φ̂_𝐱)  ]  ,
    ⟨π_𝐱π_𝐲⟩_W_t   = .δ^2 Z_t[J]/δ J^π_𝐱δ J^π_𝐲|_J=0 = 1/2Tr[ ρ̂_t (Π̂_𝐱Π̂_𝐲 + Π̂_𝐲Π̂_𝐱)  ]  .


This is explicitly verified in appendix  <ref>, where also the higher-order case is discussed. 

For a general quantum many-body system, Z_t[J] may involve more than one pair of canonically conjugated fields. These can be incorporated by adding appropriate sources J and essentially does not affect the general discussion. 


In the following, we consider only correlators of φ to lighten the notation. Nevertheless, φ may stand for either of the two fields and π is only written explicitly when necessary to avoid confusion. 
In general, we then denote all Weyl-ordered correlators as

    G^(n)_𝐱_1, …,𝐱_n(t) = ⟨φ_𝐱_1⋯φ_𝐱_n⟩_W_t = .δ^n Z_t[J]/δ J_𝐱_1⋯δ J_𝐱_n|_J=0 .

We refer to the appendix <ref> for a summary of all notational conventions used throughout this paper. In Eq. (<ref>), we have assumed a proper normalisation, ρ̂_t = 1, which implies Z_t[0] = 1.

To make use of established, powerful QFT tools, we seek a representation of Z_t in terms of functional integrals.
As shown in the appendix <ref>, the expression (<ref>) can be rewritten as a

    Z_t[J] = ∫𝒟φ 𝒟π  
    W_t[φ,π] exp(J^φ_𝐱φ_𝐱+J^π_𝐱π_𝐱)  .

The integration kernel W_t can be interpreted as a quasi-probability distribution,

    W_t[φ, π] = ∫𝒟φ̃⟨φ - φ̃/2| ρ̂_t| φ + φ̃/2⟩exp(i/ħφ̃_𝐱π_𝐱)  ,

the so-called Wigner functional. Here |φ⟩ denotes an eigenstate of the field operator Φ̂ with eigenvalue φ. The Wigner function itself has previously been applied successfully in the context of quantum optics <cit.> and also plays a prominent role in the semi-classical description of non-equilibrium quantum dynamics <cit.>. Eq. (<ref>) is the basis of the equal-time formulation of QFT and allows us to apply established procedures in the following.



 §.§ Connected correlation functions

In QFT (and analogously in classical probability theory) it is well known that the correlations encoded in Z_t are largely redundant <cit.>. The first step is the removal of additive redundancies, by the introduction of another generating functional,

    E_t[J] = log Z_t[J]  .


We denote the corresponding correlations, called connected correlators, as

    G^(n)_c,𝐱_1, …, 𝐱_n(t) = .δ^n E_t[J]/δ J_𝐱_1⋯δ J_𝐱_n|_J=0 .

Explicitly, as shown in the appendix <ref>, up to fourth order they are given by




	
    G^(1)_c,𝐱_1   = G^(1)_𝐱_1 , 
    
    	G^(2)_c,𝐱_1,𝐱_2   = G^(2)_𝐱_1,𝐱_2 - G^(1)_c,𝐱_1 G^(1)_c,𝐱_2 ,
    
    	G^(3)_c,𝐱_1,𝐱_2,𝐱_3   = G^(3)_𝐱_1,𝐱_2,𝐱_3 - (G^(2)_c,𝐱_1,𝐱_2 G^(1)_c,𝐱_3  + G^(2)_c,𝐱_2,𝐱_3 G^(1)_c,𝐱_1 + G^(2)_c,𝐱_3,𝐱_1 G^(1)_c,𝐱_2) - G^(1)_c,𝐱_1 G^(1)_c,𝐱_2 G^(1)_c,𝐱_3 ,
    
    	G^(4)_c,𝐱_1,𝐱_2,𝐱_3,𝐱_4   = 
    	G^(4)_𝐱_1,𝐱_2,𝐱_3,𝐱_4 - ( G^(3)_c,𝐱_1,𝐱_2,𝐱_3G^(1)_c,𝐱_4 +  G^(3)_c,𝐱_2,𝐱_3,𝐱_4G^(1)_c,𝐱_1 +  G^(3)_c,𝐱_3,𝐱_4,𝐱_1G^(1)_c,𝐱_2 +  G^(3)_c,𝐱_4,𝐱_1,𝐱_2G^(1)_c,𝐱_3)
       - (G^(2)_c,𝐱_1,𝐱_2G^(2)_c,𝐱_3,𝐱_4 + G^(2)_c,𝐱_1,𝐱_3 G^(2)_c,𝐱_2,𝐱_4 + G^(2)_c,𝐱_1,𝐱_4G^(2)_c,𝐱_2,𝐱_3 + ) - 
    	( G^(2)_c,𝐱_1,𝐱_2 G^(1)_c,𝐱_3G^(1)_c,𝐱_4 + G^(2)_c,𝐱_1,𝐱_3 G^(1)_c,𝐱_2G^(1)_c,𝐱_4. 
       + . G^(2)_c,𝐱_1,𝐱_4 G^(1)_c,𝐱_2G^(1)_c,𝐱_3 +G^(2)_c,𝐱_2,𝐱_3 G^(1)_c,𝐱_1G^(1)_c,𝐱_4 +G^(2)_c,𝐱_2,𝐱_4 G^(1)_c,𝐱_1G^(1)_c,𝐱_3 +G^(2)_c,𝐱_3,𝐱_4 G^(1)_c,𝐱_1G^(1)_c,𝐱_2) 
       - G^(1)_c,𝐱_1G^(1)_c,𝐱_2G^(1)_c,𝐱_3G^(1)_c,𝐱_4 ,





where we suppressed the overall time-dependence for brevity.
For every order n, the connected part G^(n)_c, is obtained by subtracting the information already given by lower-order functions G^(m<n)_c.

This can be visualised by a graphical representation in terms of Feynman diagrams, which is very helpful to organise the underlying combinatorics. By careful examination of this reorganisation, exemplified in Fig. <ref>, one learns that only connected graphs contribute to the correlators generated by E_t, hence the name connected correlations. In short, taking the logarithm of Z_t in Eq. (<ref>) removes all disconnected diagrams.



Physically, this means that by inspecting the factorisation of higher-order correlation functions, one can determine whether or not the quantum system is described by a gaussian density operator ρ̂_t. Since gaussian distributions correspond to free (non-interacting) QFTs, this in principle allows to determine the basis of conjugate fields which diagonalises the quantum many-body Hamiltonian that governs the system at hand.



 §.§ One-particle irreducible vertices

The connected correlators of order higher than two still contain redundant information. In order to access the irreducible vertices, we define the effective action

    Γ_t[Φ] = -E_t[J(Φ)] + J_𝐱(Φ) Φ_𝐱

as a Legendre transform of E_t. In Eq. (<ref>) the relation Φ_𝐱(J) = (δ E_t[J])/(δ J_𝐱) thus has to be inverted to obtain J_𝐱(Φ). We emphasise that the above notation is an abbreviation for a double Legendre transform in both of the conjugate fields Φ and Π. Accordingly, equations in this section implicitly include appropriate sums over the two fields.

Expanding the effective action in a functional Taylor series we have

    Γ_t[Φ] = ∑_n=2^∞1/n!Γ^(n)_𝐱_1, …,𝐱_n(t) ∏_j=1^n (Φ_𝐱_j -  Φ̅_𝐱_j(t))  .

Here Φ̅_𝐱(t) = ⟨φ_𝐱⟩_W_t is the mean value at time t, for which the effective action is stationary, i.e. (δΓ_t[Φ])/(δΦ)|_Φ = Φ̅ = 0.
The 1PI vertices,

    Γ^(n)_𝐱_1, …, 𝐱_n(t) = .δ^nΓ_t[Φ]/δΦ_𝐱_1⋯δΦ_𝐱_n|_Φ=Φ̅ ,

are the expansion coefficients in this series.
In Eq. (<ref>), the sum starts at n=2 because we have omitted an irrelevant constant Γ^(0) and the first order, Γ^(1), vanishes by construction due to the expansion around Φ̅. Physically, Φ̅ can take a non-vanishing value, which plays a crucial role, e.g., in the case of spontaneous symmetry breaking or the false vacuum decay <cit.>.

As shown in the appendix <ref>, the 1PI vertices up to fourth order are related to the connected correlation functions as follows:



	
		
    Γ^(1)_𝐱_1   = 0 , 
    Γ^(2)_𝐱_1,𝐱_2   =  [G_c^(2)]^-1_𝐱_1,𝐱_2 , 
    Γ^(3)_𝐱_1,𝐱_2,𝐱_3   = -Γ^(2)_𝐱_1,𝐲_1Γ^(2)_𝐱_2,𝐲_2Γ^(2)_𝐱_3,𝐲_3G^(3)_c,𝐲_1,𝐲_2,𝐲_3 , 
    Γ^(4)_𝐱_1,𝐱_2,𝐱_3,𝐱_4   = - Γ_𝐱_1,𝐲_1Γ_𝐱_2,𝐲_2Γ_𝐱_3,𝐲_3Γ_𝐱_4,𝐲_4 G^(4)_c,𝐲_1,𝐲_2,𝐲_3,𝐲_4
       + Γ^(2)_𝐱_1,𝐲_1Γ^(2)_𝐱_2,𝐲_2Γ^(2)_𝐱_3,𝐲_3Γ^(2)_𝐱_4,𝐲_4Γ^(2)_𝐳_1,𝐳_2(G^(3)_c,𝐲_1,𝐲_2,𝐳_1 G^(3)_c,𝐳_2,𝐲_3,𝐲_4+ G^(3)_c,𝐲_1,𝐲_3,𝐳_1 G^(3)_c,𝐳_2,𝐲_2,𝐲_4  + G^(3)_c,𝐲_1,𝐲_4,𝐳_1 G^(3)_c,𝐳_2,𝐲_2,𝐲_3)  .

	

We again emphasise that the explicit equations should be understood including appropriate sums over φ and π correlators (see appendix <ref>).
For higher orders these relations become more complicated and calculations are conveniently performed with the graphical notation exemplified in Fig. <ref>. These diagrams also explain the attribute 1PI: The diagrams representing the vertices can not be disconnected by cutting a single line. In this sense they are the irreducible structures from which all correlation functions and thus all physical observables can be recovered.




This also justifies the name effective action: Γ_t is the quantum generalisation of a classical action including all corrections due to quantum-statistical fluctuations. However, in contrast to the `standard' (unequal-time) action, there is one stationarity condition for each time t. Together, they do not give a time evolution equation for the one-point function in the usual sense, but one differential equation for each time t. In that sense the time t is treated as a label in the equal-time formulation of QFT.



 §.§ Measuring the effective Hamiltonian


So far, we have equivalently rewritten the quantum-statistical information of a system described by a density operator ρ̂_t in terms of generating functionals Z_t, E_t and Γ_t which encode full, connected and 1PI correlation functions, respectively. While the entries of the density matrix are typically inaccessible and less intuitive, the equal-time correlators can be measured in experiments and are directly related to relevant observables and structural information, such as occupation numbers and couplings.  Next, we employ the equal-time formalism to relate parameters of an Hamiltonian to the 1PI correlators.

As a generic example we consider a relativistic scalar field theory with potential V described by the Hamiltonian

    Ĥ =∫_𝐱[1/2 Π̂_𝐱^2 + 1/2(∇_𝐱Φ̂_𝐱)^2  + V (Φ̂_𝐱)]   .
	
Given the Hamiltonian Ĥ, it is possible to derive an evolution equation for Γ_t <cit.>.
For simplicity, we focus on the case of thermal equilibrium, which is a stationary solution Γ_β, described by the density operator ρ̂_β∼exp(- βĤ) with the prefactor fixed by normalisation.  

In order to obtain the generating functional Eq. (<ref>) we need to calculate the Wigner functional Eq. (<ref>). In the interacting case, the involved functional integration can only be performed approximately. It is however possible to derive an exact equation for the thermal Wigner functional (see appendix <ref>). It takes the form of a functional flow equation, ∂_β W_β = - (H_0 + ħ^2 H_1 + … )W_β, with
	
    H_0    = ∫_𝐱[1/2π_𝐱^2 + 1/2(∇_𝐱φ_𝐱)^2 +V_𝐱(φ)]   .

It is straightforward to solve the equation for W_β perturbatively by a semi-classical expansion in powers of ħ. 

The leading order is the classical field theory limit, where we obtain W_β∼exp(-β H_0) with the classical Hamiltonian H_0. Then the generating functional Eq. (<ref>) becomes

    Z_β[J] ∼∫𝒟π𝒟φ  e^-β H_0 + J^φ_𝐱φ_𝐱 + J^π_𝐱π_𝐱  .

Thus β H_0 plays the role of a classical action for the fluctuating fields φ and π. This allows us to calculate the effective action Γ_β in the equal-time formalism using established QFT methods, such as a the background field method employed below. 

We note that the two conjugate fields φ and π decouple in the present limit, which implies that the effective action separates as

    Γ_β[Φ,Π] = Γ_β[Φ] + Γ_β [Π]

with Γ_β[Π] = β/2∫_𝐱Π_𝐱^2 + const., as shown in appendix <ref>.
The separation Eq. (<ref>) is a property of the classical field approximation. In general, the quantum effective action of the full quantum theory requires knowledge of all equal-time correlators of φ and π, including mixed terms. However, symmetries such as time translation invariance simplify the discussion, see appendix <ref>.

By means of the background field method, we can calculate the effective action in a loop expansion. In the present formalism, we split

    Γ_β[Φ] = β H[Φ] + Γ'_β[Φ]  ,

where H[Φ] = H_0[φ = Φ, π = 0]. As shown in the appendix <ref>, the `rest' Γ_β' obeys the following functional integro-differential equation,

    e^-Γ'_β[Φ] = ∫𝒟φ exp(-β K[φ,Φ] + δΓ'_β[Φ]/δΦ_𝐱φ_𝐱) ,

where we abbreviated
	
    K[φ,Φ]    = H[Φ + φ] - H[Φ] - ∫_𝐱δ H[Φ]/δΦ_𝐱 φ_𝐱  .

The solution of this equation is organised diagrammatically as an expansion in the number of loops. At leading order (tree-level) in this expansion Γ'_β = 0 and thus the equal-time effective action is directly related to the microscopic Hamiltonian. Consequently, the 1PI vertices correspond to the interaction constants of the underlying system. Beyond the leading-order approximation, the notion of the microscopic Hamiltonian becomes a less useful concept. The effective action then plays the role of an effective Hamiltonian, with all corrections from quantum-statistical fluctuations taken into account.

Returning to the leading order approximation, which gives rise to the tree-level 1PI vertices, we explicitly have


    Γ^(2)_𝐱,𝐲   =  ∇_𝐱^2δ(𝐱-𝐲)   +∫_𝐳. δ^2 V_𝐳(Φ)/δΦ_𝐱δΦ_𝐲|_Φ=Φ̅ ,
    Γ^(n)_𝐱_1, …, 𝐱_𝐧   = ∫_𝐳. δ^n V_𝐳(Φ)/δΦ_𝐱_1⋯δΦ_𝐱_n|_Φ=Φ̅ ,


where n≥3.
Eq. (<ref>) or more explicitly Eq. (<ref>) directly show the relation between the 1PI correlation functions and the parameters of the microscopic (or more generally an effective) Hamiltonian. 
Together with the procedure to obtain the 1PI correlators, outlined below, they provide an experimental prescription for measuring a quantum many-body Hamiltonian. 



 §.§ Recipe to extract 1PI correlators


The extraction of 1PI correlators, which are the fundamental irreducible building blocks for the QFT description, from equal-time data proceeds in the following steps.  


	
  * Identify the degrees of freedom of interest which constitute the elementary fields φ of the QFT. 
	
  * Obtain many realisations (i = 1, …, N) of the desired field φ_i(𝐱) at the times t of interest.
	
  * Estimate the full correlators up to order n by averaging G^(n)_𝐱_1, …, 𝐱_n ≈1/N∑_i φ_i(𝐱_1) ⋯φ_i(𝐱_n).
	
  * Obtain the connected correlators G_c^(n) by subtracting the disconnected contributions according to Eq. (<ref>).
	
  * Calculate the 1PI correlators Γ^(n) by reducing the connected correlators according to Eq. (<ref>).


This procedure corresponds to a shift of representation from the density operator ρ̂ to Γ, the generating functional for 1PI correlators.
In the following two sections we will illustrate and verify the method in the case of the sine-Gordon model with numerically simulated data (Section <ref>) and with experimental measurements (Section <ref>).



§ EXAMPLE: SINE-GORDON MODEL
 

As an explicit example, we consider the sine-Gordon model  <cit.>  in thermal equilibrium.  It is an interacting relativistic scalar field theory described by

    βĤ_SG = ∫_x {β gΠ̂_x^2 +  λ_T/4[  1/2(∂_xΦ̂_x)^2 - 1/ℓ_J^2cos( Φ̂_x)] }  ,

where β = (k_B T)^-1 is the inverse temperature. The specific form of the Hamiltonian Ĥ_SG given above is motivated by the recent progress to quantum simulate the SG model by two tunnel coupled 1D superfluids <cit.>. See section <ref> for the physical origin of the fields Φ̂ and Π̂, the microscopic parameter g and the length scales λ_T and ℓ_J. 

The semiclassical approximation Eq. (<ref>) is valid for 

    √(4 γ)≪min[1,4/Q]   ,

where the dimensionless parameters are γ = 16 g β/λ_T and Q = λ_T / ℓ_J. In the semi-classical limit the loop expansion is controlled by Q with the tree-level approximation valid for 1/Q ≪ 1 (see appendix <ref> and <ref> for details). We therefore consider in the following λ_T = 17.35  μm and vary ℓ_J such that 1 ≲ Q ≲ 20.

Following the general discussion of the previous section, the tree-level vertices corresponding to Eq. (<ref>) are 

	
    Γ^(2),tree_p   =  λ_T/4(p^2 + 1/ℓ_J^2)  ,
    Γ^(2n),tree_p_1, … p_2n-1   = - λ_T/4ℓ_J^2  (-1)^n  , 
    Γ^(2n-1),tree_p_1, … p_2n-2   = 0 ,


where n>2 and we switched to momentum space correlators. Here and in the following we always consider the diagonal part in momentum space, i.e. , which removes the volume factors arising from translation invariance. Note that Γ^(2n-1) = 0   (∀ n ≥ 1) remains valid beyond the tree-level approximation due to the symmetries of the SG Hamiltonian.



Employing a stochastic process based on a transfer matrix formalism <cit.>, we numerically obtain thermal profiles of the field φ_x, corresponding to the operator Φ̂_x. These are exact solutions of the SG model within the semi-classical approximation and hence have contributions up to arbitrary order in the above loop expansion. With these statistical samples, we carry out the procedure described in section <ref> and calculate the 1PI correlators up to fourth order.

So far, we have implicitly assumed in Eq. (<ref>) that the correlations are obtained for an infinite system with periodic boundary conditions. The employed numerics, however, yield correlators from a finite subsystem, which is better described by open boundary conditions. We therefore employ a cosine transform and translate the results to momentum space, i.e. Fourier momenta (for details see appendix <ref>). 

Figure <ref> shows the calculated 1PI vertices in momentum space for a large value of Q ≈ 11.5 and different volumes L. Note that due to the periodicity of the sine-Gordon potential[The SG model is invariant under the shift φ→φ + ℤ× 2 π. This leads to an undefined offset for the numerical profiles φ_x and hence an undefined value of the momentum correlators for p=0.] the value of p=0 is not defined for the correlations considered. Therefore the fact that the correlation function is diagonal is crucial to be able to perform the inversion of the connected second order correlation function in order to obtain Γ^(2). The diagonal form allows to do the inversion for p ≠ 0 without knowing the value for p=0.

We find excellent agreement with the tree-level predictions for the momentum diagonal of the two and four vertex. This demonstrates the possibility to carry out the procedure described in the previous section, which allows to directly measure the microscopic parameters through equal-time 1PI correlation functions if higher-order loop corrections can be neglected.



However, when the system becomes strongly correlated the microscopic details become irrelevant and replaced by effective, momentum-dependent (so-called running) couplings. This behaviour is precisely captured by the 1PI correlation functions, which effectively replace the microscopic coupling parameters appearing in the Hamiltonian H_SG. We therefore adjust the parameters away from the weakly-coupled limit Q ≫ 1. In general, we observe an increasing deviation from the tree-level approximation, which is expected because loop corrections due to increasing fluctuations modify the physics.

Quantitatively, the corrections to the 1PI two-point function are summarized in the self-energy Σ,  defined via

    Γ^(2)_p = Γ^(2),tree_p+ Σ_p  .
 
In the appendix <ref>, we calculate the leading correction, 

    Σ_p^one-loop = -1/(4ℓ_J)   .


In Fig. <ref>, the 1PI two-point function and the self-energy are plotted as a function of p. Generically, tree-level dominates in the ultraviolet (i.e. at high momenta), which we also observe numerically. The 1PI two-point function approaches the power-law ∝ p^2 in this limit and the (normalised) self-energy vanishes.

In the infrared (low momenta), however,  loop corrections are important. It is this regime where collective macroscopic phenomena emerge and the microscopic details are washed out. We observe a negative self-energy and hence a reduction of the 1PI two-point function, which agrees with the one-loop result over an intermediate range of momenta (and Q). Physically, this result implies stronger fluctuations as Q decreases, consistent with the expectation for a strongly correlated regime of the sine-Gordon model.



Similarly, the loop corrections to the 1PI vertices lead to the notion of running couplings, i.e. momentum-dependent interaction vertices that deviate from the constant microscopic values. In the appendix <ref>. we calculate the one-loop vertex

    Γ^(4),one-loop_p = - λ_T/4 ℓ_J^2 - 1/8 ℓ_J^31/p^2 + 1/ℓ_J^2  .
 
Again, it is expected that loop corrections vanish for high momenta and the 1PI four-vertex converges to the tree-level result, i.e. the microscopic parameters of the Hamiltonian, which is confirmed by our numerical simulations.

This is demonstrated in Fig. <ref>, where the 1PI four-vertex is shown for the same values of Q as in Fig. <ref>. For very high momenta, we are again limited by finite statistics. In the infrared, we clearly observe the momentum-dependent, i.e. running, coupling. The increased values indicate stronger interactions, qualitatively consistent with the one-loop calculation. The effect is again more pronounced for smaller values of Q, as expected in the strongly correlated regime of the sine-Gordon model. 

We observe a qualitative difference between large and small values of Q. For Q ≳ 4, the magnitude of the 1PI four-vertex is increased in the infrared and shows a running coupling towards the smaller tree-level value. For Q ≲ 3, the magnitude of the vertex decreases in the infrared as compared to the tree-level value at higher momenta. This behaviour is also clearly visible in Fig. <ref>, where we show the four-vertex as a function of Q for fixed momenta.







§ EXPERIMENTAL RESULTS: PROOF-OF-PRINCIPLE





As a proof of principle to extract the 1PI vertices from experimentally measured correlations we apply the formalisms discussed above to the physical system of two tunnel-coupled one-dimensional superfluids in a DW potential on an atomchip.
Such a system can be seen as a quantum simulator of the sine-Gordon model <cit.>. The relative phase φ_x between the superfluids corresponds to Φ̂ in Eq. (<ref>) in section <ref> while the relative density fluctuations correspond to the conjugate field Π̂.

A schematic of the experimental system is given in Fig. <ref>.
The parameters in H_SG (<ref>) are related to the experimental parameters via , , and . Here the 1D effective interaction strength  is calculated from the s-wave scattering length a_s and the frequency ω_⊥ of the radial confinement; n_1D is the 1D density and m is the mass of the 87Rb atoms which the superfluids consist of. The single particle tunneling rate between the wells is denoted by J.

In the experiment, the two superfluids are prepared by slow evaporative cooling in the DW potential (the same way the slow cooled data presented in <cit.> was prepared). However, in contrast to <cit.>, the data used here was taken for a box-like longitudinal confinement <cit.> of 75 length.
Matter-wave interferometry <cit.> gives access to the spatially resolved relative phase fluctuations φ_x between the two superfluids.
More details about the experimental procedure and the data analysis can be found in <cit.>.



Starting from the measured phase profiles, we can calculate the 1PI vertices in the same way as was done for the numerics (see section <ref> and appendix <ref>). 
For box like potentials one naturally gets Neumann boundary conditions (BC) for the phase from the condition of vanishing particle current on the edges <cit.>.
From the cosine transform (compatible with the Neumann BC) of the complete system we therefore simply get the 1PI vertices of the Hamiltonian with this BC.
Acknowledging that our system is still too short to get results free from finite size effects, we nevertheless apply the conversion factors to Fourier momentum space given in (<ref>) and (<ref>) for consistency when presenting Γ^(n).

Let us start the discussion of the experimental results with the cosine transformed second-order correlation function

    G̃^(2)_c,p,p^' = 2/L( ⟨φ̃_p φ̃_p'⟩ - ⟨φ̃_p ⟩⟨φ̃_p'⟩).
 
Here φ̃_p represents the cosine transform (<ref>) over the finite interval with length L and we chose the prefactors for later convenience. The factor 2 comes from the identity (<ref>) and the factor 1/L from the delta function.  
We see from Fig. <ref> that the correlations are approximately diagonal. 
Further, note that density-phase two-point correlations, ⟨πφ⟩_W_t, vanish due to time-translation invariance of the thermal state, even beyond the semi-classical approximation Eq. (<ref>). Together, this enables us to calculate the 1PI two-point correlator as

    Γ^(2)_p = 1/G̃^(2)_c,p,p  ,

where we neglected the small off-diagonal elements of G̃^(2)_c,p,p^'.
The results are presented in Fig. <ref>.







All experimental results presented in this paper are corrected for the expected influence of the finite imaging resolution.
In our simple model, the imaging process leads to a convolution of the true phase profiles with a Gaussian function with σ_psf = 3 <cit.>.
In momentum space this simply leads to a multiplication with exp(-p^2σ_psf^2/2), which can be easily corrected by dividing the cosine transformed relative phase φ̃(p) by this factor.

In order to connect the experimental results to the theoretical model (section <ref>) we estimate λ_T = 11 for all the different measurements. The values for Q = λ_T/ℓ_J are then self consistently fitted from  <cit.>. We see good agreement between experiment and thermal sine-Gordon theory for the 1PI two-point function in Fig. <ref>.



Having obtained the two-point function, and using that the third-order correlation functions vanish for symmetry reasons, we can calculated the diagonal part of the 4-vertex as

    Γ^(4)_p = - 8/31/L⟨φ̃_p^4 ⟩_c ×(Γ^(2)_p)^4.
 
Here ⟨φ̃_p^4 ⟩_c stands for the diagonal elements of the cosine transformed fourth-order connected correlation function.
The factor 8/3 comes from the identity (<ref>), the factor 1/L again comes from the delta function. 
The results for the three lowest lying momentum modes are presented in Fig. <ref> as a function of Q.
We find qualitative agreement between experiment and theory as well as the expected approach towards the tree-level result for higher momenta.
The momentum dependence for the measurements with large experimental sample size (Fig. <ref>) reveal a running coupling with a qualitative agreement between the experiment and thermal sine-Gordon theory.



§ CONCLUSION


The presented method provides a general framework to extract and test the effective or emergent quantum field theoretical description of generic quantum many-body systems from experiments. For the example of the sine-Gordon model, which is quantum simulated with two tunnel-coupled superfluids, we have demonstrated how to experimentally obtain the irreducible vertices in thermal equilibrium and compared to theoretical expectations. This represents an essential step in the verification of the approach, which opens a new pathway to study fundamental questions of QFT through large-scale (analog) quantum simulators.
 
This becomes especially interesting for strongly correlated systems and in non-equilibrium situations, where it is often not possible to solve the theory using classical computational techniques. Extracting the irreducible building blocks of quantum many-body systems, and how they change with time, promises to provide detailed insights into the dynamics for these cases. The recent study of a spin-1 Bose condensate far from equilibrium <cit.>, which has been performed in parallel to this work and employed similar methods, presents an example where currently no theoretical solution is available. In turn, the insight from experimental measurement can support  theoretical developments in devising new approximation schemes and effective field theory descriptions.





§ ACKNOWLEDGMENTS
 
We thank Thomas Gasenzer, Philipp Kunkel, Igor Mazets, Markus Oberthaler, Maximilian Prüfer, Helmut Strobel and Christof Wetterich for discussions and Federica Cataldini, Sicong Ji, Bernhard Rauer and Mohammadamin Tajik for help with the experiment.
This work is part of and supported by the DFG Collaborative Research Centre “SFB 1225 (ISOQUANT)”, the ERC Advanced Grant “QuantumRelax”, and  the JTF project “The Nature of Quantum Networks” (ID 60478). TS acknowledges support by the Austrian Science Fund (FWF) in the framework of the Doctoral School Complex Quantum Systems (CoQuS). SE acknowledges support through the EPSRC Project Grant (EP/P00637X/1).
JS, JB and SE acknowledge the hospitality of the Erwin Schrödinger Institut in the framework of their thematic program Quantum Paths which enabled many discussions defining the foundational ideas shaping this article.








§ NOTATIONAL CONVENTIONS


In table <ref>, we have summarized the different generating functionals and the involved fields which appear throughout this paper.














We use the following notations: Operators are always indicated by a hat. Tr[…] indicates a trace over the full Hilbert space. The absence of a hat implies a c-number (i.e. commuting objects). In the whole formalism, the time t is treated as a label and often left out for brevity. Repeated spatial indices are integrated over, e.g. ; we write explicit integrals if there is room for confusion. 

In the appendix, we also use collective latin indices, a = (φ/π, 𝐱); then repeated indices are integrated or summed over as appropriate, e.g. J_a Φ̂_a = ∫ d^dx [J^φ(𝐱) Φ̂(𝐱) + J^π(𝐱) Π̂(𝐱)]. It is useful to think of correlation functions as tensors, e.g. G^(2)_c,ab and its inverse Γ^(2)_ab are the components of 2-tensors that fulfill G^(2)_c,abΓ^(2)_bd = δ_ad. Here δ_ad is the product of a discrete Kronecker delta and a continuous Dirac delta distribution.




§ TECHNICAL DETAILS OF SECTION <REF>

In this appendix we show explicit calculations and detailed dicussions that we left out in the discussion of equal-time correlation functions.




 §.§ Operator ordering at equal-time

There are three obvious choices of different orderings,
often referred to as symmetric (Weyl), normal (P) and anti-normal (Q). In terms of creation and annihiliation operators â_p and â^†_p, which fulfill [â_p, â^†_q] ∼δ_pq, they can be realized by the definitions

	
    Z^(W)_t[J]    =  Tr[ ρ̂_t exp(J_p â^†_p - J^*_p â_p )]   ,
    
    	Z^(P)_t[J]    =  Tr[ ρ̂_t exp(J_p â^†_p ) exp(  - J^*_p â_p )]   ,
    
    	Z^(Q)_t[J]    =  Tr[ ρ̂_t  exp(  - J^*_p â_p ) exp(J_p â^†_p ) ]   .


In general, there is a continuum of other choices that smoothly connect these three cases.
However, all different choices are fully equivalent in the sense that they contain all measurable information and the main difference lies in the associated quasi-probability distributions. For more details we refer to <cit.>.

Explicitly, for the choice of the main text, the ordering is resolved as


    exp(J^φ_𝐱Φ̂_𝐱 + J^π_𝐱Π̂_𝐱)    = e^J^φ_𝐱Φ̂_𝐱 e^J^π_𝐱Π̂_𝐱  e^-i/2 J^φ_𝐱J^π_𝐱
       = e^J^π_𝐱Π̂_𝐱e^J^φ_𝐱Φ̂_𝐱   e^i/2 J^φ_𝐱J^π_𝐱


where we used the BCH formula in the form

    e^Â+B̂ = e^Â e^B̂ e^-1/2[Â,B̂]  ,

which is valid for [[Â,B̂],B̂]=[[Â,B̂],Â]= 0. Thus, derivatives acting on exp(J^φ_𝐱Φ̂_𝐱 + J^π_𝐱Π̂_𝐱) from the left result in operators and additional sources according to

    δ/δ J^φ_𝐱→Φ̂_𝐱 - iJ^π_𝐱/2 ,       δ/δ J^π_𝐱→Π̂_𝐱 + iJ^φ_𝐱/2  .

Using this correspondence, it is straightforward to generate explicit expression for all correlators. For instance, at second order, we have

    δ/δ J^φ_𝐱δ/δ J^π_𝐲   →δ/δ J^φ_𝐱(  Π̂_𝐲 + iJ^φ_𝐲/2) 
       →(  Π̂_𝐲 + iJ^φ_𝐲/2)  (Φ̂_𝐱 - iJ^π_𝐱/2)+ i/2δ(𝐱- 𝐲)   .

Setting the sources to zero proves that

    ⟨φ_𝐱π_𝐲⟩_W_t   = Tr[ ρ̂_t Π̂_𝐲Φ̂_𝐱] + i/2δ(𝐱- 𝐲) 
       =  1/2Tr[ ρ̂_t (Φ̂_𝐱Π̂_𝐲 + Π̂_𝐲Φ̂_𝐱)  ]  ,

where we used the canonical commutation relations and the normalization of ρ̂_t.
In the following we drop the label t for brevity.



 §.§ Correlations and the density operator

The density operator ρ̂ can formally be recovered from Z as <cit.>

    ρ̂   = ∫𝒟 J^φ𝒟 J^π  Z[J][exp(J^φ_𝐱Φ̂_𝐱 + J^π_𝐱Π̂_𝐱) ]^-1  .

Furthermore, the mappings between the different functionals are invertible (under appropriate mathematical assumptions): W and Z are related by Fourier transforms, Z and E by an exponential (or logarithmic) map, E and Γ by Legendre transforms.
Thus, it is completely equivalent to work with the Wigner functional W or any of the generating functionals Z, E, Γ instead of the density operator ρ̂.



 §.§ Functional integral representation, Eq. (<ref>)


We seek a representation of Z in terms of classical (commuting) instead of operator-valued fields. To this end, we evaluate the trace as



	
		
    Z[J,ρ(t)] 
    		   = ∫𝒟φ^+ 𝒟φ^-  
    		φ^+ρ(t)φ^-φ^-e^J^φ_𝐱Φ_𝐱 + J^π_𝐱Π_𝐱φ^+
       = ∫𝒟φ^+ 𝒟φ^- 𝒟π̃ φ^+ρ(t)φ^-φ^-e^J^φ_𝐱Φ_𝐱π̃π̃e^J^π_𝐱Π_𝐱φ^+e^-i/2J^φ_𝐱J^π_𝐱
       = ∫𝒟φ 𝒟φ̃ 𝒟π̃ φ + φ̃/2ρ(t)φ - φ̃/2
    		e^J^φ_𝐱(φ_𝐱-φ̃_𝐱/2)+i (φ_𝐱-φ̃_𝐱/2)π̃_𝐱+J^π_𝐱π̃_𝐱-i (φ_𝐱+φ̃_𝐱/2)π̃_𝐱-i/2J^φ_𝐱J^π_𝐱
       = ∫𝒟φ 𝒟φ̃ 𝒟π̃ 𝒟π 
    		W_t[φ,π]
    		e^iπ_𝐱φ̃_𝐱+J^φ_𝐱(φ_𝐱-φ̃_𝐱/2)+J^π_𝐱π̃_𝐱-iφ̃_𝐱π̃_𝐱-i/2J^φ_𝐱J^π_𝐱
       = ∫𝒟φ 𝒟π  
    		W_t[φ,π] exp[J^φ_𝐱φ_𝐱+J^π_𝐱π_𝐱]   .

	



In the above calculation, we have again used the BCH formula, performed a change of variables φ^±≡φ±φ̃/2, and employed the definition of the Wigner functional.



 §.§ The connected correlators, Eq. (<ref>)

To discuss the explicit form of the connected and 1PI correlators, we use the short-hand notation, where sources J_a have a single index indicating space 𝐱, as well as any of the two fields φ and π. Repeated indices are summed and integrated over. Similarly, we abbreviate the fields as φ_a.

At first order, the connected correlators are directly related to the full Weyl-ordered one-point function,

    δ E/δ J_a = δlog Z/δ J_a = 1/Zδ Z/δ J_a  .

Setting the sources to zero, we have Z[J=0]=1, which proves Eq. (<ref>).

At second order, we calculate

    δ^2 E/δ J_aδ J_b = δ/δ J_a(1/Zδ Z/δ J_b) = 1/Zδ^2 Z/δ J_a δ J_b  -  1/Z^2δ Z/δ J_aδ Z/δ J_b  ,

which proves Eq. (<ref>).

The higher orders follow analogously by the combinatorics of the derivatives. E.g., the third order, Eq. (<ref>), is obtained by

    δ^3 E/δ J_aδ J_b δ J_c   = 1/Zδ^3 Z/δ J_a δ J_b δ J_c
       - (1/Z^2δ^2 Z/δ J_a δ J_bδ Z/δ J_c + 2 perm.) 
       +2  1/Z^3δ Z/δ J_aδ Z/δ J_bδ Z/δ J_c

and using Eqs. (<ref>) and (<ref>).




 §.§ The 1PI vertices, Eq. (<ref>)

The expression for the 1PI two-point function is central for the construction of the higher orders. It follows by considering a derivative of the stationarity condition,

    δ^2 Γ/δΦ_a δΦ_b = δ J_b/δΦ_a .

This is the matrix inverse of the derivative of the one-point function (in the presence of sources)

    δΦ_a(J)/δ J_b =δ/δ J_b(1/Zδ Z/δ J_a) =  δ^2 E/δ J_a δ J_b  .

Thus, we find Eq. (<ref>), or

    δ^2 Γ/δΦ_aδΦ_b = [(δ^2 E/δ J δ J)^-1]_ab  ,

which also holds without setting the sources to zero.

The higher orders follow by taking derivatives of this equation (with non-zero sources). To this end, we replace a derivative by

    δ/δΦ_c = 	δ J_c'/δΦ_cδ/δ J_c' = 	δ^2 Γ/δΦ_c δΦ_c'δ/δ J_c'

and calculate the derivative of the inverse of a matrix M(y) depending on a parameter y according to

    d/dy( M^-1) = - M^-1·dM/dy· M^-1  .

This results for the third order in

    δ^3 Γ/δΦ_aδΦ_b δΦ_c = - δ^2 Γ/δΦ_a δΦ_a'δ^2 Γ/δΦ_b δΦ_b'δ^2 Γ/δΦ_c δΦ_c'δ^3 E/δ J_a'δ J_b'δ J_c' ,

which proves Eq. (<ref>).

Similarly, the fourth order, Eq. (<ref>) , follows by the combinatorics of taking further derivatives,

	
    δ^4 Γ/δΦ_aδΦ_b δΦ_c δΦ_d =    - δ^2 Γ/δΦ_a δΦ_a'δ^2 Γ/δΦ_b δΦ_b'δ^2 Γ/δΦ_c δΦ_c'δ^2 Γ/δΦ_d δΦ_d'δ^4 E/δ J_a'δ J_b'δ J_c'δ J_d'
       + δ^2 Γ/δΦ_a δΦ_a'δ^2 Γ/δΦ_b δΦ_b'δ^2 Γ/δΦ_c δΦ_c'δ^2 Γ/δΦ_d δΦ_d'( δ^3 E/δ J_a'δ J_b'δ J_eδ^2 Γ/δΦ_e δΦ_fδ^3 E/δ J_fδ J_c'δ J_d'  + 2 perm.)   .






 §.§ The thermal case and the classical limit

In thermal equilibrium, the (unnormalized) canonical density operator ρ̂_β = e^-βĤ fulfills the equation . Employing the quasi-probability formalism <cit.>, one can show that this equation translates to an equation for W_β. It takes the form

    ∂_β W_β = -1/2[ H_W^+ + H_W^- ] W_β ,

where H^±_W = H_W[φ±iħ/2δ/δπ, π∓iħ/2δ/δφ] is a functional differential operator obtained from the Weyl-transform H_W[φ, π] of the Hamiltonian Ĥ by replacing the arguments with the given operators <cit.>. With the initial conditon W_β→∞ = const., which follows from the high-temperature limit, this functional flow equation can be solved perturbatively by exanding  in powers of ħ and comparing the coefficients. 

The first order in this expansion is the classical field theory limit, where W_β∼ e^-β H with the classical Hamiltonian H = H_W. Parametrically, this is a valid approximation when ħ is small.
Since ħ is a dimensionfull quantity, the precise power-counting of this expansion has to be determined for each theory separately. In general, any quantum system in thermal equilibrium will be governed by (at least) two dimensionless parameters ϵ_q and ϵ_th that control the strength of quantum and classical fluctuations, respectively. A sufficient condition for the validtiy of the classical approximation is

    ϵ_q≪min[1,ϵ_th]   .


The parameters ϵ_q and ϵ_th are obtained by rescaling the Hamiltonian and the fundamental fields to dimensionless quantities. Explicitly, for the sine-Gordon model in the form of Eq. (<ref>), we have

    βĤ   = 1/ϵ_th∫ dx' {1/2[(Π̂')_x'^2 + (∂_x'Φ̂'_x')^2 ] - cos(Φ̂'_x')} ,

together with the rescaled commutation relations

    [Φ̂'_x', Π̂'_y']    = i ϵ_qδ(x'-y')   ,

such that we find

    ϵ_q =  √(4 γ) ,       ϵ_th = 4 ℓ_J/λ_T = 4/Q  .

Here

    γ = mg/ħ^2 n_1D = γ_LL + 1/n_1D^2 ℓ_J^2  ,

is dominated by the 1D Lieb-Liniger parameter , such that the semi-classical approximation is valid in the weakly interacting regime γ_LL≪ 1, as expected.



 §.§ Derivation of Eq. (<ref>) and the Π-dependence
 
From Eq. (<ref>), the generating functional separates into a product Z_β[J] = Z^φ_β[J^φ] Z^π_β[J^π] with

	
    Z^φ_β[J^φ]    ∼∫𝒟φ  e^-β∫_𝐱[1/2(∇_𝐱φ_𝐱)^2 + V_𝐱(φ)] + ∫_𝐱 J^φ_𝐱φ_𝐱 ,
    
    	Z^π_β[J^π]    ∼∫𝒟π exp[-β/2∫_𝐱π_𝐱^2  + ∫_𝐱 J^π_𝐱π_𝐱]  
       ∼exp[1/2β∫_𝐱(J^π_𝐱)^2 ]    .


This directly implies that E[J] = log Z^φ_β[J^φ] + log Z^π_β[J^π] + const. and thus the effective action becomes Γ[Φ,Π] = Γ^φ[Φ] + Γ^π[Π] + const., which proves Eq. (<ref>). 

Carrying out the Legendre transform in J^π, we solve

    Π_𝐱(J^π) = .δ Z^π[J^π]/δ J_𝐱^π|_J=0 = J_𝐱^π/β⇒ J_𝐱^π(Π) = βΠ_𝐱 .

and finally obtain

    Γ^π[Π] = - log Z^π[J^π(Π)] + J_𝐱^π(Π) Π _𝐱= β/2∫_𝐱Π_𝐱^2   .




 §.§ Time translation invariance

For a stationary system, all observables are time-independent, ∂_t Tr[ ρ̂_t …] = 0. If additionally the Hamiltonian is of the form Ĥ[Φ̂, Π̂] = Ĥ[Π̂]+ Ĥ[Φ̂] with Ĥ[Π̂] = 1/2∫_𝐱Π_𝐱^2, then it follows that 0 = Tr[ ρ̂_t ( Φ̂_x Π̂_y + Π̂_x Φ̂_y) ]
and further ⟨φ_x π_y ⟩_W_t = 0. As a consequence the two-point function becomes block diagonal, which simplifies the inversion,

    G^(2)_c = [ ⟨φφ⟩    0;    0 ⟨ππ⟩ ]⇒Γ^(2) = [ ⟨φφ⟩^-1       0;       0 ⟨ππ⟩^-1 ]  .




 §.§ Loop expansion of the effective action

Starting from Eq. (<ref>), we calculate

	
    e^-Γ_β'[Φ]   = e^-Γ_β[Φ] + β H[Φ] = e^log Z[J(Φ)]-J_𝐱(Φ) Φ_𝐱 + β H[Φ]
       = ∫𝒟φ  e^-β H[φ] + J_𝐱(Φ) φ_𝐱-J_𝐱(Φ) Φ_𝐱 + β H[Φ]
       = ∫𝒟φ  e^-β(H[φ+Φ] - H[Φ]) +J_𝐱(Φ) φ_𝐱
       = ∫𝒟φ  e^-β(H[φ+Φ] - H[Φ] - δ H[Φ]/δΦ_𝐱φ_𝐱) +δΓ_β'[Φ]/δΦ_𝐱φ_𝐱 .


Here, we have used Eqs. (<ref>) and (<ref>), then performed a changed of variables φ→φ + Φ and finally expressed the sources as J_𝐱(Φ) = (δΓ[Φ] ) / (δΦ_𝐱) .

The first non-trivial correction (one-loop) is obtained by neglecting all terms beyond quadratic order in the fluctuating fields φ. The remaining gaussian integral can then be performed analytically, which gives

    e^-Γ'^,one-loop_β[Φ]   = ∫𝒟φ  e^-1/2φ_𝐱  G^-1_𝐱,𝐲[Φ] φ_𝐲
       = ( G^-1[Φ] )^-1/2 .

The name `one-loop' stems from the expansion in terms of the tree-level two-point function G_0. To see this, we rewrite

	
    Γ'^,one-loop_β[Φ]    =  1/2log  G^-1[Φ]  = 1/2Trlog  G^-1[Φ]  
       = 1/2Trlog(G_0 G^-1[Φ]  ) + 1/2Trlog   G_0^-1
       = -1/2∑_n=1^∞(-1)^n/nTr{ G_0 ( G^-1[Φ] - G_0^-1)}^n


where we have used the identity log A = Trlog A, employed the series expansion of the logarithm and dropped the irrelevant constant. Graphically the result can be pictured as a sum of loops consisting of lines that stand for G_0 connected by field insertions coming from ( G^-1[Φ] - G^-1_0 ).
For more details about the loop expansion, we refer to <cit.>.

Note that in the standard (unequal-time) formalism, the loop expansion is used as an expansion in weak quantum fluctuations. Here, in the context of the classical field theory limit in thermal equilibrium, it is employed as an expansion in weak thermal fluctuations. In terms of the dimensionless parameters introduced in section <ref>, the loop expansion is applicable for ϵ_th≪ 1. Thus, the tree-level approximation (i.e. the leading order in the loop expansion), which is used in the main text to extract the microscopic Hamiltonian parameters, is applicable when the following separation of scales holds:

    ϵ_q≪ϵ_th≪ 1   .

Colloquially speaking, this is the limit of weak thermal flucuations and even weaker quantum fluctuations.





§ DETAILS ABOUT THE DATA ANALYSIS



This section contains more technical details concerning the practical example, which is discussed in the main text.



 §.§ Cosine vs. Fourier transform and the boundary conditions

For an infinite system with translation invariance, the correlation functions in (Fourier) momentum space are directly related to the correlators obtained by a cosine transform. Explicitly, with the transforms

	
    φ_p    = ∫ dx   e^-ipxφ_x  , 
    φ̃_p    = ∫ dx  cos(px)φ_x = 1/2(φ_p + φ_-p)  ,


the two-point functions are related as

    ⟨φ̃_p^2 ⟩ = 1/2⟨φ_p φ_-p⟩ ,

where we assumed translation invariance, thus ⟨φ_p^2 ⟩ = 0.
Similarly, for the four-point functions, we have

    ⟨φ̃_p^4 ⟩ = 3/8⟨φ_p φ_p φ_-pφ_-p⟩  ,

where the prefactors arises from the 6 nonvanishing contributions out of 2^4 = 16 combinations. 

In practice, we deal with a finite system without periodic boundary conditions. A discrete Fourier transform is then not appropriate as it yields numerical artifacts. Since the calculation of 1PI correlators simplifies tremendously in Fourier space, we still prefer to work in a Fourier basis. Therefore, we calculate the correlators with a discrete cosine transform, which reduced the artifacts from the boundary conditions. Then we translate the correlators using the factors of 1/2 and 3/8 to Fourier-space correlators and subsequently calculate the 1PI correlation functions. For sufficiently large system sizes, this procedure yields the desired results and reduces numerical artifacts in a controlled way.




 §.§ The 1PI vertices from the numerical data

In practice, the numerical and experimental profiles live on a spatial lattice with lattice spacing Δ x and a finite number of lattice sites N, i.e. we have φ_x for x/Δ x∈{0, 1, …, N - 1}. We employ a discrete cosine transform of the individual realizations to obtain profiles φ̃_p_lat and later translate the results to the Fourier transform. The lattice momentum takes the values p_lat = j 2π/Δ x N with j ∈{ -N/2, …, N/2-1 }. We correct for some artifacts of the discrete transform at large momenta by considering physical momenta

    p_phys = 2/Δ xsin(p_latΔ x/2)   .

The two- and four-vertex densities (normalized to have units of 1/L with L=NΔ x) from the main text are obtained by


    Γ^(2)_p    = L/⟨| φ_p_phys|^2 ⟩_c , 
    Γ^(4)_p    = - L^3⟨| φ_p_phys|^4 ⟩_c/⟨| φ_p_phys|^2 ⟩_c^4  .


Here, the expectation values are


    ⟨| φ_p|^2 ⟩_c   = ⟨φ_p φ_-p⟩_c = 2⟨φ̃^2_p  ⟩_c , 
    ⟨| φ_p|^4 ⟩_c   = ⟨φ_p φ_p φ_-pφ_-p⟩_c = 8/3⟨φ̃^2_p ⟩_c ,


where the index c indicates connected correlators according to Eq. (<ref>).




 §.§ One-loop corrections



To obtain the one-loop correction to the effective action, we approximate


	
    β K[φ, Φ]    = β(H[φ+Φ] - H[Φ] - δ H[Φ]/δΦ_xφ_x )
       = λ_T/4∫_x {1/2(∂_x φ_x)^2 - 1/ℓ_J^2cosΦ_x [-1 + cosφ_x] - 1/ℓ_J^2sinΦ_x [φ_x - sinφ_x] }
    	=  1/2φ_𝐱 G^-1_𝐱,𝐲[Φ] φ_𝐲 + 𝒪(φ^3)



with G^-1_xy [Φ] = λ_T/4[-∂_x^2 + 1/ℓ_J^2cos( Φ_x)] δ(x-y). From Eq. (<ref>), we then have

    Γ'^,one-loop_β[Φ] = -1/2∑_n=1^∞1/n(λ_T/42/ℓ_J^2)^n Tr{ G_0  sin^2(Φ/2)}^n  ,

where 

    G_0,xy = ∫dp/2 π  e^ip(x-y) G_0,p
 
with . 

The one-loop corrections to the 1PI vertices are now obtained by taking derivatives of Eq. (<ref>), evaluated at Φ=0 in the symmetric case. Explicitly, we find at second order

    ΔΓ^(2),one-loop_xy   = -1/2ℓ_J^2  G_0,xxδ(x-y)

and at fourth order

    ΔΓ^(4),one-loop_xyzw   = 1/2ℓ_J^2 G_0,xxδ(x-y) δ(x-z) δ(x-w) 
       - 1/2ℓ_J^4[G_0,xy G_0,zwδ(x-z) δ(y-w) . 
               .+ 2 perm.]  .

The involved loop integrals are given by


    ∫_-∞^∞dq/2π1/q^2 + 1/ℓ_J^2   = ℓ_J/2 ,
    ∫_-∞^∞dq/2π1/q^2 + 1/ℓ_J^21/(p-q)^2 + 1/ℓ_J^2   = ℓ_J/(2/ℓ_J)^2 + p^2 ,


which results in the expressions given in the main text


    Γ^(2),one-loop_p   =  λ_T/4(p^2 + 1/ℓ_J^2)-1/4ℓ_J  , 
    Γ^(4),one-loop_p    = - λ_T/4 ℓ_J^2 - 1/8 ℓ_J^31/p^2 + 1/ℓ_J^2  .


Note that the one-loop correction is indeed of order 𝒪(ϵ_th) = 𝒪(1/Q) = 𝒪(ℓ_J /λ_T) compared to the tree-level approximation.





