Random matrix model with external source and a constrained vector equilibrium problem

Introduction

The random matrix model with external source is the probability measure

The model (1.1) was first studied by Brézin and Hikami and P. Zinn-Justin who showed that the eigenvalue correlations are determinantal. In it was observed that the correlation kernel can be expressed in terms of multiple orthogonal polynomials. Due to the Riemann-Hilbert problem for multiple orthogonal polynomials this opened up a new way for asymptotic analysis. For the quadratic case

with two eigenvalues $\pm a$ of equal multiplicity (thus $n$ is even), this was done in great detail in the three papers .

The quadratic case is of special interest because it has an equivalent formulation in terms of non-intersecting Brownian motions that start at one value and end at certain prescribed values which is a variation on Dyson’s Brownian motion . The quadratic model with external source (1.3) exhibits a phase transition, since for small $a>0$ the eigenvalues accumulate on one interval while for larger $a$ the eigenvalues accumulate on two disjoint intervals. At the critical value of $a$ the local eigenvalue correlations are given in terms of Pearcey integrals .

In this paper we study the external source model (1.1) with a more general potential $V$. We assume that $V$ is an even polynomial

For $a=0$ the external source model (1.1) reduces to the usual unitary matrix model

which is one of the most studied models in random matrix theory in both mathematics and physics, see e.g. for rigorous study using the Riemann-Hilbert approach. A basic fact is that for $n\to\infty$, the limiting mean eigenvalue distribution of the matrix $M$ in (1.5) minimizes the energy functional

It is an open problem to find an analogue for the equilibrium problem (1.6) in the general context of the random matrix model with external source. This paper contains a first result in this direction. We consider the external source model (1.1) in case where the potential $V$ is an even polynomial (1.4). The external source $A$ is again given by (1.3) with two eigenvalues $\pm a$ of equal multiplicity. We show that under these assumptions, the limiting mean eigenvalue distribution of the matrix $M$ in (1.1) exists, and that it arises as the first component of a pair of measures $(\mu_{1},\mu_{2})$ solving a certain vector equilibrium problem, see Section 2.

We will illustrate our results in detail for a particular case of a non-convex potential, namely the quartic double well potential

For the quartic model (1.5), (1.7) (without external source) it is known that the eigenvalues accumulate on either one or two intervals. The local eigenvalue correlations for the critical value of $t=t_{cr}=2$ are given in terms of $\Psi$-functions associated with the Hastings-McLeod solution of the Painlevé II equation .

So in the quartic model with external source there exist at least two mechanisms by which a transition from one to two intervals can occur: namely a Pearcey transition and a Painlevé II transition. It will be one of the outcomes of the present paper that we can determine precisely the location of the phase transitions in the $ta$-plane.

Statement of results

The main ingredient in our analysis is a new vector equilibrium problem associated with the random matrix model (1.1) with external source. We emphasize that it only applies in the setting we are considering, namely an even polynomial potential $V$ as in (1.4) and an external source (1.3) with two eigenvalues of equal multiplicity. This setting gives a symmetry with respect to the origin, which we use in an essential way.

The equilibrium problem is as follows. We minimize the energy functional

with respect to all pairs of measures $(\mu_{1},\mu_{2})$ satisfying

$\mu_{1}$ and $\mu_{2}$ have finite logarithmic energy,

Standard references on potential theory in the complex plane are .

The equilibrium problem (2.1) has both an external field $V(x)-a|x|$ acting on $\mu_{1}$, and an upper constraint $\sigma$ acting on $\mu_{2}$. The interaction between $\mu_{1}$ and $\mu_{2}$ is of Nikishin type . This type of vector equilibrium problem also appeared recently in a model of non-intersecting squared Bessel paths and in the two-matrix model with quartic potential .

Our first result concerns the structure of the minimizer of the equilibrium problem.

There is a unique minimizer $(\mu_{1},\mu_{2})$ which satisfies

The support of $\mu_{1}$ is bounded and consists of a finite union of intervals

The measure $\mu_{1}$ is absolutely continuous with density

where $h$ is a nonnegative function on $\operatorname*{supp}(\mu_{1})=\bigcup_{j=1}^{N}[a_{j},b_{j}]$ that is real analytic, except possibly at zero.

The support of $\mu_{2}$ is the full imaginary axis and there exists $c\geq 0$ such that

and in that case, $\mu_{2}$ has the density

If (2.7) is not satisfied then $c>0$ is determined by the condition

Both $\mu_{1}$ and $\mu_{2}$ are symmetric with respect to the origin.

Theorem 2.1 will be proved in Sections 3.1 and 3.2. Note that (2.4)–(2.8) are similar to statements proved in , while (2.9) and (2.11)–(2.12) have apparently not been stated before.

2 Variational conditions

The minimizer $(\mu_{1},\mu_{2})$ to the equilibrium problem in Section 2.1 is characterized by the following Euler-Lagrange variational conditions. We write

for the logarithmic potential of a measure $\mu$.

These relations follow directly from the variational conditions of the equilibrium problem.

3 Regular and singular cases

For the analysis in this paper, we will assume that both $\mu_{1}$ and $\mu_{2}$ are regular. The following lemma follows immediately from this assumption and is stated only for further reference.

The measures $\mu_{1}$ and $\sigma-\mu_{2}$ satisfy the following square root behavior near their endpoints $a_{j}$, $b_{j}$ and $\pm ic$:

If the measure $\mu_{1}$ is regular then it has a density of the form (2.5) with $h$ strictly positive on $\bigcup_{j=1}^{N}[a_{j},b_{j}]$.

If $c>0$ then the measure $\sigma-\mu_{2}$ has a density of the form

where $k$ is an analytic, strictly positive function on $(-i\infty,-ic]\cup[ic,i\infty)$.

Lemma 2.3(a) follows immediately from the definition of $\mu_{1}$ being regular. Lemma 2.3(b) follows from (2.11).

4 Limiting eigenvalue distribution

Our main result deals with the global distribution of eigenvalues as $n\to\infty$.

Let $V$ be an even polynomial, and let $A$ be a diagonal matrix with two eigenvalues $\pm a$ of equal multiplicity. Let $(\mu_{1},\mu_{2})$ be the solution of the equilibrium problem in Section 2.1, and assume that both $\mu_{1}$ and $\mu_{2}$ are regular in the sense explained above. Then the mean eigenvalue distribution of a matrix $M$ from the random matrix model

We strongly expect that the conclusion of Theorem 2.4 remains valid in the case where $\mu_{1}$ and/or $\mu_{2}$ is singular.

Theorem 2.4 will be proved in Section 5.7.

5 About the proof

The proof of Theorem 2.4 is based on the Riemann-Hilbert problem for multiple orthogonal polynomials and its connection with the external source model (1.1).

The multiple orthogonal polynomials $P_{n}(x)$ in question are orthogonal with respect to the weights

More precisely, $P_{n}(x)$ is a monic polynomial of degree $n$ that is characterized by the multiple orthogonality conditions (we assume $n$ is even)

The polynomial $P_{n}(x)$ is also the average characteristic polynomial

The RH problem (2.20)–(2.21) has a unique solution. The $(1,1)$-entry of $Y(z)$ is the multiple orthogonal polynomial $P_{n}(z)$ characterized by (2.18).

It is known that the eigenvalues of the random matrix model with external source (1.3) form a determinantal point process with correlation kernel

Theorem 2.4 then comes down to the following statement about the limiting behavior of the kernels $K_{n}$:

We will establish (2.23) in Section 5.7, thereby proving Theorem 2.4.

From the RH analysis it is possible to obtain universality results for the local eigenvalue correlations as well. In the regular cases, this leads to the usual sine kernel in the bulk and Airy kernel at the edge points of the spectrum. We will not discuss this any further and refer to the papers , among others, for a detailed analysis in a similar context.

6 Organization of the paper

The rest of the paper is organized as follows. In Section 3 we discuss the structure of the equilibrium measures and we prove Theorem 2.1. In Section 4 we introduce the Riemann surface built from the solution of the equilibrium problem. Section 5 contains the steepest descent analysis of the RH problem for $Y(z)$, leading to the proof of Theorem 2.4. In Section 6 we make some general remarks on the expected phase transitions of our model, and in Section 7 we study this in detail for the case of a quartic potential.

The equilibrium problem

In this section we prove the existence of the minimizer $(\mu_{1},\mu_{2})$ of the equilibrium problem in Section 2.1. To this end we follow [20, Section 4].

The energy functional (2.1) can be written as

denotes the logarithmic energy of a signed measure $\nu$. Occasionally we will also write

to denote the mixed energy of a pair of measures $\nu_{1}$ and $\nu_{2}$.

Since $I(\nu)\geq 0$ if $\nu$ is a signed measure with $\int d\nu=0$, we find from (3.1) that

where the last inequality follows from standard logarithmic potential theory with external fields, see e.g. . Thus the energy functional is bounded from below.

The extra term $U^{\mu_{2}}(x)$ comes from the interaction between $\mu_{1}$ and $\mu_{2}$. It is a term that attracts the $\mu_{1}$ mass towards the origin. It can indeed be proved (as in ) that if the minimizer in external field $V(x)-a|x|$ is contained in $[-X,X]$, then the minimizer in external field $V(x)-a|x|-U^{\mu_{2}}(x)$ is also contained in $[-X,X]$ (and so $X$ is independent of $\mu_{2}$).

If we fix $\mu_{1}$ on $[-X,X]$ then the problem for $\mu_{2}$ is to minimize

among all $\mu_{2}\leq\sigma$ with total mass $1/2$. As in , equality in the constraint is attained precisely on an interval of the form $[-ic,ic]$ for certain $c\geq 0$. We will show further that the minimizer $\mu_{2}$ satisfying this constraint is given explicitly by (2.8)–(2.12). From these explicit formulas it follows immediately that for a measure $\mu_{1}$ on $[-X,X]$, the corresponding minimizer $\mu_{2}$ satisfies

with a constant $K$ that only depends on $X$.

As shown above, we may assume in addition that

Then it follows as in that the sequences $(\mu_{1,n})_{n=1}^{\infty}$ and $(\mu_{2,n})_{n=1}^{\infty}$ are tight. There is a convergent subsequence of $(\mu_{1,n},\mu_{2,n})_{n=1}^{\infty}$ and the limit is the vector of minimizing measures, see also .

Summarizing, we have now proved the existence of the solution $(\mu_{1},\mu_{2})$ to the equilibrium problem. The uniqueness of the solution follows in a standard way from the convexity of the energy functional, see e.g. (3.1) and . ∎

2 Proof of Theorem 2.1

The proof of Theorem 2.1(a) follows as in , while Part (c) is evident from the symmetry of the problem.

It remains to prove (2.7)–(2.12) in Theorem 2.1(b). For $k=1,2$, define the Cauchy transforms

By differentiating the variational condition (2.15) we find that

Here we assume that the imaginary axis is oriented from bottom to top, so that the $+$-side is on the left, and the $-$-side is on the right, as usual. If $c>0$ then $\frac{d\mu_{2}}{|dz|}=\frac{d\sigma(z)}{|dz|}=\frac{a}{\pi}$ on $(-ic,ic)$, from which it follows that

We can solve equations (3.3)–(3.4) for $F_{2}$. We consider the two cases: $c=0$ and $c>0$.

Define the Cauchy transforms of the restrictions of the measure $\mu_{1}$ to the positive and negative half-axes,

and, due to the uniqueness of the solution of the scalar Riemann-Hilbert problem (3.5), we have

which is equivalent to (2.8). The case $c=0$ is valid if and only if the density (2.8) is bounded by $a/\pi$. Since (2.8) assumes its maximum for $|z|=0$, this happens if and only if (2.7) holds.

where $\sqrt{R(z)}=\sqrt{z^{2}+c^{2}}$ is defined with a cut $(-i\infty,-ic]\cup[ic,i\infty)$, and $\sqrt{R(0)}=c$. From equation (3.3) we obtain that

where $+$ again denotes the limiting value from the left half plane. Observe that for $y>c$,

In addition, from equation (3.4) we obtain that

By (3.6) and (3.7) the Sokhotski-Plemelj formula implies that

The first term in the right-hand side of (3.8) can be written as

where the contour $\Gamma$ is depicted in Figure 1. From (3.2) and Fubini’s theorem,

By contour deformation and Cauchy’s theorem we have that

For the second term in the right-hand side of (3.8) we have

By inserting (3.11) and (3.12) in (3.8), we obtain that

Note that by taking $z\to+\infty$ in (3.13), we find the relation (2.9) between $c$ and $a$. Now the density of $\mu_{2}$ is equal to

which is equivalent to (2.11). Then (2.12) follows from this and (2.9). ∎

3 Structure of the equilibrium measures in the regular case

Theorem 2.1 implies that in the regular case, the structure of the equilibrium measures near the origin is described by one of the following three cases. This distinction will be important at several places of our RH steepest descent analysis.

For the quadratic potential $V(x)=x^{2}/2$, it turns out that we are in Case I (with $N=2$) for large values of $a$ and in Case III (with $N=1$) for small values of $a$. The Case II does not occur.

Let $\mu_{1}^{L}$ and $\mu_{1}^{R}$ denote the restrictions of $\mu_{1}$ to the negative and positive real axis, respectively. By symmetry, we then have in Case I that $\mu_{2}$ is the balayage of either $\mu_{1}^{L}$ or $\mu_{1}^{R}$ onto the imaginary axis, and moreover

Note that the equality is valid not only on the imaginary axis, but also in a full half-plane. This follows from an easy application of the minimum and maximum principles for harmonic functions [35, Chapter 0].

Then the following string of equations is easy to verify:

it then follows from (3.16)–(3.18) that the energy functional (2.1) can be rewritten in Case I as

Therefore the equilibrium problem is equivalent to the following equilibrium problem of Angelesco type for $\mu_{1}^{L}$ and $\mu_{1}^{R}$: Minimize

with respect to all pairs of measures $(\mu_{1}^{L},\mu_{1}^{R})$ satisfying

For the quadratic potential $V(x)=\frac{1}{2}x^{2}$, this equilibrium problem was described in .

Riemann surface

From the minimizer $(\mu_{1},\mu_{2})$ of the vector equilibrium problem we construct a three sheeted Riemann surface $\mathcal{R}$, whose three sheets are given as follows.

The sheet $\mathcal{R}_{1}$ is connected with $\mathcal{R}_{2}$ via the intervals $[a_{j},b_{j}]$ on the positive real line, $\mathcal{R}_{1}$ is connected with $\mathcal{R}_{3}$ via the intervals $[a_{j},b_{j}]$ on the negative real line, and (in Case II and Case III) $\mathcal{R}_{2}$ is connected to $\mathcal{R}_{3}$ via the interval $[-ic,ic]$ on the imaginary axis. The connections are in the usual crosswise manner. The Riemann surface is compact and has genus

Here the Cases I, II and III were defined in Section 3.3. An illustration of the Riemann surface for each of these three cases is shown in Figures 2–4.

Recall the functions $F_{1}$ and $F_{2}$ in (3.2). These functions are used to define a meromorphic function on the Riemann surface, compare with [20, Lemma 5.1]:

For $k=1,2,3$, let $\xi_{k}(z)$ on the sheet $\mathcal{R}_{k}$ be defined by

Then these functions have an analytic continuation to a meromorphic function (denoted by $\xi(z)$) on the Riemann surface whose only pole is at the point at infinity on the first sheet.

Let us first check that $\xi_{k}$ is analytic on the sheet $\mathcal{R}_{k}$, $k=2,3$. This reduces to showing the equality

which is a direct consequence of the variational condition (2.15), see also (3.3).

which is a direct consequence of the variational condition (2.13). The other equalities are checked similarly, see also (3.4). ∎

It follows from Proposition 4.1 that the function $\xi(z)$ is an algebraic function satisfying an equation of the third degree in $\xi$, known as the spectral curve:

where $p_{0}$, $p_{1}$, $p_{2}$ are polynomials. Here

is known, but the determination of the polynomials

cannot be done in general. We can only certify that (we use that $V(z)$ is a polynomial of degree $2d$ and $F_{1}(z)=1/z+O(1/z^{3})$, $F_{2}(z)=1/(2z)+o(1/z)$ as $z\to\infty$)

For $d=1$ we are in the quadratic case. If

then $p_{1}(z)=1-a^{2}$, $p_{0}(z)=a^{2}z$ so that the spectral curve is

This is known as Pastur’s equation . It plays an important role in .

For $d=2$ we are in the quartic case. Let’s take

This is McLaughlin’s equation, named after K.T-R McLaughlin who derived it first for the case $t=0$, see also . We will analyze this case in more detail in Section 7 below.

Proof of Theorem 2.4

Recall the RH problem for $Y(z)$ in (2.20)–(2.21). In Subsections 5.1–5.6 we will perform a Deift-Zhou steepest descent analysis of this RH problem. This will then lead to the proof of Theorem 2.4 in Section 5.7.

In the first transformation we open up an unbounded lens around $\operatorname*{supp}(\sigma-\mu_{2})$ which is bounded by a contour $\Gamma=\Gamma^{+}\cup\Gamma^{-}$. We choose the lens so that it is symmetric under reflection with respect to both the real and the imaginary axis. The construction of the lens depends on whether we are in Case I or in one of the other cases (Case II or Case III).

Next, consider the Cases II and III. Then we take $\Gamma=\Gamma^{+}\cup\Gamma^{-}$ as in Figure 6. The part of $\Gamma^{+}$ in the upper half plane is a Jordan curve going from $\infty$ at an angle $\theta\in(0,\pi/2)$ to the point $ic$. The other part of $\Gamma^{+}$ is its reflection in the real axis, and $\Gamma^{-}$ is obtained from $\Gamma^{+}$ by reflection in the imaginary axis. We orient $\Gamma$ as shown in Figure 6.

The precise way to choose the contour $\Gamma$ will be described further on.

The contour $\Gamma$ divides the complex plane into an inner and an outer part. By definition, we say that $\operatorname*{supp}(\sigma-\mu_{2})$ is inside the lens, while $\operatorname*{supp}(\mu_{1})$ is outside the lens. Note that our definitions are such that the outside of the lens is always on the left when traversing $\Gamma$ according to its orientation.

We define a new $3\times 3$ matrix valued function $X$ by

Recall that $q$ is used to denote the intersection of $\Gamma^{+}$ with the real axis (in Case I). In Cases II and III we put $q=0$.

Then $X$ satisfies the following RH problem.

on the interval $(-ic,ic)$ oriented upwards, and

2 Second transformation X↦Tmaps-to𝑋𝑇X\mapsto T

In the second transformation $X\mapsto T$ we use the minimizers $\mu_{1}$ and $\mu_{2}$ of the vector equilibrium problem, and the associated $g$-functions

The behavior of the real and imaginary parts of the $g$-functions is described in the following lemma.

The equalities follow immediately from the definitions, where we have to be careful with the choice of branches of the logarithm as discussed above.

Only the equality in (5.11) for $z\in[-ic,ic]$ needs some extra comment. We have that $\mu_{2}=\sigma$ on $[-ic,ic]$ and so if $z\in[0,ic]$,

Note that $\lambda_{2}(z)$ and $\lambda_{3}(z)$ are defined with a cut along the entire imaginary axis. In fact, from (5.11)–(5.13) it follows that

and a similar formula holds for $\lambda_{2}$.

We can now reformulate Lemma 5.1 in terms of the $\lambda$-functions. This leads to the following two lemmas.

These are reformulations of the variational conditions (2.13)–(2.16) associated with the equilibrium problem, taking into account (5.10)–(5.13). The fact that we have strict inequalities follows from our assumption that the measure $\mu_{1}$ is regular. ∎

The last expression is purely imaginary and its imaginary part is strictly increasing in terms of $\operatorname{Im}(z)$ as $z\in(-i\infty,-ic)\cup(ic,i\infty)$.

This is a straightforward calculation using (5.10)–(5.13). ∎

We define the new matrix valued function $T$ as

Then $T$ satisfies the following RH problem.

But from the jump relations in (5.14) we see that

by the fact that $n$ is even. In a similar way one shows that $e^{n(\lambda_{2,+}(z)-\lambda_{2,-}(z))}=1$, and so the jump matrix in (5.24) is just the identity matrix.

where $b_{0}=-\infty$, $a_{N+1}=+\infty$ and

We would like the jump matrices on $\Gamma$ to be exponentially close to the identity matrix as $n\to\infty$. From the above RH problem, we see that this is achieved provided $\Gamma$ lies in the region where $\operatorname{Re}(\lambda_{3}-\lambda_{2})<0$ if $\operatorname{Re}z>0$ and $\operatorname{Re}(\lambda_{3}-\lambda_{2})>0$ if $\operatorname{Re}z<0$. The fact that $\Gamma$ can indeed be chosen in this way, follows by applying the Cauchy-Riemann equations to the last equality in Lemma 5.3, and using the last line in the statement of that lemma.

Around each of the intervals $[a_{j},b_{j}]$ we open up a small lens to transform the oscillatory entries of the jump matrix into exponentially decaying entries. Since the non-trivial part of the jump matrix is locally of size $2\times 2$ only, this can be done in the standard way .

More precisely, we take Jordan curves $\Gamma_{j}^{+}$ and $\Gamma_{j}^{-}$ surrounding the interval $[a_{j},b_{j}]$ as in Figure 7. The region between these curves is called the lens, and $\Gamma_{j}^{+}$ and $\Gamma_{j}^{-}$ are the upper and lower lip of the lens, respectively. We choose them sufficiently close to the real axis so that

The fact that this is possible, follows from applying the Cauchy-Riemann equations to the first equation of (5.16), cf. .

Then $S$ satisfies the following RH problem.

For $x\in\operatorname*{supp}(\mu_{1})$ we have that

On the lips $\Gamma_{j}$ of the lenses we have

The jumps of $S(z)$ on the other contours are the same as those for $T(z)$.

From (5.15) and (5.30)–(5.31), it can be checked that all the non-constant entries in the jump matrices for $S(z)$ tend to as $n\to\infty$, uniformly for $z$ bounded away from the branch points $a_{j}$, $b_{j}$, $\pm ic$. The only case that requires more explanation is the $(3,1)$ entry in the jump matrix in (5.33). In that case, one can factorize

and observe from (5.15) and (5.31) that for $z\in[-iz_{0},iz_{0}]$, the leftmost factor is uniformly bounded by $1$ while the rightmost factor is uniformly exponentially decaying as $n\to\infty$.

4 Global parametrix

The global parametrix we look for is a $3\times 3$ matrix valued function $M$ with jumps (obtained from the jumps of $S$ by ignoring all entries which are exponentially small for $n\to\infty$)

$M$ has at most fourth-root singularities at the branch points $a_{1},b_{1},\ldots,a_{N}$, $b_{N}$, $\pm ic$.

We can solve this problem with the help of meromorphic differentials on the Riemann surface. Such a construction was first used in and later developed further in .

To the Riemann surface we associate a canonical homology basis $\{A_{1},\ldots,A_{g}$, $B_{1},\ldots,B_{g}\}$ where $g$ is the genus. The details of the construction depend on whether we are in Case I, II or III, see Figures 8–10.

For brevity, we give a detailed description only for Cases II and III. Then the genus is $g=N-1$. The cycles $B_{j}$ are on the first sheet and $B_{j}$ encircles $[a_{1},b_{j}]$ once in the counterclockwise direction. The cycles $A_{j}$ are partly in the upper half-plane on the first sheet and partly in the lower half-plane on the second or third sheet. $A_{j}$ passes through $[a_{j},b_{j}]$ and $[a_{j+1},b_{j+1}]$.

The anti-holomorphic involution $\phi$ is defined by mapping $z$ to $\overline{z}$ on the same sheet. The fixed point set of $\phi$ is the disjoint union of $N=g+1$ closed curves $\bigcup_{j=0}^{g}\Sigma_{j}$ on the Riemann surface. Here $\Sigma_{j}$ is homotopic to $A_{j}$ as a closed curve, $j=1,\ldots,g$, while $\Sigma_{0}$ is the unbounded component.

If $P_{j}\in\Sigma_{j}$ for $j=1,\ldots,g$, then the divisor

The proof is based on the following result which can be found e.g. in [19, Theorem 2.4.2]. The divisor $D$ is non-special if and only if $\theta(u(P)-u(D)-\vec{K})$ does not vanish identically for $P$ on the Riemann surface.

Using the antiholomorphic involution $\phi$ it can be shown that the Riemann period matrix is purely imaginary. Finally, one then shows that $u(D)+\vec{K}$ has a real representative modulo the lattice $L$ . By taking into account the results in the last two paragraphs, the desired result then follows. ∎

We now basically follow . Given $(P_{1},\ldots,P_{g})$ with $P_{j}\in\Sigma_{j}$, we define a meromorphic differential $\omega_{P}$ so that

$\omega_{P}$ has simple poles in $a_{1},b_{1},\ldots,a_{N},b_{N}$, $\pm ic$, $P_{1},\ldots,P_{g}$, $\infty_{2}$ and $\infty_{3}$ with residues

is a bijection. These claims follow in the same way as in .

Thus there exist $P_{j}^{(1)}\in\Sigma_{j}$ so that

Let $\omega_{P}^{(1)}$ be the corresponding meromorphic differential.

We take the base point $P_{0}=\infty_{1}$ and define three functions $v_{1}(z)$, $v_{2}(z)$ and $v_{3}(z)$ of a complex variable $z$ as follows. We have

where $z$ is considered as a point on the $k$th sheet of $\mathcal{R}$, and where the path of integration is as follows

for $k=1$, the path of integration is on the first sheet and does not intersect the real line,

where the jump matrices on the real line are

and all $v_{j}$ functions have fourth-root singular behavior at the branch points $a_{1},b_{1},\ldots,a_{N},b_{N}$, $\pm ic$.

Now define the first row $(M_{11},M_{12},M_{13})$ as follows

This gives the correct jumps for $M$. The other rows of $M$ can be constructed in a similar way, or by a simple transformation of the first row .

5 Local parametrices

In the regular case (the one we are considering) we construct local parametrices out of Airy functions near each of the branch points. We denote the local parametrices by $P$. The local parametrices match with the global parametrix on the boundary of a small circle around the branch points. Since the non-trivial part of the jump matrix is locally of size $2\times 2$ only, and since in the regular case we have square root behavior near all the branch points (Lemma 2.3), this construction can be done in the usual way . We omit the details.

6 Final RH problem

7 Proof of Theorem 2.4

Having performed the steepest descent analysis of the RH problem, we can now prove Theorem 2.4. The proof follows the same pattern as in the papers .

First assume that $x,y\in(a_{j},b_{j})$ with $x,y>0$. We will transform (2.22) under the series of transformations $Y\mapsto X\mapsto T\mapsto S$. By virtue of (5.1) we have

Now it follows by standard arguments (e.g. [8, Section 9]) that

uniformly in $n$. Inserting this into (5.39) and setting $h(x)=V(x)-\operatorname{Re}(\lambda_{1}(x))$ yields

uniformly in $n$. By letting $y\to x$ and using l’Hôpital’s rule we find

as $n\to\infty$. From (4.5) and the Stieltjes inversion principle we conclude that

Phase transitions: General discussion

Recall the Cases I, II and III describing the structure of the equilibrium measures in Section 3.3. Intuitively one expects the following possible behavior in terms of the parameter $a$. For large $a$ we are in Case I. The measure $\mu_{1}$ is supported on two (or more) disjoint intervals with a gap around . The constraint $\sigma$ is not active.

When $a$ decreases the gap around shrinks. Then one of two things could happen. It could happen that for a certain value of $a$ the constraint becomes active, while the gap in the support of $\mu_{1}$ around is still there. Then we are in Case II. Then if $a$ further decreases the gap may close or not. The latter depends on whether in the unitary matrix model with potential $V$ (without external source) is in the support or not. If the support closes then we are in Case III.

The other situation that could happen is that the constraint $\sigma$ remains inactive all the way until for a certain value of $a$ the gap in the support of $\mu_{1}$ is closed. Then if $a$ further decreases the constraint becomes active. The transition is then from Case I to Case III without passing through the Case II. This is precisely what happens in the quadratic case $V(x)=\frac{1}{2}x^{2}$. More generally, one expects this kind of behavior when the potential $V(x)$ is convex, or ‘nearly’ convex.

For those values of $a$ for which a transition between one of the Cases I, II, III to another takes places, one expects that the local eigenvalue correlations near the origin are described by special functions related to ODE’s. For typical cases one expects such special functions as Pearcey integrals and the Hastings-McLeod solution to the Painlevé II equation . However, our model allows for new kinds of critical and multi-critical behavior as well, but it remains an open problem to describe these new critical phenomena.

In Section 7 we will illustrate the above considerations in detail for the case of a quartic potential. See in particular Figure 11.

A case study: The quartic potential

Let us investigate the case of a quartic potential

and the associated McLaughlin equation (4.9):

The discriminant of the McLaughlin equation (w.r.t. $\xi$) is a polynomial $D_{12}(z)$ of degree $12$ in $z$. We calculated it with Maple, but it is too long and not too interesting to reproduce it here in full. The first terms are

The branch points of the Riemann surface are among the zeros of $D_{12}(z)$. There are other zeros, and they should come with higher multiplicities.

For general $\alpha$ and $\beta$ the McLaughlin equation has genus $4$ (according to Maple). The special choices for $\alpha$ and $\beta$ that are relevant to us will lead to a reduction of the genus. The genus can be at most one, as the following lemma shows.

A similar result occurs in [20, Prop. 5.2.5], but the proof given there is incorrect. Here we give a self-contained proof which may be used for the situation in as well. The proof uses an idea due to Lun Zhang (personal communication).

Assume that $x\mapsto V(\sqrt{x})$ is convex for $x>0$. Then the support of $\mu_{1}$ is either one interval (in case $0\in\operatorname*{supp}(\mu_{1})$) or a disjoint union of two intervals (in case $0\not\in\operatorname*{supp}(\mu_{1})$), and the measure $\mu_{1}$ can have singular behavior only at zero.

By fixing $\mu_{2}\leq\sigma$ in the energy functional (2.1), we see that $\mu_{1}$ is the unique minimizer for the energy functional

over all probability measures $\mu$ on $[0,\infty)$.

We are going to show that the external field $V(\sqrt{x})-a\sqrt{x}-U^{\mu_{2}}(\sqrt{x})$ in (7.3) is convex for $x>0$. Since

which due to the constraint $\mu_{2}\leq\sigma$ can be bounded by

It follows that $-a\sqrt{x}-U^{\mu_{2}}(\sqrt{x})$ is convex. Due to the assumption in the lemma on the convexity of $V(\sqrt{x})$, it then follows that $\left(V(\sqrt{x})-a\sqrt{x}-U^{\mu_{2}}(\sqrt{x})\right)$ is indeed convex for $x>0$.

The three cases in Section 3.3 now specialize as follows.

In the first case there are four real branch points $\pm b_{1}$, $\pm b_{2}$ (with $b_{1}>b_{2}>0$) and no other branch points. The remaining zeros of the discriminant (7.2) come as four double zeros. The Riemann surface has genus zero.

Case II:

In the second case there are four real branch points $\pm b_{1}$, $\pm b_{2}$ (with $b_{1}>b_{2}>0$) and two purely imaginary branch points $\pm ic$ (with $c>0$). The remaining zeros of the discriminant come in the form of a six-fold zero at $z=0$, see Lemma 7.2 below. The Riemann surface has genus one.

Case III:

In the third case there are two real branch points at $\pm b_{1}$ (with $b_{1}>0$) and two purely imaginary branch points at $\pm ic$ (with $c>0$). There are four double zeros of the discriminant. The Riemann surface has genus zero.

The genus one region

Let us first investigate Case II. This is the case of genus 1. It turns out that in this case, the parameters $\alpha$ and $\beta$ in the McLaughlin equation (7.1) take on a particularly simple form: they are both equal to zero. This is the content of the next lemma.

where $D_{6}$ is the degree six polynomial

From the general descriptions above we know that in the genus 1 case, the McLaughlin equation has six simple branch points $\pm b_{1}$, $\pm b_{2}$, and $\pm ic$ ($b_{1},b_{2},c>0$), which are simple zeros of the discriminant. The remaining zeros of the discriminant should come as three double zeros, possibly coalescing. By symmetry is a double zero, and this forces $\alpha=0$, cf. (7.2).

If we substitute $\alpha=0$ in (7.1) and calculate the discriminant with respect to $\xi$ we obtain the $12$th degree polynomial

In conclusion, we have shown now that $\alpha=\beta=0$. Inserting this in the McLaughlin equation and computing its discriminant by a direct calculation then leads to (7.4)–(7.5). ∎

Note that the factor $z^{6}$ in (7.4) corresponds to the six-fold zero at $z=0$, while the zeros of $D_{6}(z)$ should yield the branch points $\pm b_{1}$, $\pm b_{2}$ and $\pm ic$. In particular, four of these zeros should be real and the other two purely imaginary. The next lemma describes when this happens.

that connect the point $(t_{2},a_{2})$ with $(t_{1},a_{1})$ and $(t_{3},a_{3})$, respectively. The region $\mathcal{D}$ is shown in the bottom right part of Figure 11.

Rewrite $D_{6}(z)$ as a cubic polynomial in the variable $y=z^{2}$:

The two previous lemmas show that the genus of the McLaughlin equation can only be 1 if $(t,a)$ lies in the region $\mathcal{D}$. Outside $\mathcal{D}$ the genus must necessarily be zero.

It remains to show that inside $\mathcal{D}$ the genus is exactly 1 (and not 0). This is taken care of by the next lemma.

Inside $\mathcal{D}$ the genus is either identically $1$ or identically .

There exists at least one point $(t,a)$ in $\mathcal{D}$ for which the genus is $1$.

Now let $\mathcal{D}_{1}\subset\mathcal{D}$ be the region formed by those $(t,a)\in\mathcal{D}$ for which the genus is 1. We show that $\mathcal{D}_{1}$ is both open and closed in $\mathcal{D}$. To show that it is open, let $(t,a)\in\mathcal{D}_{1}$. Then the discriminant has six simple zeros and by continuity the same holds in an open neighborhood of $(t,a)$. To show that $\mathcal{D}_{1}$ is closed in $\mathcal{D}$, we take a sequence of points $(t_{k},a_{k})\in\mathcal{D}_{1}$, $k=1,2,\ldots$, which converge to a limit point $(t,a)\in\mathcal{D}$. By Lemma 7.2 we have $\alpha=\beta=0$ for each $(t_{k},a_{k})$ so by continuity the same must hold for the limit point $(t,a)$. But then Lemma 7.3 shows that the discriminant has six distinct simple zeros, which implies that the genus is 1. Hence $(t,a)\in\mathcal{D}_{1}$.

For Part (b), we only outline a proof. The idea is to show that for any fixed $t>2$, we have $(t,a)\in\mathcal{D}_{1}$ for all $a$ small enough. This relies on the fact that for $t>2$ the eigenvalues in the unitary matrix model with potential $\frac{1}{4}x^{4}-\frac{t}{2}x^{2}$ (without external source) are supported on two intervals . The claim then follows from a continuity argument for $a\to 0$; we do not go into the details.

An alternative approach to prove Part (b) would be to pick a numerical point $(t,a)\in\mathcal{D}$ and show by direct means (using the McLaughlin equation with $\alpha=\beta=0$) that this algebraic curve makes the RH steepest descent analysis work, in a similar vein as in . ∎

We summarize our findings with the following

The genus zero region

According to , $c$ is actually the largest positive root of equation (7.11), but we will not need this in what follows.

Thus the equations (7.11)–(7.12) have a common root $c$. In other words, the resultant of these two equations with respect to the variable $c$ should be zero. Computing this resultant with Maple yields the following condition on $(t,a)$:

Next we consider the case $c^{4}+2c^{4}u-u^{2}=0$. From (7.8)–(7.9) this implies that $\alpha=\beta=0$ and we know from earlier considerations (or from a similar resultant calculation as above) that this is only possible if $(t,a)$ is such that (7.7) holds. The phase transition on this curve will be discussed in the next section. The phase transition on the curve (7.13) will be discussed in the section thereafter.

Painlevé II transition

At the two curved boundaries of the region $\mathcal{D}$ in Figure 11, we have a transition from genus 0 to genus 1. Recall that these boundaries are described by the relevant branches of the equation

More precisely, these branches are given by

We have $A_{1}>0$ for every $t\geq\sqrt{3}$ whereas $A_{2}>0$ only for $\sqrt{3}\leq t<2$ and $0<A_{2}\leq A_{1}$ for these $t$-values. See Figure 11.

On the above curves we have a transition from genus 0 to genus 1, and we expect that the phase transition is of Painlevé II type . More precisely, we expect that the following happens. If one lets $a$ decrease towards the curve $a^{2}=A_{1}$ (with $t>\sqrt{3}$), then the constraint $\sigma$ on the imaginary axis becomes active, and we have a transition from Case I to Case II. If one further decreases $a$ towards the curve $a^{2}=A_{2}$ ($\sqrt{3}<t<2$), then the gap in the support of $\mu_{1}$ closes and hence we have a transition from Case II to Case III.

On the curve $a^{2}=A_{2}$ the phase transition involves the eigenvalue measure $\mu_{1}$ and therefore we expect Painlevé II behavior in the local eigenvalue correlations at the origin . On the curve $a^{2}=A_{1}$, however, the phase transition takes place on the ‘non-physical’ sheets of the Riemann surface and therefore it is not felt in the eigenvalue statistics. But then we expect Painlevé II behavior in the recurrence coefficients for the associated multiple orthogonal polynomials, as in .

Note that for $t>2$ the transition at $a^{2}=A_{2}$ does not occur. This is consistent with the fact that for $t>2$ the eigenvalues in the unitary matrix model with potential $\frac{1}{4}x^{4}-\frac{t}{2}x^{2}$ (without external source) are supported on two intervals , as mentioned before.

Pearcey transition

The branch $A_{4}$ is negative, and so is irrelevant for us. The other branch $A_{3}$ is positive and we expect that for $t<\sqrt{3}$ a phase transition of the Pearcey type takes place for $a^{2}=A_{3}$. See Figure 11.

Acknowledgements

The first author is supported in part by the National Science Foundation (NSF) Grant DMS-0652005.

The second author is a Postdoctoral Fellow of the Fund for Scientific Research - Flanders (Belgium).

The third author is supported in part by FWO-Flanders project G.0427.09, by K.U. Leuven research grant OT/08/33, by the Belgian Interuniversity Attraction Pole P06/02, by the European Science Foundation Program MISGAM, and by grant MTM2008-06689-C02-01 of the Spanish Ministry of Science and Innovation.

References