In this section, we present the experimental results obtained with pure electron plasmas. Electrons are confined in a cylindrical vacuum chamber radially by a strong magnetic field along the axis of the chamber and axially by the negatively biased electrostatic potentials at both ends. When the electrostatic potential at one end of the chamber is turned off, electrons are released along the axis and strike a phosphor screen positioned perpendicular to the magnetic field. The resulting luminosity distribution is recorded by a Charge Coupled Device (CCD) camera. As this measurement is destructive, the experiment requires high reproducibility. Details of the experimental configuration are found in Refs. 9, 24, 25. The electron motion perpendicular to the magnetic field is described by the 2D Euler equation, ^{6, 7} in which the electron density is proportional to the vorticity.

The time evolution of the electron distribution is shown in Figure 2. The brightness in this figure is proportional to the electron density, i.e., the vorticity. In this experiment, we first established an uniform equilibrium electron distribution as the background vortex, and then injected an electron population with a density higher than that of the background. A region where electrons are concentrated at a density higher than the background distribution is referred to as a clump. In this case, the vorticity inside the clump is about 180 times higher than that of the background.

The clump migrates toward the center of the background vortex through the interaction of its self-induced swirling flow with the background vortex. As the clump moves, it entrains low-vorticity regions from outside the background vortex, creating a ring hole of reduced vorticity around it. We will explain the mechanism of ring hole formation using the two-point correlation function of vorticity.

The density depression around the clump is depicted in Figure 3. Figure 3 plots the density profiles along lines A-A’, B-B’, and C-C’ shown in Figures 2 (c), (d), and (e), respectively. Points A, A’, B, B’, C, and C’ are located in the background vortex, and regions with density lower than the background are highlighted in green. The green regions indicate areas exhibiting negative correlation.

In this section, we define a two-dimensional (2D) system considered and introduce prerequisite for calculating the explicit formula of a two-body correlation function for the system.

1. 2D point vortex system

Let us consider a collection of N singular point vortices, which defines a 2D vorticity field ω ̂ ( r , t ) : ω ̂ ( r , t ) ≔ ∑ i = 1 N Ω δ ( r − r i ( t ) ) . (1)

The position vector of the i -th point vortex is given by r i . Each vortex has a constant circulation Ω . The notation δ ( r ) is the 2D Dirac delta function. A quantity depending on the positions r i of the singular point vortices is indicated by ⋅ ̂ . We call a variable with ⋅ ̂ as the microscopic one. The total circulation N Ω is assumed to be finite and constant, namely Ω ∼ 1 N . (2)

A stream function is defined by ψ ̂ ( r ) ≔ ∑ i Ω G ( r , r i ) (3) where G ( r , r ′ ) is the 2D Green function for 2D Laplacian G ( r , r ′ ) ≔ − 1 2 π log | r − r ′ | . (4)

The stream function (3) and the vorticity (1) are related by the Poisson equation − Δ ψ ̂ ( r ) = ω ̂ ( r ) (5)

2. Hamiltonian

The following integral corresponding to energy of the system diverges at r = r ′ 1 2 ∫ R 2 ∫ R 2 G ( r , r ′ ) ω ̂ ( r ) ω ̂ ( r ′ ) d r d r ′ = 1 2 ∫ R 2 ψ ̂ ( r ) ω ̂ ( r ) d r ≕ H ̂ 0 + E inf , (6) where E inf is a constant including a diverging part. Notation H ̂ 0 is the Hamiltonian of the point vortex system in the usual sense. H ̂ 0 ≔ Ω 2 2 ∑ i ≠ j G ( r i , r j ) (7)

A quantity without an external perturbation (see below) is indicated by the suffix ⋅ 0 .

3. External perturbation

A macroscopic stream function ψ e ( r ) without ⋅ ̂ is introduced as an external arbitrary perturbation, which defines a microsopic stream function ψ ̂ e Ω ψ ̂ e ≔ ∫ R 2 ψ e ( r ) ω ̂ ( r ) d r = Ω ∑ i ψ e ( r i ) (8) H ̂ e ≔ H ̂ 0 + Ω ψ ̂ e (9)

The quantity H ̂ e denotes Hamiltonian with an external perturbation ψ ̂ e .

For convenience, we introduce a suffix ⋅ ∗ where the asterisk matchs either 0 or e . In an equation which contais two or more asterisks, all the asteriscs are unifiedly replaced by either 0 or e .

4. Partition function

A partition function is defined by Z ∗ ≔ ∫ R 2 N exp ( − β H ̂ ∗ ) d r 1 ⋯ d r N . (10)

According to the assumption (2), the inverse temperature β is chosen to be β ∼ O ( N ) . (11)

Under the above condition and the definition (7), the dependence of Hamiltonian H ̂ ∗ on N is estimated as β H ̂ ∗ ∼ β Ω 2 ( N − 1 ) G ∼ N ⋅ 1 N 2 ⋅ N ∼ O ( 1 ) (12) with respect to each r i . This ensures that the orders of the partition function Z ∗ and the free energy (see Eq. (21)) are estimated as Z ∗ ∼ O ( C N ) , F ∗ = − 1 β ln Z ∗ ∼ O ( 1 ) . (13)

This shows the free energy is well-defined regardless of the magnitude of N .

5. Canonical average

The canonical average converts a microscopic quantity into a macroscopic one. There are two kinds of canonical average, ⟨ ⋅ ⟩ 0 and ⟨ ⋅ ⟩ e . The notation ⟨ ⋅ ⟩ 0 represents the average with H 0 , and the notation ⟨ ⋅ ⟩ e represents the average with H e , namely, ω ∗ ( r ) ≔ ⟨ ω ̂ ( r ) ⟩ ∗ = 1 Z ∗ ∫ R 2 N ω ̂ ( r ) exp ( − β H ̂ ∗ ) d r 1 ⋯ d r N . (14)

6. Fluctuation

The fluctuation δ ω ̂ ∗ ( r ) is defined by the subtraction of the canonical averaged quantity from the microscopic one which depends on the particle position r i : δ ω ̂ ∗ ( r ) ≔ ω ̂ ( r ) − ω ∗ ( r ) (15)

7. One-body distribution function

The one-body distribution function P 1 , ∗ N ( r ) is defined by P 1 , ∗ N ( r ) ≔ ⟨ δ ( r − r i ) ⟩ ∗ . (16)

Due to the symmetry of H ̂ ∗ under the permutation of the i -th and the j -th particles, it satisfies P 1 , ∗ N ( r ) = 1 Z ∗ ∫ R 2 ( N − 1 ) exp ( − β H ̂ ∗ ( r , r 2 , ⋯ , r N ) ) d r 2 ⋯ d r N (17) for any i = 1 , ⋯ , N . Thus, the following relation also holds: ω ∗ ( r ) = ⟨ ω ̂ ( r ) ⟩ ∗ = N Ω P 1 , ∗ N . (18)

8. Two-body distribution function

The two-body distribution function P 2 , ∗ N ( r , r ′ ) is defined by P 2 , ∗ N ( r , r ′ ) ≔ ⟨ δ ( r − r i ) δ ( r ′ − r j ) ⟩ ∗ . (19) and satisfies P 2 , ∗ N ( r , r ′ ) = 1 Z ∗ ∫ R 2 ( N − 2 ) exp ( − β H ̂ ∗ ( r , r ′ , r 3 , ⋯ , r N ) ) d r 3 ⋯ d r N . (20)

Note that P 2 , ∗ N ( r , r ′ ) is not defined on r = r ′ as H ̂ 0 ( r , r ′ ) is not defined on r i = r j ( i , j = 1 , ⋯ , N ) .

The goal of this section is to derive a linear response formula for the N -point vortex system to a fluctuation (path (1) in Figure 1).

1. Linear response formula

The free energy is defined by F ∗ ≔ − 1 β ln Z ∗ . (21)

As a general procedure, we calculate a functional derivative two times to obtain a correlation of the fluctuation δ ω ̂ ( r ) _. ²⁶ The functional derivative is carried out by introducing an arbitrary external perturbation ϕ e ( r ) ≔ s δϕ e ( s ∈ R ) and differentiating with respect to s twice. d ds F e = 1 Z e ∫ R 2 N Ω ∑ i δϕ e ( r i ) exp ( − β H ̂ e ) d r 1 ⋯ d r N = ∫ R 2 ( 1 Z e ∫ R 2 N ω ̂ ( r ) exp ( − β H ̂ e ) d r 1 ⋯ d r N ) δϕ e ( r ) d r ⁽²²⁾

We represent the integrand as the functional delivative δ F e δϕ e ( r ) = 1 Z e ∫ R 2 N ω ̂ ( r ) exp ( − β H ̂ e ) d r 1 ⋯ d r N = ⟨ ω ̂ ( r ) ⟩ e = ω e ( r ) . ⁽²³⁾

Similarly, the second-order delivative is obtained: − 1 β ⋅ d 2 d s 2 F e = − 1 β ⋅ d ds ( 1 Z e ∫ R 2 N Ω ∑ i δ ϕ e ( r i ) exp ( − β H ̂ e ) d r 1 ⋯ d r N ) = − 1 Z e 2 ( ∫ R 2 N Ω ∑ i δ ϕ e ( r i ) exp ( − β H ̂ e ) d r 1 ⋯ d r N ) 2 + 1 Z e ∫ R 2 N ( Ω ∑ i δ ϕ e ( r i ) ) 2 exp ( − β H ̂ e ) d r 1 ⋯ d r N = − ∫ R 2 ( 1 Z e ∫ R 2 N ω ̂ ( r ) exp ( − β H ̂ e ) d r 1 ⋯ d r N ) δ ϕ e ( r ) d r × ∫ R 2 ( 1 Z e ∫ R 2 N ω ̂ ( r ′ ) exp ( − β H ̂ e ) d r 1 ⋯ d r N ) δ ϕ e ( r ′ ) d r ′ + ∫ R 2 ∫ R 2 ( 1 Z e ∫ R 2 N ω ̂ ( r ) ω ̂ ( r ′ ) exp ( − β H ̂ e ) d r 1 ⋯ d r N ) δ ϕ e ( r ) δ ϕ e ( r ′ ) d r d r ′ ⁽²⁴⁾

We represent this result as − 1 β δ 2 F e δϕ e ( r ′ ) δϕ e ( r ) = − 1 β δω e ( r ) δϕ e ( r ′ ) = − ⟨ ω ̂ ( r ) ⟩ e ⟨ ω ̂ ( r ′ ) ⟩ e + ⟨ ω ̂ ( r ) ω ̂ ( r ′ ) ⟩ e = ⟨ δ ω ̂ e ( r ) δ ω ̂ e ( r ′ ) ⟩ e . (25)

In the limit s → 0 , i.e., ϕ e = 0 ⋅ δϕ e = 0 , Eq. (25) is written as − 1 β δω e ( r ) δϕ e ( r ′ ) | ϕ e ( r ′ ) = 0 = ⟨ δ ω ̂ 0 ( r ) δ ω ̂ 0 ( r ′ ) ⟩ 0 (26)

As an important result in relation to Eq. (26), a response to an external perturbation defined by δω e ( r ) ≔ ⟨ ω ̂ ( r ) ⟩ e − ⟨ ω ̂ ( r ) ⟩ 0 = ω e ( r ) − ω 0 ( r ) (27) is linearly approximated by δω e ( r ) ∼ − β ∫ ⟨ δ ω ̂ 0 ( r ) δ ω ̂ 0 ( r ′ ) ⟩ 0 ϕ e ( r ′ ) d r ′ . (28)

It should be noted that Eq. (28) has the same form as the general linear-response formulae. ²⁶

2. Two-body correlation function

The two-body correlation funcntion g 0 ( r , r ′ ) in a thermal equilibrium state is defined by P 2 , 0 ( r , r ′ ) ≔ P 1 , 0 ( r ) P 1 , 0 ( r ′ ) ( 1 + g 0 ( r , r ′ ) ) . (29)

The aim of the present analysis is to determine g 0 ( r , r ′ ) in the limit of vanishing external perturbation. Characteristics of the correlation function for long-range interacting particle system was discussed by Ornstein and Zernike. ^{27, 28} A basic idea introduced by them is to split an effect of the correlation function into a self-correlation and a mutual correlation, which is called Ornstein-Zernike formula. For the point vortex system, the formula is rewritten as ⟨ δ ω ̂ 0 ( r ) δ ω ̂ 0 ( r ′ ) ⟩ 0 = Ω ω 0 ( r ) δ ( r − r ′ ) + ω 0 ( r ) ω 0 ( r ′ ) g 0 ( r , r ′ ) . (30)

The first term of Eq. (30) in the right hand side corresponds to the self-correlation and the second term the mutual correlation. Substituting Eq. (30) into (28), we obtain the formula for the linear response in the point vortex system. δω e ( r ) = − β Ω ω 0 ( r ) ϕ e ( r ) − βω 0 ( r ) ∫ g 0 ( r , r ′ ) ω 0 ( r ′ ) ϕ e ( r ′ ) d r ′ . (31)

This formula is the main conclusion, representing the linear response of the system derived from the two-body correlation function.

Mean-field theory is an approximation technique in which terms proportional to square of the fluctuation are neglected. The goal of this section is to derive a linear response formula for the mean-field approximated point vortex system.

1. Hamiltonian in the mean-field approximation

In the mean-field theory, physical quantities are represented by the sum of the mean-filed quantity and a fluctuation. There are two possible representations, with (suffix “e”) and without (suffix “0”) the external field: ω ̂ ( r ) ≔ ω ¯ ∗ ( r ) + δ ω ¯ ̂ ∗ ( r ) , (32) where the notation ⋅ ¯ indicates the quantity measured with the mean-field approximation. Namely, ω ¯ ∗ ( r ) is the vorticity in the mean-filed approximation. Mean-field vorticity is related to the stream function by the following Poisson equation. − Δ ψ ¯ ∗ ( r ) = ω ¯ ∗ ( r ) . (33)

Note that the second term in the right hand side of Eq. (32) depends on the particle coordinate r i . Thus, the term has the hat ⋅ ̂ .

The Hamiltonians H ¯ ̂ 0 and H ¯ ̂ e in the mean-filed approximation are given by H ̂ 0 + E inf = 1 2 ∫ R 2 ∫ R 2 G ( r , r ′ ) ω ̂ ( r ) ω ̂ ( r ′ ) d r d r ′ ≈ 1 2 ∬ G ( r , r ′ ) ω ¯ 0 ( r ) ω ¯ 0 ( r ′ ) d r d r ′ + ∬ G ( r , r ′ ) δ ω ¯ ̂ 0 ( r ) ω ¯ 0 ( r ′ ) d r d r ′ = − 1 2 ∬ G ( r , r ′ ) ω ¯ 0 ( r ) ω ¯ 0 ( r ′ ) d r d r ′ + Ω ∑ i ψ ¯ 0 ( r i ) ≕ H ¯ ̂ 0 ( r 1 , ⋯ , r N ) (34) H ̂ e + E inf = H ̂ 0 + E inf + ∫ R 2 ϕ e ( r ) ω ̂ ( r ) d r ≈ − 1 2 ∬ G ( r , r ′ ) ω ¯ e ( r ) ω ¯ e ( r ′ ) d r d r ′ + Ω ∑ i { ψ ¯ e ( r i ) + ϕ e ( r i ) } ≕ H ¯ ̂ e ( r 1 , ⋯ , r N ) (35)

2. Partition function

The partition functions Z ¯ e and Z ¯ 0 are defined by Z ¯ ∗ ≔ ∫ R 2 N exp ( − β H ¯ ̂ ∗ ) d r 1 ⋯ d r N (36) Z ¯ e = exp ( β 2 ∬ G ( r , r ′ ) ω ¯ e ( r ) ω ¯ e ( r ′ ) d r d r ′ ) × [ ∫ exp ( − β Ω ( ψ ¯ e ( r ′ ) + ϕ e ( r ′ ) ) d r ′ ] N (37) Z ¯ 0 = exp ( β 2 ∬ G ( r , r ′ ) ω ¯ 0 ( r ) ω ¯ 0 ( r ′ ) d r d r ′ ) × [ ∫ exp ( − β Ω ( ψ ¯ 0 ( r ′ ) ) d r ′ ] N . (38)

3. Free energy

The free energy F ¯ e and F ¯ 0 are defined by F ¯ ∗ ≔ − 1 β ln Z ¯ ∗ (39) F ¯ e = − 1 2 ∬ G ( r , r ′ ) ω ¯ e ( r ) ω ¯ e ( r ′ ) d r d r ′ − N β ln ∫ exp ( − β Ω ( ψ ¯ e ( r 1 ) + ϕ e ( r 1 ) ) ) d r 1 (40) F ¯ 0 = − 1 2 ∬ G ( r , r ′ ) ω ¯ 0 ( r ) ω ¯ 0 ( r ′ ) d r d r ′ − N β ln ∫ exp ( − β Ω ( ψ ¯ 0 ( r 1 ) ) ) d r 1 (41)

4. Thermal equilibrium state

A thermal equilibrium state is defined by minimization of the free energy. The condition is given by δ F ¯ ∗ δ ω ¯ ∗ ( r ) = 0 . (42)

Upon explicit calculation, we obtain: δ F ¯ e δ ω ¯ e ( r ) = − ∫ G ( r , r ′ ) ω ¯ e ( r ′ ) d r ′ − N β ( − β Ω ) ∫ exp [ − β Ω ( ψ ¯ e ( r ′ ) + ϕ e ( r ′ ) ) ] δ ψ ¯ e ( r ′ ) δ ω ¯ e ( r ) d r ′ ∫ exp [ − β Ω ( ψ ¯ e ( r ′ ) + ϕ e ( r ′ ) ) ] d r ′ = − ψ ¯ e ( r ) + N Ω ∫ exp [ − β Ω ( ψ ¯ e ( r ′ ) + ϕ e ( r ′ ) ) ] G ( r , r ′ ) d r ′ ∫ exp [ − β Ω ( ψ ¯ e ( r ′ ) + ϕ e ( r ′ ) ) ] d r ′ = 0 . (43)

By operating − Δ on Eq. (43), the following formula is obtained. − Δ ψ ¯ e ( r ) = ω ¯ e ( r ) = N Ω exp ( − β Ω ( ψ ¯ e ( r ) + ϕ e ( r ) ) ) ∫ exp ( − β Ω ( ψ ¯ e ( r ′ ) + ϕ e ( r ′ ) ) ) d r ′ . (44)

Similarly, the formula without a perturbation is also obtained. − Δ ψ ¯ 0 ( r ) = ω ¯ 0 ( r ) = N Ω exp ( − β Ω ψ ¯ 0 ( r ) ) ∫ exp ( − β Ω ψ ¯ 0 ( r ′ ) ) d r ′ (45)

The solution to this Poisson equation is derived here, as it will be required in subsequent sections. Substituting u 0 = − β Ω ψ ¯ 0 ( r ) (46) into Eq. (45), we get − Δ u 0 = λ e u 0 ∫ e u 0 ( r ′ ) d r ′ ( = − β Ω ω ¯ 0 ) , (47) where λ ≔ − βN Ω 2 . (48)

The solution to Eq. (47) in R 2 was discussed by Chen and Li. ²⁹ The solution to Eq. (47) exists only when λ = 8 π and the solution is given by u 0 ( r ) = ln 8 ϵ 2 ( 1 + | ϵ r | 2 ) 2 − ln 8 π Z (49) and its parallel translation with the parameter Z = ∫ e u 0 ( r ′ ) d r ′ = ∫ e − β Ω ψ ¯ 0 ( r ′ ) d r ′ , ϵ > 0 . (50)

Thus, the solution of Eq. (45) is given by the following formulae and those of the parallel translation: ψ ¯ 0 ( r ) = − 1 β Ω u 0 ( r ) = − 1 β Ω ln 8 ϵ 2 ( 1 + | ϵ r | 2 ) 2 + 1 β Ω ln 8 π Z (51) ω ¯ 0 ( r ) = − Δ ψ ¯ 0 ( r ) = ( − 1 β Ω ) ( − Δ u 0 ) = − 1 β Ω 8 ϵ 2 ( 1 + | ϵ r | 2 ) 2 (52)

5. Response in the mean-field theory

The responses in the mean-field theory is defined by δ ω ¯ e ( r ) ≔ ω ¯ e ( r ) − ω ¯ 0 ( r ) , (53) δ ψ ¯ e ( r ) ≔ ψ ¯ e ( r ) − ψ ¯ 0 ( r ) (54)

Using Eqs. (53) and (54), Eq. (44) is linearlized. δ ω ¯ e ( r ) = N Ω exp ( − β Ω ( ψ ¯ e ( r ) + ϕ e ( r ) ) ) ∫ exp ( − β Ω ( ψ ¯ e ( r ′ ) + ϕ e ( r ′ ) ) ) d r ′ − N Ω exp ( − β Ω ψ ¯ 0 ( r ) ) ∫ exp ( − β Ω ψ ¯ 0 ( r ′ ) ) d r ′ ≈ N Ω exp ( − β Ω ψ ¯ 0 ( r ) ) ( − β Ω ) ( δ ψ ¯ e ( r ) + ϕ e ( r ) ) ∫ exp ( − β Ω ψ ¯ 0 ( r ′ ) ) d r ′ − N Ω exp ( − β Ω ψ ¯ 0 ( r ) ) { ∫ exp ( − β Ω ψ ¯ 0 ( r ′ ) ) d r ′ } 2 × ∫ exp ( − β Ω ψ ¯ 0 ( r ′ ) ) ( − β Ω ) ( δ ψ ¯ e ( r ′ ) + ϕ e ( r ′ ) ) d r ′ = − β Ω ω ¯ 0 ( r ) ( δ ψ ¯ e ( r ) + ϕ e ( r ) − 1 N Ω ∫ ω ¯ 0 ( r ′ ) ( δ ψ ¯ e ( r ′ ) + ϕ e ( r ′ ) ) d r ′ ) (55)

Equation (55) is the final formula of the response to the external field in the mean-field theory.

We have otaintained two formulae for the response, Eq. (31) and Eq. (55). We assume that δω e ( r ) , ω 0 ( r ) in Eq. (31) coincde with δ ω ¯ e ( r ) , ω ¯ 0 ( r ) in Eq. (55), respectively: − β Ω ω ¯ 0 ( r ) ϕ e ( r ) − β ω ¯ 0 ( r ) ∫ g 0 ( r , r ′ ) ω ¯ 0 ( r ′ ) ϕ e ( r ′ ) d r ′ = − β Ω ω ¯ 0 ( r ) ( δ ψ ¯ e ( r ) + ϕ e ( r ) − 1 N Ω ∫ ω ¯ 0 ( r ′ ) { δ ψ ¯ e ( r ′ ) + ϕ e ( r ′ ) } d r ′ ) , ↔ 1 Ω ∫ g 0 ( r , r ′ ) ω ¯ 0 ( r ′ ) ϕ e ( r ′ ) d r ′ = δ ψ ¯ e ( r ) − 1 N Ω ∫ ω ¯ 0 ( r ′ ) { δ ψ ¯ e ( r ′ ) + ϕ e ( r ′ ) } d r ′ (56)

Applying the operator − Δ r + β Ω ω ¯ 0 ( r ) to Eq. (56), where Δ r denotes the Laplacian with respect to r , we obtain 1 Ω ∫ { − Δ r + β Ω ω ¯ 0 ( r ) } g 0 ( r , r ′ ) ω ¯ 0 ( r ′ ) ϕ e ( r ′ ) d r ′ (57) = δ ω ¯ e ( r ) + β Ω ω ¯ 0 ( r ) ( δ ψ ¯ e ( r ) − 1 N Ω ∫ ω ¯ 0 ( r ′ ) { δ ψ ¯ e ( r ′ ) + ϕ e ( r ′ ) } d r ′ ) (58)

Substituting Eq. (55) into Eq. (58), the following equation is obtained. ∫ { − Δ r + β Ω ω ¯ 0 ( r ) } g 0 ( r , r ′ ) − β Ω 2 ω ¯ 0 ( r ′ ) ϕ e ( r ′ ) d r ′ = ω ¯ 0 ( r ) ϕ e ( r ) (59)

As ϕ e ( r ) is an arbitrary function, Eq. (59) is reduced to − ( Δ r − β Ω ω ¯ 0 ( r ) ) g 0 ( r , r ′ ) − β Ω 2 = δ ( r − r ′ ) . (60)

In the next section, we will solve Eq. (60) and obtain an explicit formula for g 0 ( r , r ′ ) .

1. Homogeneous solution

Now, we are ready to solve Eq. (60). At first, we will obtain the homogeneous solution to Eq. (60), namely, the solution to − ( Δ r − β Ω ω ¯ 0 ( r ) ) g 0 ( r , r ′ ) − β Ω 2 = 0 (61) where the explicit formula for ω ¯ 0 ( r ) is given in Eq. (52). The partial derivative of Eq. (47) with respect to ϵ at fixed Z is given by − Δ ∂ u 0 ∂ϵ = ∂ ∂ϵ ( λ e u 0 ∫ e u 0 ) = λ e u 0 ∫ e u 0 ∂ u 0 ∂ϵ = − β Ω ω ¯ 0 ∂ u 0 ∂ϵ (62)

Here, we have used the following relation with fixed Z . 0 = ∂ Z ∂ϵ = ∂ ∂ϵ ∫ e u 0 = ∫ e u 0 ∂ u 0 ∂ϵ (63)

Equation (62) is reduced to − ( Δ − β Ω ω ¯ 0 ( r ) ) ∂ u 0 ∂ϵ = 0 (64) and the homogeneous solution G h ( r ) to Eq. (61) is obtained G h ( r ) − β Ω 2 = ∂ u 0 ∂ϵ = 2 ϵ 1 − | ϵ r | 2 1 + | ϵ r | 2 (65) where the explicit formula for u 0 is given by Eq. (49).

2. Special solution

In the next we will obtain a special solution to Eq. (60) with r ′ = 0 .

Consider Eq. (47) with a perturbation cδ ( r ) whose solution is denoted by u c ( r ) : − Δ u c ( r ) = λ e u c ( r ) ∫ e u c ( r ′ ) d r ′ + cδ ( r ) (66) where c is a parameter. Equation (66) is reduced to (47) in the limit of c → 0 . A solution of (66) has been discussed by J. Prajapat and G. Tarantello. ³⁰ In the case of sufficiently small | c | , the solution of Eq. (66) in R 2 exists only when λ = 8 πn (67) n ≔ 1 − c 4 π (68) and is given by u c ( r ) = c 2 π ln | r | − 1 + ln 8 ϵ 2 n ( 1 + | ϵ r | 2 n ) 2 − ln 8 π nZ . (69) where Z = ∫ e u c ( r ′ ) d r ′ , ϵ > 0 . (70)

The limit c → 0 of the partial derivative of Eq. (66) with respect to c at fixed Z is given by − Δ ∂ u c ∂ c | c = 0 = ∂ ∂ c ( λ e u c ∫ e u c ( r ′ ) d r ′ ) | c = 0 + δ ( r ) = − 2 e u 0 ∫ e u 0 ( r ′ ) d r ′ + 8 π e u 0 ∫ e u 0 ( r ′ ) d r ′ ∂ u c ∂ c | c = 0 + δ ( r ) (71)

The above equation is reduced to the following form: − ( Δ − β Ω ω ¯ 0 ) ( ∂ u c ∂ c | c = 0 − 1 4 π ) = δ (72)

By employing the practical form of ∂ u c ∂ c | c = 0 − 1 4 π = 1 2 π ( 1 − | ϵ r | 2 1 + | ϵ r | 2 ln 1 | ϵ r | − 1 ) , (73) we obtain the special solution G s ( r ) to Eq. (60) G s ( r ) − β Ω 2 = ∂ u c ∂ c | c = 0 − 1 4 π = 1 2 π ( 1 − | ϵ r | 2 1 + | ϵ r | 2 ln 1 | ϵ r | − 1 ) . (74)

We conclude that the complete solution G ( r ) to Eq. (60) is given by the combination of the homogeneous solution (65) and the imhomogeneous solution (74), namely, G ( r ) = G s ( r ) + α G h ( r ) (75) where α is an arbitrary constant.

In general, the correlation function satisfies the following two conditions: g ( x 0 , 0 ) = 0 ( | x 0 | ≫ 1 ) (76) ∫ R 2 w ( x ) g ( x , 0 ) d x = 0 (77) where w ( x ) is a one-body distribution function in an equilibrium state. The reason why the second condtion (77) holds is discussed in Appendix A. In this section, we examine if the obtained distribution function (52) and correlation function (75) satisfies the above condition.

To derive the correlation function, we have used an mean field approximation. Mean-field approximation generally fails in systems with strong correlations. That is, we consider that the first condition is already fulfilled.

Let us examine the second condition. One-body distribution function in an equilibirium state is given by Eq. (52) ω ¯ 0 ( x ) = ϵ 2 π 1 ( 1 + | ϵ x | 2 ) 2 (78) where the coefficient is determined to hold the normalization condition ∫ R 2 ω ¯ 0 ( x ) d x = 1 . (79)

Integrating Eq. (75) in R 2 , we obtain ∫ R 2 ω ¯ 0 ( x ) g ( x , 0 ) d x = − β Ω 2 ∫ R 2 ϵ 2 π 1 ( 1 + | ϵ x | 2 ) 2 × { 1 2 π ( | ϵ x | 2 − 1 | ϵ x | 2 + 1 ln | ϵ x | − 1 ) + α | ϵ x | 2 − 1 | ϵ x | 2 + 1 } d x (80)

Perform a change of variables in the integral ϵ x = y , d x = 1 ϵ 2 d y (81) and continue the calculation. = − β Ω 2 ∫ R 2 1 π 1 ( 1 + | y | 2 ) 2 × { 1 2 π ( | y | 2 − 1 | y | 2 + 1 ln | y | − 1 ) + α | y | 2 − 1 | y | 2 + 1 } d y (82)

Axial symmetry is assumed. = − β Ω 2 ∫ 0 ∞ 1 π 1 ( 1 + r 2 ) 2 × { 1 2 π ( r 2 − 1 r 2 + 1 ln r − 1 ) + α r 2 − 1 r 2 + 1 } 2 πr d r (83)

Perform a change of variables in the integral r 2 = s , 2 rdr = ds (84) and continue the calculation. = − β Ω 2 ∫ 0 ∞ 1 ( 1 + s ) 2 { 1 2 π ( s − 1 s + 1 1 2 ln s − 1 ) + α s − 1 s + 1 } ds = − β Ω 2 { 1 4 π ∫ 0 ∞ s − 1 ( s + 1 ) 3 ln sds − 1 2 π ∫ 0 ∞ ds 1 + s + α ∫ 0 ∞ s − 1 ( s + 1 ) 3 ds } = β Ω 2 4 π (85)

We find that the α -dependence does not remain in the integral result. To overcome this difficulty, cutoff radius R is introduced. According to this change, we replace the integral domain R 2 with B R ( 0 ) and search the value of α which satisfies G R ( x ) ≔ G s ( x ) + α G h ( x ) (86) ∫ B R ( 0 ) ω ¯ 0 ( x ) G R ( x ) d x = 0 (87)

Returning to Eq. (85), we re-evaluate the above integral. ∫ B R ( 0 ) ω ¯ 0 ( x ) G R ( x ) d x = − β Ω 2 { 1 4 π ∫ 0 ( ϵR ) 2 s − 1 ( s + 1 ) 3 ln sds − 1 2 π ∫ 0 ( ϵR ) 2 ds ( s + 1 ) 2 + α ∫ 0 ( ϵR ) 2 s − 1 ( s + 1 ) 3 ds } (88) = − β Ω 2 { 1 4 π ( − ( ϵR ) 2 ( ( ϵR ) 2 + 1 ) 2 ln | εR | + ( ϵR ) 2 ( ϵR ) 2 + 1 ) − 1 2 π ( ϵ R ) 2 ( ϵ R ) 2 + 1 − α ( ϵ R ) 2 ( ( ϵ R ) 2 + 1 ) 2 } = 0 (89)

Thus, we obtain α satisfying condition (87). α = − 1 4 π ln ( ϵ R ) 2 − 1 4 π ( ( ϵ R ) 2 + 1 ) (90)

By introducing a conditional probability that assumes the presence of a particle at the origin as the source of perturbation, the two-body distribution function effectively reduces to a one-body distribution function. Applying the above condition to Eq. (29) yields P 2 , 0 ( 0 , r ′ ) = 1 ⋅ P 1 , 0 ( r ′ ) ( 1 + g 0 ( 0 , r ′ ) ) = P 1 , 0 ( r ′ ) ( 1 + G ( r ′ ) ) (91)

Among the terms in the above expression, P 1 , 0 ( r ′ ) G ( r ′ ) corresponds to the linear response part, where G ( r ′ ) is defined in Eq. (75). Figure 4 shows a plot of this function. As anticipated, a region with negative correlation forms in the vicinity of the origin in response to the perturbation imposed at the origin. It should be noted that since the term G ( r ′ ) is proportional to Ω ∼ O ( 1 / N ) , the magnitute of the response is small compared to the leading-order term of 1. Thus, although we successfully obtained the two-body correlation function G ( r ) analytically, its effect is small and only qualitative agreement is achieved. In other words, additional quantitative analysis may be required to refine the explanation for the formation of the depleted region observed in the experiments.

Further efforts are necessary to quantitatively explain the experimental results. The two-body correlation function that we analytically derived in this work remains within the linear response regime and does not address the longer-time dynamics. Since analytical treatment of the time evolution of the system is expected to be difficult, we are preparing large-scale numerical simulations using a supercomputer equipped with the latest PEZY-SC4s processor, a platform that one of the authors (Y.Y.) has been utilizing. Numerical simulations will enable us to capture the detailed particle transport behavior, allowing us to investigate in detail aspects such as the destination of particles during the formation of the depleted region.