Collapse of n Point Vortices, Formation of the Vortex Sheets and Transport of Passive Markers

: In this paper, the motion of the n-vortex system as it collapses to a point in ﬁnite time is studied. The motion of vortices is described by the set of ordinary differential equations that we are able to solve analytically. The explicit formula for the solution demands the initial location of collapsing vortices. To ﬁnd the collapsing locations of vortices, the algebraic, nonlinear system of equations was built. The solution of that algebraic system was obtained using Newton’s procedure. A good initial iterate needs to be provided to succeed in the application of Newton’s procedure. optimization procedure was used to ﬁnd such a good initial iterate. The numerical studies were conducted, and numerical evidence was presented that if in a collapsing system n = 50 point vortices include a few vortices with much greater intensities than the others in the set, the vortices with weaker intensities organize themselves onto the vortex sheet. The collapsing locations depend on the value of the Hamiltonian. By changing the Hamiltonian values in a speciﬁc interval, the collapsing curves can be obtained. All points on the collapse curves with the same Hamiltonian value represent one collapsing system of vortices. To show the properties of vortex sheets created by vortices, the passive tracers were used. Advection of tracers by the velocity induced by vortices was calculated by solving the proper differential equations. The vortex sheets are an impermeable barrier to inward and outward ﬂuxes of tracers. Arising vortex structures are able to transport the passive tracers. In this paper, several examples showing the diversity of collapsing structures with the vortex sheet are presented. The collapsing phenomenon of many vortices, their ability to self organize and the transportation of the passive tracers are novelties in the context of point vortex dynamics.


Introduction
Studies of the discrete vortex systems began in 1858 with the publication of the Helmholtz paper [1], where the interaction of two linear vortices was investigated.
Nowadays, the dynamics of the point vortex system still rouse interest among scientists. Good introductory reference material to the subject can be found in [2][3][4][5][6][7][8]. The concept of replacing the continuous field of vorticity by the distribution of δ-Dirac functions, called point vortices, constitutes the foundation of the vortex methods and provides a very useful numerical methodology to study the inviscid or viscous flow problems [9][10][11]. The set of point vortices that was obtained by replacing the continuous vorticity led to the generation of a velocity field that approximates the solution of the Euler equation [12] (refer to Appendix A). Essential ingredients of modern modeling of any fluid flow are studies of evolution of vorticity using the point vortices [6,10,13,14]. Sometimes mathematical analysis of a collection of a few point vortices' dynamics can shed light on interesting features of fluid motion, which is the reason why scientists remain highly interested in the point vortex system's dynamics. The most fascinating phenomenon related to the collapse of the system of vortices is the point in finite time. This phenomenon is known for a long time since the dissertation of Gröbli 1877 [15]. It was rediscovered once again in the 70s by Aref [16] and Novikov [17,18]. Due to the fact that the distance between vortices changes during the collapse motion, different scales of motion arise, which is the elementary feature of turbulence. The collapse of vortices is regarded as an elementary act in two-dimensional turbulence kinetics [19]. Since Aref's and Novikov's papers were published, many scientists have studied the collapse of the three-vortices problem and its applications; for example, [8,18,[20][21][22][23][24]. The number of papers devoted to the collapsing of three vortices indicates the importance of the collapsing phenomenon, even though only three vortices were taken into account. The author's papers showed numerically by several examples (references [25,26] that the collapse of vortices is possible for any number of vortices. In this paper, numerical results related to some special feature of collapse system of vortices are presented. It was discovered that when one or a few vortices in the collapse system are much larger than others (such an inhomogeneity typical in turbulent motion [27,28]), then the vortices with the smallest intensity always organize themselves into the vorticity threads. Weaker vortices are gathered along the curve segments, which can be regarded as the vortex sheets. The sheets behave as a solid barrier and are able to transport the passive markers en route to the singularity point. As stated in [27], the coherent vortex structures that appeared in the turbulent rotating flows are able to trap the passive tracer. The boundary of the coherent vortex structure is impermeable for tracers from the inside out and vice versa. It is evident from the numerical experiments that a collapse system of vortices that possess strong vortices organizes itself into a coherent vortex structure.
In [29] the numerical algorithm for the collapse of the vortex sheets accompanied by a strong vortex was given. In that paper, the shape of vortex sheet was initially presumed and then it was replaced by identical weak point vortices. To find the collapse locations, the algebraic system resulting from the point vortex approximation of the Rott-Birkhoff equation for vortex sheet, together with some invariant of motion was solved.
In present paper, it is shown numerically that a collapse system with one or a few strong vortices organizes itself into a vortex sheet. The self-organization of vortices generates interesting vortex structures that are able to transport the passive tracers. The vortex sheet essentially modifies the advection of passive tracers. Passive particle advection is a crucial component in many geophysical processes, such as the dispersion of pollutants in the atmosphere and in the ocean [7].
The numerical results here are a novelty in the field of the point vortex system dynamics. In this paper, the description of the numerical algorithm is provided, followed by several examples of different shapes of collapsing structures of vortices with one, two and four strong vortices. Passive tracers were used to demonstrate the transporting features of collapsing vortex structures.

Equations of Motion for n Vortices
Equations of motion for the system of n-point vortices were first provided by Helmholtz in 1858 [1] (refer to Appendix A). Since then, his equations played a key role in vortex methods and were used in multiple applications of fluid motion modeling. Here, Helmholtz equations are used to study the dynamics of point vortices that will collapse to a point in a finite time. Details of assumptions and methodology are provided below and in Sections 3 and 4.
First, let us assume that there are n point vortices on the plane with distinct positions z k , and circulations (intensities) Γ k , each Γ k ∈ R \ 0. Each of distinct positions can be expressed in terms of complex number z k = x k + iy k . Using the complex variables representation, Helmholtz equations of motion for an n-point vortex system can be expressed as [30][31][32] (see Appendix A): where the prime on the summation indicates the omission of the term with j = k, and the overline denotes complex conjugation; v k is the velocity of the kth vortex. System (1) is associated with the useful identities [3,30] : Using the above identities, it is not hard to obtain the integrals of motion for Equation (1) [3,32]: where (5) is a Hamiltonian, (6) and (8) are a linear and an angular impulse respectively and invariant V (Equation (7)) is a virial. Invariant (5) is expressed through the mutual distance between the vortices. It permits one to express the equation of motion (1) as a Hamiltonian system [32]: By substituting z k in Equation (1) with new variableẑ k = az k + b, it can be observed that the shape of trajectories stays unchanged apart from time, which is scaled: t/|a| 2 . Without lost of generality, one can shift the center of a vorticity's mass to the origin of the coordinate system, (0, 0), assuming that σ = ∑ k Γ k = const. = 0 [26]. Consequently, we have M x = 0, M y = 0 and S = 0 [3].

Self Similar Motions of n-Vortices
In this article, the following definition of the self-similar collapse motion of the nvortices was introduced: [25,26]. Definition 1. The system of the n-vortices is in a self-similar collapse motion if there exist a complex function λ(t), λ(t) = λ r (t) + iλ i (t), Re(λ) = λ r (t) < 0 and λ i (t) = 0, and the equation of motion takes the form Let us assume that the solution of Equation (1) allows self-similar motion, and let us introduce the new variables (r(t), ϕ(t)), such that the solution takes form: By inserting Equation (11) into (1), the variables (r(t), ϕ(t)) can be easily separated, allowing the following representation: Integration of (12) provides: which allows us to represent the solution of Equation (11) as [25,26] (refer also to [33]): Solution (14) represents the logarithmic spiral. The collapse time is given by Using the expression for collapse time (15) one can rewrite the (14) as follows: It is easy to check that when t → T c ; then z k (t) → 0, r(t) → 0 and ϕ(t) → +∞ when λ i (0) > 0 or ϕ(t) → −∞ when λ i (0) < 0. When λ r (t) > 0 then the system expands. To change collapsing system it is enough to change the signs of all vortex circulations to the opposite ones [32,34]. When the real part of λ(t) equals zero, λ r = 0, then the vortices are in relative or absolute equilibrium. The entire system can rotate like a solid body or take stationary positions.
In the collapsing system of vortices, distances between vortices go to zero, seemingly contradicting the Hamiltonian constant. However, by inserting Equation (11) into to Equation (5) we obtain From Equation (17), it can be concluded that Hamiltonian (5) during the self-similar motion is constant, when the virial V = 0.

Algebraic Equations for Collapse Locations of Vortices
A configuration space of n vortices requires 2n real numbers. From the definition of the self-similar motion (10), we have v 1 z k = v k z 1 , which provides only n − 3 independent relations [3,26]. This leads to 2n − 6 algebraic equations of the form: f j = v j z j+1 − v j+1 z j = 0, for real and imaginary parts of f j . In addition, we have three more equations from the invariants: f n−5 = M x = 0, f n−4 = M y = 0, and f n−3 = S = 0 and one more equation from property of identity (2), it is f n−2 = Re(∑ n i=1 Γ i v i ) = 0 [3]. This leads to the nonlinear system of n − 2 algebraic equations ( f 1 , f 2 , . . . , f n−2 ) T = 0 where each has the form f j (x 1 , y 1 , . . . , x n−1 , y n−1 ) = 0.
To further reduce the number of equations of the vortex, z n = (x n , y n ) was fixed. The nonlinear system of 2n − 2 equations was solved using the Newton method (refer to [25,26]). To succeed in a numerical solution of an algebraic nonlinear system of equations, one should provide a good guess for the initial iterate [35]. To obtain such an initial iterate, an unconstrained Levenberg-Marquart optimization procedure was used to minimize the sum: The vector (x, y) min that minimizes (18) was then used as an initial iterate for Newton's procedure (refer to [25,26]). All calculations were performed using the Mathematica © v.12. To succeed in the minimization of (18), one should perform the calculations with high digital precision. In Mathematica © the digital precision is fixed by including parameter called "WorkingPrecision" in the procedure. By setting WorkigPrecision → ω, calculations were perform with ω-digits of internal computation [36] (p. 409). In all case studies, ω was set to 500. The minimum threshold of function F was ≈ 10 −700 or smaller. All invariants (M x , M y , S, V) were on the level ≈ 100 −500 . Typically, the Newton procedure demands only one iteration.

Prelude Example
First, we present an example that highlights the difference when the collapsing system of 50 vortices has and does not have a strong vortex. Both systems have a fixed vortex in the same place, marked by navy-blue dot with a dashed circle. Both systems start from a randomly selected initial location. Figure 1a,b present calculated locations of a collapsing set of vortices with intensities Γ 1−25 = 1 and Γ 26−50 = −4/3. Vortices with positive circulation (red points) and negative circulation (blue points) are laid without visible order. In Figure 1c, vortices were in the same initial starting locations, as in Figure 1a, but the intensities of vortices were different. It was assumed that vortices from 1 to 49 have equal intensity Γ 1−49 = 1, and that there is one strong vortex with intensity Γ 50 = −24. Red dots represent the positive vortices with intensity equal to one, and blue, thick dots with dashed circles around the represent the negative strong vortex.
The locations of collapsing vortices were obtained by first applying the unconstrained Leverber-Marquart optimization procedure (18) and then Newton's procedure. Computed locations are presented in Figure 1c. It can be observed that a large number of vortices with intensity equal to one lay very closely to each other on the arc of a curve. This arc can be treated as a vortex sheet. As will be shown in the following examples, it is possible to obtain a solution for which all vortices with small intensities (small in absolute value) lay on the same curve. The appearance of the a vortex sheet is a generic feature, which can be observed in a collapsing vortex set if one or several strong vortices are present in the collapsing set.
(a) In the following sections only selected examples out of many numerical case studies are presented. These examples were chosen to present interesting properties of collapsing vortices. Further, to make presented results easily reproducible, the initial locations of vortices were placed on a well defined curve. In Sections 5.1.1 and 5.1.2, additional details regarding the dynamics of the passive tracer and the selection of intensity values for strong vortices are provided.

Dynamics of Passive Tracers
The passive tracers were used to demonstrate how the vortex sheets impact the motion of fluid and passive tracers in the collapse of vortices. Passive tracers do not affect the motion vortices [37]. They move according to the velocity induced by vortices. Due to the fact that the locations of vortices are known in time (refer to (14)), the vortex-induced velocity field (u(x, y, t), v(x, y, t)) can be computed. The trajectory of motion of passive tracers x(t) and y(t) can be then obtained by solving the following differential equations: where ψ(x, y) is a stream function that with the intensity Γ j can be expressed as The differential equations were solved using embedded Bogacki-Shampaine 5(4) Runge-Kutta method implemented in Mathematica [38].

Remarks on the Determination of the Values of Strong Vortices
In this paper, the calculations for one, two and four strong vortices inside the collapsing vortices with uniform intensities are presented. Here, we explain how to determine what value of intensity the strong vortex should have to fulfill the condition for virial V = 0 (refer to (7) and (17)). For simplicity, steps taken are described using the example. Here, the intensities of 49 vortices were assumed to be Γ 1−49 = 1 and the intensity of a strong vortex was assumed be Γ 50 = a. By inserting information about intensities of vortices into Equation (7) and using the ability for symbolic calculations in Mathematica © , the following expression for viral V can be obtained V = 1176 + 49a = 0. It is easy to see that in this case a = −24. In the same way, we can calculate the values of intensity for two, three and four strong vortices. For example, for two strong vortices with the same intensity a and the remaining vortex intensities Γ 1−48 = 1, the virial calculated by Mathematica can be written as V = 1128 + 96a + a 2 = 0. Such a representation yields two possible solutions for intensities Γ 49 = Γ 50 = a, a = 2(−24 − 7 √ 6) or a = 2(−24 + 7 √ 6). If the system has four strong vortices instead of two, the equation for virial takes form V = 1035 + 184a + 6a 2 = 0, which yields two possible solutions for intensities Γ 47 = Γ 48 = Γ 49 = Γ 50 = a, i.e., a = 1/6(−92 − 7 √ 46) or a = 1/6(−92 + 7 √ 46).

Numerical Results for the Collapse of 50 Vortices with One Strong Vortex
In this example, the 49 vortices with unit intensity Γ 1 = · · · = Γ 49 = 1 were assumed to be initially located on the arc of the circle. For vortex i, i = 1, . . . 49, the initial position on the arc was defined as (y and R = 62, x r = 10, x l = −10, y 0 = 60 (green dots in Figure 2a). The initial location of a fixed strong vortex with intensity Γ 50 = −24 was (x 50 , y 50 ) = (−6, 1/2)-navy-blue dot in dashed circle in Figure 2a,b. As pointed out in Section 5.1, the initial locations of 50 vortices were selected so that anyone can reproduce the example without difficulties. The collapsing set of vortices are marked by black dots, and one fixed vortex was marked by a navy-blue dot with a dashed circle in Figure 2a. Due to the fact that the Hamiltonian is expressed by mutual distance between the vortices, replacing the identity (refer to Section 4)       Figure 5 shows the evolution of passive tracers that were put around the strong vortex (blue color) and that were put on the other side of the vortex sheet (green color). Evidently, the blue tracers approached the vortex sheet but did not penetrate to the other side of the sheet. It can be observed that blue tracers rotate violently around the strong vortex and they are carried over time with the whole system to the center of vorticity. By contrast, the green tracers are calmly moved by advection. The strong vortex together with vortex sheet rotate quicker than green tracers. Due to the process of collapse, the vortex sheet shrinks and becomes shorter, and it gradually loses the blue tracers. It can be also be observed that over time, the vortex sheet catches the green patch of tracers. The vortex sheet initially only sliced the green patch, but over time the patch of green tracers was stroked and smashed.

Collapse of 50 Vortices with Two Strong Vortices
The behaviour of the system of collapsing vortices with more than one strong vortex is more complex and more interesting. In Sections 6.1-6.4, several selected case studies of the system of vortices with two strong vortices are presented. In Section 6.5 a case study for the system with four strong vortices is presented.
6.1. The Case A 1 Figure 6a displays an example of the initial locations (green dots and two thick navyblue dots) and locations of the collapsing set of vortices (red dots together with two thick navy-blue points). Figure 6b displays the selected trajectories of collapse vortices. The 48 vortices with uniform intensity (Γ 1−48 = 1) were laid along the same arc, as in the example presented in Figure 2a. The two strong vortices with intensities Γ 49 = Γ 50 = 2(24 + 7 (6)) were located at (x 49 , y 49 ) = (8, 0) and (x 50 , y 50 ) = (−8, 0) (navy-blue, thick dots) (Figure 6a). It can be noticed that the collapsing locations of vortices resemble sinusoidal-like functions with antisymmetric positions of vortices. Peculiar to note is the strong vortex at location (x 49 , y 49 ) = (8, 0) did not move afterthe application of the minimization of (18) and Newton's procedure, as if that vortex was in the correct location to begin with. Figure 7 shows the evolution of vortex locations when the Hamiltonian values change in the interval H ∈ [3.59, 99.65]. The collapse time and Hamiltonian value for each example are provided at the top of each frame. Noticeably, the wavy sinusoidal-like function of vortices' locations at H ≈ 3.59 straightens up gradually as the Hamiltonian value changes to H ≈ 99.65 and vortices take the positions along the vertical axis y. The system seems to be near the relative equilibrium at collapse time, T c ≈ 25170.   (6)). It can be observed that vortices lay in a nearly straight line that creates a vortex sheet, and two sets of tracers are located on opposite sides. Over time, when the system is collapsing, the length of vortex sheet gets shorter, and the tracers are gradually released into the surroundings. Collapse time is T c ≈ 104.22.

The Case A 2
In this case study, the initial locations of vortices with intestines Γ 1−48 = 1 were the same as in case study (Section 6.1), i.e., along the arc, but the locations of two strong vortices were changed. Here, the two strong vortices were placed at (x 49 , y 49 ) = (4/5, 7/2) and (x 50 , y 50 ) = (−8, −1) (fixed vortex position). Figure 9a displays starting locations of weak vortices (black dots) and two strong vortices were marked by navy-blue and light-blue dots. Figure 9a displays also the collapsing locations: the weaker, positive vortices (red dots) and two navy-blue dots as strong vortices with negative intensities. It can be observed that the second strong vortex moved farther away from its initial position (x 49 , y 49 ) = (4/5, 7/2) to (x 49 , y 49 ) ≈ (−1.684, 9.927) and that two disjoint pieces of vortex sheet were formed. Each piece of vortex sheet surrounds different strong vortices. Figure 9b presents the collapse trajectories of the system. The collapse time in this case is T c ≈ 12.559; the Hamiltonian value is H ≈ 23.145.. Figure 10 presents   Figure 11 presents the evolution of three patch of tracers in the velocity fields induced by the collapsing vortices. Two patches of tracers (green and orange) were placed around strong vortices and one patch of tracers was placed on the right site of vortex sheets. The initial locations of vortices for t = 0 correspond to the third frame in Figure 10 for the Hamiltonian value H ≈ 37.77. Red dots, represent the positive, weak vortices, Γ 1−48 = 1, and navy-blue dots represent two strong vortices with intensities Γ 49 = Γ 50 = 2(24 + 7 (6)). It can be observed that over time, green and orange tracers rotate quickly around the strong vortices and they are carried with the whole system toward the center of vorticity (0, 0). By contrast, the gray tracers are slowly convected around. The strong vortices together with the vortex sheets rotate faster than gray tracers. However, over time the vortex sheets catch the gray tracers. Initially, the vortex sheets only sliced the patch of gray tracers into two, but in the next turns of the system, patch of the gray tracers were torn off. The branch shaped vortex sheet behaves as a solid barrier that limit the influence of strong vortices. Just behind that barrier, in the direction opposite to the locations of strong vortices, the values induced by all vortices velocities are much smaller than in front of the sheets. At the vortex sheets, normal velocity is zero. It means that flux of mass from inside in and inside out is also zero. The bar of weaker vortices between the strong vortices separates the green and orange tracers, creating two compartments. As the vortices collapse, the vortex sheet shrinks and becomes shorter, gradually losing the green and orange tracers. The collapse time for this system is T c = 75.51.

The Case A 3
In this case study presented in Figure 12a, the initial starting locations of 48 vortices with weaker intensities Γ 1−48 = 1 were laid along the same arc as in the example presented in Figures 6a and 9a. In this example, two strong vortices with intensities Γ 49 = Γ 50 = 2(24 + 7 (6)) were located in (x 49 , y 49 ) = (5,4) and (x 50 , y 50 ) = (−6, 0) (navy-blue dots in Figure 12a). Figure 12b displays the collapse locations that were calculated using Newton's procedures (red dots, two navy-blue dots). It can be observed that some of weaker vortices took approximately the locations on the perimeter of the ellipse (dashed line in Figure 12b) and the remaining ones created a short line on y axis, almost symmetric around the origin of the coordinate system that divides the inner space of ellipse into two compartments. The exemplary collapse trajectory is shown in Figure 12b. Figure 13 presents the shift in collapsing locations when the Hamiltonian values were changed from H ≈ −192.57 (see left figure in Figure 13) to H ≈ −202.81 (right ellipse in Figure 13). The collapse times T c =≈ 26.45 and T c =≈ 26131.9, respectively. The vortices that are at the perimeter of the ellipse on the right side of Figure 13 are close to relative equilibrium, and the distribution of the vortices is nearly symmetrical in terms of y axis. Noticeably, one of the vortices lost its connection with right part of the vortex sheet and took a position beneath the vortex bar, near the y-axis, (x, y) = (−0.187, −8.681).

The Case A 4
In this case study, it was assumed that the 48 vortices with Γ 1−48 = 1 were distributed initially on the circle (green dots in Figure 15a with radius R = 3 and Cartesian coordinates defined as x j = R cos((j − 1)∆θ), y j = R sin((j − 1)∆θ); and ∆θ = (2π/(n − 2), j = 1, 2 . . . n − 2). Two vortices with strong intensities Γ 49−50 = 2(24 + 7 (6)) marked by navy-blue dots were located at (x 49 , y 49 ) = (1, −0.5) and (x 50 , y 50 ) = (3.4, −2), respectively. The vortex in location (x 50 , y 50 ) = (3.4, −2) was fixed (navy-blue dot with dashed circle). Figure 15a presents both the initial locations of vortices (green dots and navy-blue dots as fixed vortices, and a light-blue dot as another strong vortex) and their collapse locations (red dots and two navy-blue). Figure 15b displays the selected collapsing trajectories.    Figure 17 corresponds to the first frame of Figure 16. Instantaneous locations of a set of 48 weaker vortices (red dots), two strong vortices (two navy-blue dots) and the passive green and blue tracers located around two strong vortices can be observed. Time is provided at the top of each frame. The collapse time and Hamiltonian value are T c ≈ 20.85 and H ≈ 14.74. It can be noted that two sets of traces are separated by the vortex sheet, which surrounds the strong vortices. The vortex sheet takes a tear-shape structure. It also can be noted that two strong vortices move as a vortex pair. Over time, the passive markers are pressed out from the loop of the vortices because of fluid motion being incompressible and the vortex system's shrinkage to the singularity.  Figure 18 presents how two patches of passive tracers and the system of collapsing vortices presented in the last frame in Figure 16, behave over time. Similarly, in Figure 16, instantaneous locations of a set of 48 weaker vortices, two strong vortices and the passive green and blue tracers located around two strong vortices are presented. It can be observed that weaker vortices are organized into the tail, which as time pass by, behaves as a solid wall separating two strong vortices. The structure created by the vortices presented in Figure 18 resembles the structure of a tadpole galaxy. The tracers pressed out from the loop of vortices flowed up along the vortex tail. The tail acts as a solid arm of the large mixer.

The Case A 5
In this case study, four strong vortices were put in the system. The initial locations of 46 vortices with weaker intensity Γ 1−46 = 1 were placed at the peripheral of ellipse with axes a = 12, b = 2 and the cartesian coordinates of four strong vortices were (x 47 , y 47 ) = (−3, 5, 3.5), (x 48 , y 48 ) = (9.5, 2), (x 49 , y 49 ) = (1, −3.5) and the last, fixed vortex was located at (x 50 , y 50 ) = (−6, 0.5) (see Figure 19a). Figure 19b presents the calculated collapsing locations of vortices using the initial location of vortices presented in Figure 19a. As can be observed in Figure 19b, all vortices lay inside of a circle marked by dashed gray lines. Inside that circle, the strong vortices took positions at distinct cells. The boundaries of cells are created by the weak vortices that arrange themselves loosely along the arc of an ellipse and straight lines which make the walls of cells. The vortices are not closely located, but they are still able to keep passive tracers together and transport them together with strong vortices over long time. Figure 20 presents how the system of collapsing vortices presented in Figure 19b behaves over time. Passive color tracers were put around strong vortices (navy-blue dots). The tracers rotate quickly around vortices and they stay in each individual cell for a relatively long time. In a collapsing process, the distances between the vortices become smaller, the boundaries of the cells become shorter and the tracers are pressed out through cracks in cell walls. The time is given at the top of each frame.   The red dots correspond to the left end of the interval H ≈ −77.34, and the red dots to the right end, H ≈ −76.86]. The paths of the vortices from red to green locations are marked by gray dotted lines. It can be noticed that the Hamiltonian interval is relatively narrow, and the changes of the vortex locations are small. Figure 21b presents the dependence between collapse time and the value of Hamiltonian. It can be observed that as the Hamiltonian value increases, the collapse time T c increases from 148.35 to 37, 251.45. For Hamiltonian values near the right edge of the Hamiltonian interval, the collapse time sharply increases, which suggests that the system of vortices (green points in Figure 21a) is close to the relative equilibrium.

Conclusions
In this paper, a numerical approach to find the collapsing system of vortices was provided. First the algorithm was provided to find the collapsing locations from which vortices collapse. Then through series of selected, reproducible case studies very interesting properties and behavior of the system of collapsing vortices with at least one strong vortex were highlighted and examined. It was demonstrated numerically that the system of collapsing vortices in a self-similar motion can organize itself in vortex structures that are able to transport passive markers caught inside of that structures. The collapsing system of vortices should fulfill condition for the integral of motion V (viral), i.e., V = 0 (Equation (7)). The viral V is expressed as double sum mutual products of vortex intensities. The condition V = 0 can be realized in many ways. If we assume that the majority of the vortices in the collapsing system have the intensity equal to one then a few vortices with negative intensities will have intensities much greater than one. Presented numerical results confirmed that such a situation always results in a tight gathering of vortices with weaker intensities along some curves. These densely located vortices along some curves can be regarded as the vortex sheets. The vortex sheets behave as a solid barrier, that are impermeable to tracers from inside out and vice versa. The self-organization of point vortices into vortex structures that are capable of transporting passive tracers is a novel element of point vortex dynamics. This can be regarded as the prototype of the coherent vortex [27] and may provide a new mechanism of transport in geophysical processes like the dispersion of pollutant in the atmosphere and in the ocean by advection [7].
In this article, only selected examples emphasizing the properties and behavior of the system of collapsing vortices were presented so that the reader can reproduce the results independently. Examples were selected out of dozens of collapse system vortices obtained by author. Based on author experience, it is possible to find solutions describing the location from which vortices collapse for huge variety of vortex systems with different initial conditions.
The author chose certain examples because they seemed interesting in form, but at no time were forms of vortex structures exhausted. It seems that the number of different collapse structures of vortices is almost infinite. The change in the initial starting locations nearly always leads to a different collapse configuration. Further investigation using field experiments or observational studies of atmospheric events would be helpful to reinforce the discovery of the collapse phenomenon with the presence of a vortex sheet. From the dynamical systems theory point of view, the analytical proof that the collapsing system of vortices can self-organized in vortex sheets would be desirable. The numerical results here are a great stepping stone to further understand the collapse phenomena of many vortices and their ability to self-organize and transport the passive tracers.

Acknowledgments:
The author would like to thank Grzegorz Karch from the Institute of Mathematics of Wroclaw University for his helpful comments on the first draft of this paper. The author would also like to thank Maria Kudela for a thorough review of this paper.

Conflicts of Interest:
The author declares no conflict of interest.

Appendix A. Derivation of Helmholtz Equations
The equations of incompressible inviscid fluid motion (Euler equations) are ∂v ∂t where v = (u, v, w) is velocity vector, ρ is fluid density, p is pressure and f is mass force, for which we assume that there is potential Φ, f = ∇Φ. By using the vector to identify (v · ∇)v = ∇ 1 2 v 2 + ω × v, Equation (A1) can be rewritten as ∂v ∂t where B = v 2 2 + Φ + p ρ . Taking rotation from both sides of Equation (A3) results in: Clearly, the motion of the fluid describing by the Euler Equation (A1) with the potential mass field f is possible only when the above equation is fulfilled. Further, we will regard only a flow in 2D space. Hence, the velocity has only two components v = (u, v, 0), while the vorticity takes the form ω = (0, 0, ω). The transport vorticity equation takes the form: The transport vorticity Equation (A5) results that along the particle paths, vorticity ω is constant. The trajectory of α-particle x(t, α) is defined by the flow map Φ : Ω 0 → Ω t (see Figure A1). Figure A1. Trajectory of the α particle.
The motion of particles along the path is described by differential equations: Incompressibility ∇ · v = 0 guarantees the existence of the vector potential A, such that v = rotA. In 2D space, this vector potential has only one component: A = (0, 0, ψ), where ψ is a known stream function. The stream function ψ and velocity can be computed by solving the Poisson equations: The solution of the Poisson equation in the entire 2D space can be expressed as ψ(x, y) = − 1 2π D ln r ω(x , y )dx dy , r = (x − x ) 2 + (y − y ) 2 (A8) and velocity as v(x) = (u(x, y), v(x, y)) = R 2 ((y − y ), −(x − x )) 2πr ω(x )dx dy (A9) Now we approximate the vorticity by delta Dirac measure distribution: where the value of Γ k is referred to as the intensity of the point vortex which is equal to the circulation around the point k. In practice one can regard Γ k as a mean value of vorticity over a small area A k = h 2 .
By substituting (A10) into the Equation (A9), we obtain the system of differential equations: , y(0) = α 2 , j = 1, . . . n (A13) where r jk is the distance between the vortices i and j and the prime on the summation indicates omission of the term with j = k. The system (A12) for two-point vortices was used the first time by Helmholtz (1858) [1].