Equations of motion for weakly compressible point vortices

Equations of motion for compressible point vortices in the plane are obtained in the limit of small Mach number, M, using a Rayleigh–Jansen expansion and the method of Matched Asymptotic Expansions. The solution in the region between vortices is matched to solutions around each vortex core. The motion of the vortices is modified over long time scales O(M2log⁡M) and O(M2). Examples are given for co-rotating and co-propagating vortex pairs. The former show a correction to the rotation rate and, in general, to the centre and radius of rotation, while the latter recover the known result that the steady propagation velocity is unchanged. For unsteady configurations, the vortex solution matches to a far field in which acoustic waves are radiated. This article is part of the theme issue ‘Mathematical problems in physical fluid dynamics (part 2)’.


Introduction
Vorticity is a key aspect of fluid mechanics, in particular in high-Reynolds-number and turbulent flows. being governed by a set of ODEs. A discussion of the justification for the equations of motion of point vortices and generalizations is given in [1], including a review of the momentum argument set out by [2].
The majority of the work to date on point vortices has been for plane incompressible flows [2]. Attempting to extend the notion of a point vortex to plane compressible flows is a daunting task in the general case, but for the case of low-Mach number flows, a Rayleigh-Jansen expansion in Mach number provides one approach. With the Mach number used in the expansion defined by the velocities induced by the vortices' motion and the speed of sound, the O(1) incompressible velocity field increases so as to become supersonic near the location of a point vortex. One needs to consider further physics near the vortex location, i.e. the vortices have a small core region. Barsony-Nagy, Er-El & Yungster (hereafter BNEEY) [3] showed how to obtain steady point vortex configurations in this manner, relating the core behaviour to a solution obtained by Taylor [4].
Since then, there have been a few similar studies. These have examined the translating vortex pair [5][6][7], for which it was found in [7] (hereafter L06) that the speed of propagation was unchanged at O(M 2 ), and the von Kármán vortex street [8], for which the speed of propagation for both staggered and unstaggered streets can either increase or decrease depending on parameters of the flow. (There have also been works on steady weakly compressible hollow vortices, as in [9][10][11], but these do not consider point vortices.) As pointed out by [5], the existence of a family of continuous shock-free transonic compressible flows with embedded vortices is of intrinsic interest, given that similar flows for transonic aerofoils do not persist under small perturbations. Our goal is to extend the work on weakly compressible point vortices to the unsteady case. We extend the approach of BNEEY to obtain equations of motion for the positions of the vortices up to O(M 2 ). Our approach is based on conservation of momentum, which has been used for incompressible constant-density flows and which we now review (see [1]). We compute the rate of change of momentum inside a moving closed contour C from Newton's Second Law in complex notation, where the contour C is described in the positive sense, the complex momentum inside C is given by the area integral P = S ρw dS, the complex position and velocity are z = x + iy and w = u − iv respectively, and the velocity of C is given by W = U − iV. Using Bernoulli's equation, substituting local expansions for the variables into (1.1), and taking the limit as the contour shrinks down to the vortex, gives ζ t = W =w, since the contour moves with the vortex. Here,w is the desingularized velocity at the vortex: physically, a point vortex moves with the local desingularized flow. In this incompressible argument, point vortices have no internal structures, but for the compressible case, we will need to consider the flow in the vortex cores on scales of O(M) smaller than the distance between vortices. In addition, at large distances of O(M −1 ) from the vortical region, as pointed out by L06, the Rayleigh-Jansen expansion will become disordered, indicating the presence of a wavelike far field. This feature was already present in previous work on sound generation by vortical flows in aeroacoustics using matched asymptotic expansions (MAE; see [12,13] for an overview).
Terms in the Rayleigh-Jansen expansion will evolve on slow time scales to allow the position of vortices to change. Figure 1 illustrates the different regions of the flow and gives some notation. In §2, we present the governing equations for the vortical region, and then examine the core region to understand the nature of the expansion. In §3, we obtain the equations of motion for a vortex to O(M 2 ) by solving at successive orders in the Rayleigh-Jansen expansion, and then discuss the form of the global solution. We examine the case of two compressible vortices in §4. We discuss the far field in §5 and relate it to previous work. Finally, §6 concludes the paper. Electronic supplementary material includes details of algebra, as well as an account of the formal matching near the vortex cores.

Problem formulation
We consider irrotational adiabatic compressible flow in the plane. The adiabatic relation between pressure, p * , and density, ρ * , takes the form p * /p * 0 = (ρ * /ρ * 0 ) γ , where γ is the constant ratio of specific heats, and p * 0 and ρ * 0 are reference values for pressure and density, taken to be the values at infinity where the flow is at rest. The momentum equation can be transformed into the unsteady Bernoulli equation, where φ * is the velocity potential and the speed of sound (squared) is given by c 2 = dp * /dρ * | s * = γ p * /ρ * with constant value c 2 0 = γ p * 0 /ρ * 0 at infinity. It is convenient to combine the above equation and the continuity equation into a single equation for the velocity potential, the Blokhintsev equation [7,14].
We non-dimensionalize using a length characteristic of the distance between vortices L, a typical velocity V induced by one vortex on another, the resulting time scale L/V, as well as the value of density at large distances, ρ 0 and the dynamic pressure scale ρ 0 V 2 . Then the Mach number is M = V/c 0 , and the Blokhintsev equation becomes, dropping the stars and using the summation convention with subscripts running from 1 to 2, . (2.4) The above equations are valid in the region of length scale L between vortices, which we call the vortex region. They break down near the vortex cores, as pointed out by BNEEY and also examined by L06. To understand the flow behaviour in a vortex core and its impact on the subsequent matching process, we work in a reference frame co-moving with the vortex, so that for a vortex at location X moving with velocity U, one has x =x + X, u =û + U, t =t and φ =Φ + U ·x, (2.5) wherex andû are the position and velocity with respect to the vortex core in the moving frame, respectively. Here,t is used to emphasize that partial time-and space-derivatives in the core frame are taken with constantx andt, respectively. The velocity potential in the core frame isΦ. Following previous authors, we now define an appropriately scaled variable in the core region usingx = Ms, so that the radial coordinate measured from the vortex core is r = Ms. In terms of these variables, the Blokhintsev equation becomes where dots above U j indicate time derivatives. The Rayleigh-Jansen expansion is an expansion in small Mach number. Before writing it down in §3, we consider the momentum equation in complex form, as in [1], taken over the small circle with radius e. We take e M, as we are interested in the momentum balance over circles that are asymptotically small with respect to the region between vortices but much larger than the vortex cores, i.e. e is an intermediate variable in the terminology of MAE. This means that one can use either the core solution or the vortex solution when evaluating the right-hand side of (1.1), since the terms on the right-hand side are all contour integrals evaluated at radius e. The left-hand side, however, is a surface integral that must be calculated in the inner variable. To leading order in M, we haveΦ = κθ where κ = Γ /(2π ) is the scaled circulation and θ is the polar angle, so that in vectorial form, the momentum is where t is the unit vector tangential to the circle. The lower limit s v is the smallest value of s for which the pressure and density in (2.3) and (2.4) are positive, and is obtained from the condition 2s v 2 = (γ − 1)κ 2 . We see that the term in t cancels by symmetry. In complex form, we then obtain where ψ(z) is the digamma function and we write β = (γ − 1) −1 for brevity. We see that (

Derivation of the equation of motion for a vortex (a) Rayleigh-Jansen expansion of the global solution in the vortical region and time dependence
Motivated by the discussion above, we consider a modified version of the Rayleigh-Jansen expansion in the form where the arguments of the two first terms reflect the fact that they correspond to incompressible flow in the vortical region, since from (2.2) we have ∇ 2 φ 0 = ∇ 2 φ 1 = 0. The governing equation for φ 2 , does not contain γ . Define a complex potential F 0 = φ 0 + iψ 0 , since the flow at leading order is incompressible and irrotational. Similarly, there is a complex potential F 1 = φ 1 + iψ 1 . These potentials are harmonic functions that decay far from the vortex. Since F 0 has logarithmic singularities, F 1 cannot have singularities of higher order, while logarithmic singularities in F 1 are disallowed by requiring the vorticity to be entirely at O (1). Hence, as an analytic function bounded at infinity with no singularities, F 1 is a constant that can be taken to be 0 without loss of generality. This means that the Rayleigh-Jansen expansion does not in fact have a term at O(M 2 log M). The term entering the matching from (2.8) must therefore come from taking the O(M 2 log M) time dependence into account appropriately, as is done in (3.6) below.
Using φ 0 = (F 0 + F 0 )/2, (3.2) becomes where c.c. stands for complex conjugate and we have defined a function F 2 (t, z, z) such that φ 2 = (F 2 + F 2 )/2. Note that, because the flow at O(M 2 ) is no longer incompressible, there is no streamfunction corresponding to φ 2 , so we call F 2 a potential but not a complex potential. Only the real part of F 2 matters, and the complex velocity w 2 = u 2 − iv 2 is given by w 2 = ∂ z (F 2 + F 2 ). We can integrate (3.3) and obtain a particular solution for φ 2 as the real part of where Z 0 is a time-dependent centre of vorticity that can be picked to simplify the analysis for specific cases. The functions I(z) and J(z) are defined globally by The integration limits z I and z J will also be picked depending on the global nature of the flow. It turns out that the location of the vortex, considered as an expansion in Mach number, is not uniquely defined. We can remove this ambiguity by requiring that the location not be expanded in Mach number (this is reminiscent of slaving principles as in [15]). However, we need to allow the position to evolve in time at higher order in M to allow matching of the terms in (2.8). This leads to the use of multiple time scales. Given the form of (2.8), we consider all variables to be functions of t 0 = t, t 1 = tM 2 log M and t 2 = tM 2 (and possibly further time scales). We define where we write ζ ,i = ∂ζ /∂t i . We can now write pressure and density using the full time dependence by expanding (2.3) and (2.4) written in complex notation and There is a dynamically irrelevant component p −2 = (γ M 2 ) −1 that can be ignored, while ρ 0 = 1.

(b) Local solution for the O(M 2 ) potential
We now consider the solution near the vortex located at z = ζ . The following expansions provide the terms needed to compute the equations of motion, writing = z − ζ and e = | |: 14) and We have expressed the coefficients of F 0 in the above form for later convenience. The coefficients f 0 , f 1 , . . . depend on time in general. While the coefficient f 0 appears to have no dynamical meaning, it is different from one vortex to another and hence is kept. The first step is to ensure that the velocity field obtained from F 2 in (3.4) is single-valued. While F 0 has a multi-valued logarithmic term near ζ , the resulting complex velocity, w 0 , is singlevalued. However, F 2 contains logarithmic terms of the form log multiplied by functions of z. The resulting velocity is not single-valued, but becomes single-valued when adding a homogeneous solution with the same z-dependence multiplied by log . This corresponds to including the following contribution in G(z): The square bracket defines the function l(z). The effect of (3.17) in the expansion of F 2 near a vortex is to replace log terms by log e 2 . We then include two homogeneous terms in G(z) The coefficients μ and ξ allow us to remove unacceptable singularities in F 2 near = 0. Since F 0 has logarithmic singularities, F 2 cannot have singularities of higher order, or else the expansion would be disordered near = 0. Logarithmic singularities in F 2 are removed by requiring the vorticity to be entirely at O (1). Unlike [8] we do not require a term in −2 . Finally, the local expansion of G(z) also contains terms from expanding the counterparts of the terms inside the square bracket in (3.4) due to other vortices and to other homogeneous contributions to G(z). From the forms of F 0,0 and of I(z), the former lead to terms of O(e −1 ) as well as analytic terms. The latter are denoted by K(z). The result is a contribution g −1 −1 + g 0 + We can now write down the expansion for F 2 needed to remove singular terms: The term in −2 is purely imaginary and hence can be ignored, while the term in −1 log e 2 cancels from the O(1) result, which will be rederived in the current framework and notation below. Since F 2 only enters the solution via its real part, removing the singular terms in its real part leads to the following conditions: and

(c) Conservation of momentum
We now return to the conservation of momentum, viewing it as a matching problem for an expansion in M, with e serving as the independent variable. We define Q as the left-hand side of (1.1) and express it as an expansion in the inner variable, s, and define q as the right-hand side and express it as an expansion in the vortex variable, e. We expand the time-derivative of (2. variable. Then Q 0 = Q (0,.) = 0 since there are no O(1) terms in (3.22). Expanding the right-hand side of (1.1) leads to the exact result (3.23) The (0, 0) term in Van Dyke's rule is Q (0,0) = q (0,0) , so that Using ζ ,0 = W 0 leads to i.e. the incompressible result expressed in the current notation.
We should now group the O(M 2 log M) and O(M 2 ) terms together to continue with Van Dyke's rule. We should also compute Q (2,.) . We shall avoid doing this, and instead carry out the matching informally. This approach works, but to be safe we will revisit the formal matching in the electronic supplementary material. In the vortex region, the right-hand side of (1.1) gives the O(M 2 log M) contribution To obtain this, we need to consider further terms in the local expansion of F 2 . Using the conditions (3.20) and (3.21), we find where the relation means that the equality ignores purely imaginary terms. The coefficients needed are and Substituting into (3.29) and computing the integrals leads to extensive cancellation, yielding We see from (3.22) that the log e term in q 2 cancels the log e term at O(M 2 ) in Q. Recalling that ζ ,2 = W 2 gives the equation for the slow evolution of ζ as It is useful to check the behaviour of a single point vortex. The incompressible complex potential is F 0 = −iκ log (z − ζ ). The point vortex does not move at O(1). Hence ζ ,1 = 0 and The arbitrariness of z I is irrelevant, as it is cancelled by μ 1 when removing the simple pole in F 2 .
The leading-order term is purely imaginary so it can be ignored. Hence the O(M 2 ) velocity of a single point vanishes, a necessary feature for this model.

(d) Global solution
The results above are applicable near every vortex, because neither the vortex circulation nor the vortex location has a preferred value. We can now assemble a global solution that is valid everywhere in the vortical region. The O(M 2 ) potential is given by the sum of the inhomogeneous part (3.4) and of a homogeneous part. The homogeneous part takes the form (3.38) and includes contributions from each vortex of the form (3.17) and (3.18), while K(z) includes possible further terms (e.g. to satisfy boundary conditions or to set the circulation around objects in the flow). We now consider the expansion near vortex n of the sum in (3.38) omitting term n, writing the rest of the sum as g −1 −1 n + g 0 + g 1 n + O(e 2 n ), with n = z − ζ n and e n = | n |. The calculations above show that g 0 is not needed. We find, for vortex n, and where ζ mn = ζ n − ζ m and the prime in the summation indicates that term n is omitted (the derivation is given in the electronic supplementary material).

Two vortices in the plane
We consider the simplest situation consisting of two point vortices in the infinite plane. In cases such as this, there are no other contributions to the potential F 0 beyond the vortices, which means that K(z) = 0 in (3.38). While some simple geometries can also be solved using the method of images and could hence be considered as consisting of a finite number of vortices, the dynamics of the actual and image vortices are different: the latter are not physical so that their motion is set by boundary conditions rather than matching.
Results for the leading-order potential, complex velocity, I and J are equally simple for N vortices; we will take N = 2 after presenting general results. We have The decay properties of w(z) for large |z| mean that we can take z I = ∞. It is known that the incompressible two-and three-vortex cases are integrable. The system has four real conserved quantities, two of which combine to give the complex vortex momentum N n=1 Γ n ζ n . Conservation of this quantity means that integral J(z) is convergent at infinity, so that we can take z J = 0. To calculate the integrals I and J, we use the primitives where the primed sums indicate m = n.
For the two-vortex case, the total circulation is Γ ∞ = Γ 1 + Γ 2 , so there are two different cases, corresponding to Γ ∞ = 0 and Γ ∞ = 0. In the latter case, the conservation laws show that the vortices must stay in a bounded area of the plane. The former case corresponds to a copropagating dipole pair in the incompressible limit. Expressions for f 0 , f 1 , J 0 and g −1 for the two vortices are given in the electronic supplementary material. The relation (3.21) gives ξ 1 = ξ 2 = 0 if Z 0 is taken to be on the line joining ζ 1 and ζ 2 , although the final result for the motion of the vortices is independent of Z 0 .
Calculations (see the electronic supplementary material) lead to H 1 = 0, with the Z 0 terms cancelling. Substituting into (3.36), along with f 1,0 = 0, means that ζ does not change with t 2 .
Hence we recover the result of L06 and [8]: the translation speed of the co-propagating vortex pair does not change at O(M 2 ). The current procedure is of course lengthier than that needed to obtain this result in a co-moving frame in which the pair is at rest, but we can now address the fundamentally unsteady co-rotating pair.
For the co-rotating case, κ = κ 1 = κ 2 , so that ζ 1 + ζ 2 is independent of t 0 . Without loss of generality, we take the positions of the vortices at t = 0 to be ±a 0 with a 0 real, yielding ζ 1 = ae iϕ and ζ 2 = −ae iϕ (4.5) with ϕ = θ + κt/(2a 2 ) = θ + ωt. Here, a and θ are functions of t 1 and t 2 with a(0, 0) = a 0 and θ(0, 0) = 0. Since ζ 1,00 = −ω 2 ζ 1 , we find from (3.28) that a ,1 = 0 and θ ,1 = ω 1 = κω 2 The radius is only a function of t 2 , while the rotation rate varies with t 1 . Since the θ ,1 term has the same sign as the O (1)  and Here, the constant C is given in (2.10). Once again there is a correction to the rotation speed In dimensional form, we can combine the frequencies to obtain where stars represent dimensional quantities. The circulation scale is Γ = LV. There is no unique choice of scalings, but the simplest choice is probably Γ = 2κ * and L = 2a * , so that Γ is the total circulation divided by 2π and L is the distance between vortices. Then For the general two-vortex case, the vortices rotate about their centre of vorticity at O(1). Write ζ 1 = ae iϕ and ζ 2 = −be iϕ with aκ 1 = bκ 2 so that the centre of vorticity is at the origin with Z 0 = 0. Then . (4.12) The O(M 2 log M) equations become The velocities of the vortices are the same, but their angular velocities differ, so that as they move their trajectories will no longer be circles. The centre of the vorticity ζ c = (κ 1 ζ 1 + κ 2 ζ 2 )/(κ 1 + κ 2 ) moves slowly with the O(M 2 log M) velocity (4.14) The O(M 2 ) motion is  This shows that the relative equilibrium of the co-rotating vortices evolves slowly in time due to the weak compressibility effects, although in the symmetric case, it is only the rotation rate that changes.

Outer region
As pointed out by L06, at distances of O(M −1 ) from the vortices, there is a region in which the dynamics are wave-like. We define far-field upper-case variables by X = Mx. The Blokhintsev equation becomes One can obtain expressions for the pressure and density analogous to (2.3) and (2.4). The far-field limit of the O(1) potential in the vortical region is where the total circulation is 2πκ ∞ and the g i (t) are functions of time alone. The incompressible dynamics of the vortex region imply that κ ∞ , and g 1 , which is related to the vortex impulse I [2], are independent of time. We set g 0 = 0 as the global constant of integration for F 0 . Rewritten in the outer variable, we have From (5.1), Φ 0 is a harmonic function, which must be κ ∞ θ, so that the real part matches. Since the higher terms in (5.1) are multiplied by powers of M 2 , Φ 1 is again a potential function, and must be Re(g 2 Z −2 ). This result requires checking the behaviour of F 2 for large |z|. The general question of compressible far-field waves relates to the domain of aeroacoustics. Excellent reviews of aeroacoustics are given e.g. in [12,17]. The use of MAE in solving aeroacoustic problems is reviewed in [13]. For the general case, we can use these results. We hence summarize the necessary results following [12]: wavelike solutions with dipole R −2 cos (2θ − χ ) behaviour near the origin are given by a synthesis of monochromatic solutions of the form 2 is the Hankel function of the second kind. These solutions are superposed, with amplitude coefficients A(ω) that can be related to g 2 . One obtains quadrupole radiation of flow-generated sound, as expected in a situation with no boundaries or mass sources.
For the two cases considered in §3, we can follow previous authors directly. Acoustic emission from the co-rotating vortex pair was first examined by [18], who found a change of rotation frequency equivalent to that in §3. The velocity of the co-propagating pair is found to be independent of t 2 , t 1 and t 0 , so the analysis of L06, who considered a steady dipole in a moving reference frame, is equivalent to the current one. L06 shows by using a coordinate system moving with the vortex that the O(M) solution, a propagating dipole, is efficiently represented in the coordinate system ((1 − M 2 ) 1/2 X, Y). This is equivalent to a different way of writing an asymptotic solution that gives the same result to the order obtained.

Conclusion
We have obtained equations of motion for small-Mach-number compressible point vortices in the plane, in which compressibility manifests itself as an evolution over slow time scales of O(M 2 log M) and O(M 2 ). The first correction (3.28) is quite simple and vanishes for steady configurations. The second correction is more involved.
We have examined the corrections to O(M 2 ) for the simplest case of vortex pairs. We recover the known result for the co-propagating pair that the velocity is unchanged at O(M 2 ). The symmetric co-rotating vortex pair exhibits a change of angular velocity. For the general two-vortex case, however, the centre of rotation and radius of the orbit evolve slowly, while the motion of each vortex is instantaneously perpendicular to the line of centres and the motion remains circular on the O(1) time scale.
The solution in the far-field region with spatial scale M −1 L, corresponding to the wavelength of the emitted sound, can be obtained by matching, following previous work. If the vortical flow is steady, the response is a dipole moving with the speed of the centre of vorticity, as in L06. If the vortical flow is unsteady, an expression for the quadrupole radiation is obtained in terms of the quadrupole moment g 2 (t) (presumably this could be applied to the calculation of wave radiation by chaotic point vortex evolution as in [19,20]). Following previous work, quantities such as the power radiated to infinity could be obtained.
The back-reaction of the wave field is deliberately ignored here, as is usual in aeroacoustics. This means that while radiation is present in the current formulation, the coupling of the flow in the vortical region to the far field only appears at an order higher than M 2 . In the geophysical context, the corresponding effect of gravity wave radiation on vortex dynamics has been examined [21] (see also [22] for a related discussion for scattering). Here, this would require a calculation to O(M 4 ), most likely with logarithmic terms.
A list of interesting extensions comes to mind: efforts at simplifying the equations further in special cases; possible efficient solution techniques; whether any of the other known equilibria of point vortices survive to O(M 2 ); the effects of more complicated boundaries and whether a better model for the core regions is warranted. The effect of boundaries is of particular interest, following on from BNEEY. The case of vortices inside and outside a circle is currently being examined: the O(1) solution can be obtained using the method of images, but the corrections require extensive algebra. A further example is the half-plane with a vortex moving around it considered in [23]). Finally, it could be interesting to compare the present results to numerical simulations. These would require very highly resolved aeroacoustic-type calculations with large separations of scales between vortex cores, vortical region and wave field. 1 Data accessibility. The data are provided in electronic supplementary material [24].