Evidence of singularities for a family of contour dynamics equations
 ^{*}Instituto de Matemáticas y Física Fundamental, Consejo Superior de Investigaciones Científicas, C/Serrano 123, Madrid 28006, Spain; ^{‡}Departamento de Matemáticas, Universitad Autónoma de Madrid, 28049 Madrid, Spain; and ^{§}Mathematics Department, Yale University, P.O. Box 208283, New Haven, CT 06520208283
See allHide authors and affiliations

Communicated by Charles L. Fefferman, Princeton University, Princeton, NJ, March 10, 2005 (received for review January 19, 2005)
Abstract
In this work, we show evidence of the existence of singularities developing in finite time for a class of contour dynamics equations depending on a parameter 0 < α ≤ 1. The limiting case α → 0 corresponds to 2D Euler equations, and α = 1 corresponds to the surface quasigeostrophic equation. The singularity is pointlike, and it is approached in a selfsimilar manner.
One of the most important open problems in mathematical fluid dynamics is whether the solutions to the Euler and NavierStokes equations modeling the evolution of incompressible inviscid and viscous fluids, respectively, may develop singularities in finite time. Several possible scenarios for a singularity have been proposed in the past (see ref. 1 for an account of some of them), although none of them led to a rigorous proof of the formation of singularities. One of these scenarios corresponds to the socalled vortex patch problem that we briefly describe below.
A vortex patch consist of a 2D region D(t) (simply connected and bounded) that evolves with a velocity given at each instant of time by [1]
where the stream function ψ is such that ω = Δψ, and the vorticity ω has a constant value ω_{0} over D(t). Vortex patches are, therefore, weak solutions of Euler equations. The appearance of finite time singularities at the contour of D(t) was the subject of strong debate based on numerical results showing that the curvature of the boundary might grow superexponentially in time (see refs. 24). Nevertheless, work of Chemin (5) and Bertozzi and Constantin (6) rigorously proved global existence of regular solutions for the dynamics of the vortex patch, and, at least in this case, singularities cannot appear.
A very natural singularity scenario would correspond to a vortex patch of the socalled surface quasigeostrophic equation (see ref. 7) for which the relation between the stream function and potential temperature θ (that play the role that vorticity plays in 2D Euler equations) is θ = (Δ)^{1/2}ψ. The interest of this equation lies in its strong analogies to the 3D Euler equation as it has been argued in ref. 8 and its physical relevance as a model for the formation of temperature fronts in some geophysical contexts (see refs. 9 and 10).
Hence, a natural question to ask is whether models for which θ = (Δ)^{1}(α/2)ψ representing an interpolation between 2D Euler and quasigeostrophic patches develop singularities for 0 < α ≤ 1. In this work, we provide numerical evidence showing that this scenario is indeed the case and describe the selfsimilar structure of such singularities. This result represents a previously undescribed singularity scenario for incompressible flows.
The Model
Following Zabusky et al. (4), we can invert the relation between ψ and θ to obtain the following formula for the velocity with 0 < α < 1 [2]
where is the position of points over C(t), the boundary of D(t), parameterized with γ, and . Here is the value of θ in the patch and the factor c_{α} =Γ(α/2)/2^{1}αΓ(2α/2) results from inverting (Δ)^{1}(α/2). In the limiting case α = 1, one should use the following formula for the velocity (see ref. 7) [3]
The contour dynamics equation is then
Given the fact that only the normal component of the velocity is able to deform the curve, the contour dynamics equation also can be written as
In the periodic setting, where the front is described in the form , with θ = 1 above the curve, and 0 below the curve, one arrives then at the following evolution equation: [4]
If ϕ is smooth enough, then the integrand in Eq. 4 is not singular. This formulation is used below to prove local existence. The generalization of these equations to the case of nonsimply connected domains D(t) is straightforward.
Local Existence. Here, we sketch the proof of local existence and uniqueness of a solution. For simplicity, we will restrict to the periodic case, but it is easy to see that the argument presented here extends to the case of any sufficiently regular initial contour or finite set of disconnected initial contours.
In the periodic setting, we consider the scalar function θ(x, y, t) as given by [5]
To obtain an equation involving only the curve describing the front, we start by inverting (▵)^{1}(α/2)ψ = θ and then use the relationship u = ((∂ψ/∂y), ∂ψ/∂x) between the velocity and the stream function to eliminate ψ from the system.
The expression obtained for the velocity becomes singular as we approach the front, but the singularity is in the direction of the tangent to the front. We can use the fact that θ is convected by the fluid (i.e., ) to observe that u can be modified by adding a term in the direction of because θ always satisfies . By using the adequate function h, it is possible to obtain a less singular expression for the equation, we obtain For the above equation, it is possible to set up a NashMoser argument (11) following the ideas introduced in ref. 7, which thus provides as a result the local existence and uniqueness for the front.
The Numerical Method. We have numerically solved Eq. 2 for 0 < α < 1 for the case of two vortex patches, [6]
The initial state consists of two identical circular vortex patches of unit radius whose centers are separated by 2.5 vortex radii and such that θ_{k} = 1. The evolution is calculated by using ideas developed in refs. 3 and 12 for the case of the 2D Euler equations. In particular, as described in ref. 12, contours are represented by an adaptive node spacing adjusted nonlocally at each time step and cubic spline interpolation between nodes. Contour integrals along each cubic contour segment are evaluated as follows: series expansions on small parameters are used to evaluate both singular and regular integrands. However, in the latter case, convergence is not always satisfied and then the contour integral is numerically evaluated with a fifthorder adaptive RungeKutta method. After node representation of contours, Eq. 6 is transformed in a set of coupled ordinary differential equations that are forward timeintegrated with a variable timestep fourthorder RungeKutta method. Timestep size Δt is adjusted with Δx, which is the minimum distance between nodes in contours C_{1}(t) and C_{2}(t) at a time t. In particular we choose Δt = AΔx, where A is tuned for different values of α. For the case α = 1, we simulate the two vortex patches version of Eq. 3, [7] where singularities in the integrands are removed; however, the numerical ideas presented above are still valid. In fact, Eq. 7 is also valid for 0 < α < 1.
Evidence of Singularities
We performed several numerical experiments for different values of α. In Fig. 1, we present the profiles at t = 0, 7.066, 11.046, 13.103, 15.383, and 16.515 for α = 0.5, and in Fig. 2, we present the same at t = 0, 1.962, 3.3, 4.031, 4.369, and 4.464 for α = 1 (the surface quasigeostrophic equation). As we can see, two corners with high values of curvature do develop in the last plots. It looks like a sharp front appears, in the sense that the boundary of the two patches seems to collapse along a curve. In fact, this is not the case, as one can see from Fig. 3, where the two curves appear clearly separated except for one point at which they get so close that they seem to touch each other (the corner). In Fig. 4, a magnification of the corner region at several times is represented, together with a representation of the maximum curvature of the patch contour as a function of time.
The fast growth of the curvature suggests that it might blow up at a finite time t_{0}. In fact, if we rescale the spatial variable in the form and introduce the new time variable τ = log(t_{0}  t), Eq. 2 transforms into [8]
provided that δ = 1/α. Convergence to τindependent solutions to Eq. 8 would represent convergence to solutions of Eq. 2 such that the maximum curvature (inverse of the minimum radius of curvature R) grows as [9]
The minimum distance between the two contours would follow the law
To test these laws, we represent in Fig. 5, for the case α = 1, the values of κ^{1} and d with respect to the time for the whole simulation period. Along the evolution we can identify several regimes. At approximately t ≈ 3, a crossover between tendencies does take place, leading to a linear behavior for both κ^{1} and d as we can see in Fig. 5 Inset. It is clear from the figure the tendency toward zero of both. In fact, we stopped the simulation when κ^{1} reached values of the order 10^{4}, which is indistinguishable from zero in the figure. This abrupt change between different regimes is quite characteristic of problems involving singularities. Finally, in Fig. 6, we represent the numerical profiles depicted in Fig. 4a rescaled a rate (t_{0}  t)^{(1/α)} in both spatial directions. It is clear from the plot that they converge to a stationary profile that would be a τindependent solution of Eq. 8 and a selfsimilar solution of Eq. 2.
Finally, we have repeated the simulation for several values of α to verify our scaling laws. In Table 1, we summarize the results obtained for some values of α. In particular, we have performed a leastsquares fitting to the curvature law (Eq. 9) for a set of values of the curvature κ_{i} close to the collapse time t_{0}. This technique estimates t_{0}, the exponent 1/α, and the proportionality constant C in Eq. 9. Observe that t_{0} increases very rapidly as α tends to zero. In fact, it is known (5, 6) that t_{0} =∞ for α = 0, which suggests that t_{0} should be very large for small values of α, but eventually the patches also should develop a pointlike singularity. We have run experiments for 0 < α < 0.5, which confirm this.
The computational part of this work was performed on the hidra 40 processor (AMD Athlon XP at 2 Ghz) beowulf cluster at the Centro de Apoyo Tecnológico of the Universidad Rey Juan Carlos (Madrid) and on the superordenador virtual gallego (80 processor Pentium 4 at 3.2 Ghz) at Centro de Supercomputación de Galicia.
To summarize, we have presented in this work a previously undescribed kind of singularity and supported its existence in both scaling arguments and numerical evidence. A rigorous study of the stationary or quasistationary solutions of Eq. 8 would lead to an eventual proof of the existence of these singularities.
Acknowledgments
We thank A. Córdoba and J. Eggers for their useful comments and suggestions. D.C., M.A.F., and J.L.R. were partly supported by Ministerio de Ciencia y Tecnología Grant BFM200202042. A.M.M. was partly supported by Ministerio de Educación y Ciencia Grant MTM200400797 and a Ministerio de Ciencia y Tecnología (Spanish Government) Ramón y Cajal Research Fellowship.
Footnotes

↵† To whom correspondence should be addressed. Email: dcg{at}imaff.cfmac.csic.es.

Author contributions: D.C., M.A.F., A.M.M., and J.L.R. performed research and wrote the paper.
 Received January 19, 2005.
 Copyright © 2005, The National Academy of Sciences
References
 ↵
Majda, A. & Bertozzi, A. (2002) Vorticity and Incompressible Flow, Cambridge Texts in Applied Mathematics (Cambridge Univ. Press, Cambridge, U.K.), Vol. 27.
 ↵
 ↵
 ↵
 ↵
Chemin, J. Y. (1993) Ann. Sci. École Normale Supérieure 26, 517542.
 ↵
 ↵
 ↵
 ↵
Córdoba, D., Fefferman, C. & Rodrigo, J. L. (2004) Proc. Natl. Acad. Sci. USA 101, 26872691.pmid:14978276
 ↵
Pedlosky, J. (1987) Geophysical Fluid Dynamics (Springer, New York), pp. 345368 and 653670.
 ↵
 ↵