Painlev\'e IV Critical Asymptotics for Orthogonal Polynomials in the Complex Plane

We study the asymptotic behaviour of orthogonal polynomials in the complex plane that are associated to a certain normal matrix model. The model depends on a parameter and the asymptotic distribution of the eigenvalues undergoes a transition for a special value of the parameter, where it develops a corner-type singularity. In the double scaling limit near the transition we determine the asymptotic behaviour of the orthogonal polynomials in terms of a solution of the Painlev\'e IV equation. We determine the Fredholm determinant associated to such solution and we compute it numerically on the real line, showing also that the corresponding Painlev\'e transcendent is pole-free on a semiaxis.


Introduction
In this work we consider the critical asymptotic behaviour for the orthogonal polynomials that are associated to the matrix model described here. Over the set of normal n × n matrices dA(λ 1 ) · · · dA(λ n ), (1.2) where λ j are the complex eigenvalues of the normal matrix M [11] and dA(λ) = d Re(λ)d Im(λ) is the area element in the complex plane. Associated to the above data we consider the sequence of (monic) orthogonal polynomials p n (z) = z n + · · · , C p n (λ)p m (λ)e −N W (λ) dA(λ) = h n,N δ n,m , h n,N > 0, n, m = 0, 1, . . . , (1.3) where W : C → R is the external potential. Likewise in the standard GUE, the orthogonal polynomial p n (z) is also the expectation value of the characteristic polynomial of the random normal matrix M : p n (z) = det(z1 − M ) , where the expectation is with respect to the measure (1.1) (see, e.g., [18]). The statistical quantities related to the eigenvalues of matrix models can be expressed in terms of the associated orthogonal polynomials, p n (λ). In particular, the probability measure (1.2) can be written as a determinantal point field (sometimes also called "process") with Christoffel-Darboux kernel K(λ, η) = e −N (W (λ)+W (η))/2 The average density of eigenvalues is ρ n,N (λ) = 1 n K(λ, λ) and in the limit it converges to the unique probability measure, µ * , in the plane which minimizes the functional [19,26] I(µ) = log |λ − η| −1 dµ(λ)dµ(η) + 1 T W (λ)dµ(λ).
The functional I(µ) is the Coulomb energy functional in two dimensions and the existence of a unique minimizer is a well-established fact under mild assumptions on the potential W (λ) [37]. If W (λ) is twice continuously differentiable then the equilibrium measure is given by where χ D is the characteristic function of the compact support set D = supp(µ * ) and ∆ = 4∂ λ∂λ . For example, if W (λ) = |λ| 2 , one has the Ginibre ensemble [23] and the measure dµ * (λ) is the uniform measure on the disk of radius √ T . Following [4], we consider the external potential to be of the form (1.4) The case d = 1 is reducible to Hermite polynomials and we will consider only d ≥ 2. Up to a rotation λ → e iθ λ we can assume without loss of generality that t ∈ R + , leading to the external potential and orthogonality measure W (λ) = |λ| 2d − t λ d +λ d , λ ∈ C, e −N (|λ| 2d −t(λ d +λ d )) dA(λ).
(1.5) There remains a residual discrete rotational Z d -symmetry which allows us to further reduce the problem as observed in [5] and [20]. The potential W in (1.5) can be written as W (λ) = 1 d Q λ d , and hence the equilibrium measure for W can be obtained from the equilibrium measure associated to Q by an unfolding procedure. Considering the particular case of (1.4), the potential Q can be rewritten as which corresponds to the Ginibre ensemble whose equilibrium measure is the normalized area measure of the disk centered at u = t and of radius t c |u − t| ≤ t c , t c = T d .
Pulling back to the λ-plane, we obtain the equilibrium measure of W where dA is the area measure and χ D is the characteristic function of the domain D whose boundary is a lemniscate (1.6) A plot of the domain D for different values of t is shown in Fig. 1. Depending on the radius t c , the domain D is either connected for t ≤ t c or has d connected components for t > t c . The cases t > t c and t < t c were analysed in [4]. In this manuscript we study the critical regime t = t c . Boundary points where ∆W > 0 have been classified by Sakai [38] and fall in the category of regular points, cusps or double points. Details about the application of Sakai's theory to the setting when ∆W > 0 near a singular point can be found in several papers, e.g., [33]. The classification of boundary points when ∆W = 0 is still an open problem to the best of our knowledge. In this manuscript we take a first step in considering a specific example when ∆W = 0. In particular we will study the asymptotic behaviour of the orthogonal polynomials associated to this model. The asymptotic analysis of orthogonal polynomials when the measure is of planar type, is much harder than the orthogonality on contours. When the orthogonality planar measure is not varying there are several studies like in [36] or [24]. In a very recent work [27] the asymptotic of orthogonal polynomials for a varying planar measure was obtained in the regular case. Critical regimes have been considered in [34], in [9,32] studying a normal matrix model with cut-off, and in [1] where new types of determinantal point fields, have emerged. The equilibrium problem for a class of potentials to which our case belongs, was also recently considered in [2]. In this paper we derive the asymptotic behaviour for large n of the orthogonal polynomials associated to the exponential weight (1.5) in the critical regime t c . The zeros of the polynomials accumulate on the Szegö curve that was first observed in relation to the zeros of the Taylor polynomials of the exponential function [40]. For large n we obtain the expansion of the orthogonal polynomials in terms of a special solution of a Painlevé IV transcendent. We show that such solution is pole free on the negative semi-axis. We obtain the Fredholm determinant associated to such solution and compute it numerically on the real line. Painlevé IV critical asymptotic was observed for orthogonal polynomials with deformed Laguerre weights [15]. Finally we remark that normal matrix models are strictly related to the 2-dimensional Toda lattice, see, e.g., [41,43]. It would be interesting to study the critical case considered in this paper in this perspective as done in [6] or [12,13,15] for the one-dimensional Toda lattice.

Zeros of orthogonal polynomials: statement of the result
The zeros of the orthogonal polynomials (1.3) in the critical case t = t c concentrate on the curveĈ that is described below. Let us introduce the function and consider the level curveĈ Observe that the level curveĈ consists of a closed contour contained in the set D defined in (1.6) with t = t c because it is the pullback under the map z = 1 − λ d tc of the celebrated (closed) Szegö curve ze 1−z = 1 [40] (see Fig. 2).
Let us define the measureν, associated with this family of curves, given by and supported onĈ.
Lemma 1.1. The a-priori complex measure dν in (1.9) is a probability measure on the contourĈ defined in (1.8).
When t = t c , the zeros of the orthogonal polynomials (1.3) behave in the following way. 1 d are zeros of the polynomials p kd+d−1 with multiplicity k and λ = 0 is a zero with multiplicity d − 1.
• for n = kd + , = 0, . . . , d − 2 the polynomial p n (λ) has a zero in λ = 0 with multiplicity and the remaining zeros in the limit n, N → ∞ such that accumulate on the level curveĈ defined in (1.8). The measureν in (1.9) is the weak-star limit of the normalized zero counting measure ν n of the polynomials p n for n = kd + , = 0, . . . , d − 2.
Our next result gives strong uniform asymptotics as n → ∞ for the polynomials p n (λ) in the complex plane. We consider the double scaling limit in such a way that lim k→∞, t→tc with S in compact subsets of the real line. In the description of the asymptotic behaviour of the orthogonal polynomials, p n (λ), in this double-scaling limit, the Painlevé IV transcendent (see Section 2) with Θ 0 = Θ ∞ = γ 2 and γ = d− −1 d ∈ (0, 1) plays a major role. (1) For λ in compact subsets of the exterior ofĈ one has with H(S) the Hamiltonian (2.7) of the Painlevé IV equation (2.3).
(2) For λ in the region nearĈ and away from the point λ = 0 one has (3) For λ in compact subsets of the interior region ofĈ and away from the point λ = 0 one has (4) In the neighbourhood of the point λ = 0 one has where Ψ is the solution of the Painlevé IV Riemann-Hilbert Problem 2.1 with Stokes multipliers specified in (2.9) and (2.10) and the regions Ω 0 , Ω ∞ and Ω 2 are specified in Fig. 4.
Here Ψ 11 and Ψ 12 are the entries of the matrix Ψ.
We observe that, in compact subsets of the exterior ofĈ, there are no zeros of the polynomials p n (λ). The only possible zeros are located in λ = 0 and in the region where the second term in parenthesis in the expression (1.10) is of order one. Since Reφ(λ) is negative insideĈ and positive outsideĈ, it follows that the possible zeros of p n (λ) lie insideĈ and are determined by the condition We conclude with the following proposition. This manuscript is organised as follows. In Section 2 we present the Lax pair and Riemann-Hilbert problem for the Painlevé IV equation and we show how to associate for a particular value of the Stokes matrices the solution of the Painlevé IV equation in terms of a Fredholm determinant. We compute such solution numerically for several values of the parameter γ that parametrizes the asymptotic behaviour of the orthogonal polynomials. In Section 3 we show how to associate to the orthogonal polynomials on the planar domain a Riemann-Hilbert problem and in Section 4 we perform the asymptotic analysis of this Riemann-Hilbert problem thus deriving the asymptotic behaviour of the orthogonal polynomials p n (λ) as n → ∞ in the complex plane.

Painlevé IV and Fredholm determinant
The purpose of this section is to introduce the Painlevé IV equation and study various aspects related to it that are connected with the work that we will be doing in the following sections.
We will begin by introducing and deriving the general Painlevé IV equation following the ideas outlined by Miwa, Jimbo and Ueno [30,31]. In the sequence of this, we will also discuss the associated Stokes' phenomenon and the monodromy problem (see Wasow [42]).
The second part of this section will be devoted to applying the general setting of Painlevé IV in a particular case, in order to obtain a special solution for it. This will correspond to a special Riemann-Hilbert problem whose properties will be outlined. The reason and motivation for introducing this solution is the fact that this Riemann-Hilbert problem will be seen to be useful in the following section.
The third and final part of this section will be devoted to the study of the Fredholm determinant, as it has been established by Its [28], Harnad [25], Bertola and Cafasso [8], among others. This will be used in order to study the τ -function of the special solution of Painlevé IV that we introduced before and analyze its validity with respect to the value of the parameter s.

The general Painlevé IV
We recall the Lax pair formulation of the Painlevé IV equation as in the original work of Miwa, Jimbo and Ueno [30,31], with notation adapted to our purposes.
Consider the matrix-valued function Ψ, associated with Painlevé IV, as the fundamental joint solution of the system of linear differential equations, the Lax pair, given by where σ 3 is the third Pauli matrix, Y , Z and U are functions of s whose dependence we now derive, while Θ 0 and Θ ∞ are arbitrarily chosen constants. The compatibility of (2.1) yields the zero-curvature equation which, in turn, yields the following system of nonlinear ODEs In the equations above the prime means differentiation with respect to s. By eliminating U , Z from the above equations we obtain the classical form of the fourth Painlevé equation

Stokes' phenomenon and the Monodromy problem
The review material of this section is mostly based on [21] and [30,31]. The solution of the first equation in (2.1) admits a formal sectorial expansion at λ = ∞ of the form Figure 3. Stokes' phenomenon and the corresponding jumps in the Ψ(λ) Riemann-Hilbert problem.
In each of the sectors q = (0, I, II, III, IV, V) of Fig. 3 there exists a unique solution Ψ (q) (λ) which is Poincaré asymptotic to Ψ formal (λ), even though the latter is only a non-convergent (in general) formal series. The matrices in Fig. 3 on the rays Γ q are the Stokes' matrices S q (and the formal monodromy on the negative axis) and relate the different solutions as follows: Ψ (q+1) = Ψ q S q (q modulo 6). The rays Γ 1,5 can be chosen arbitrarily within the respective quadrants, or more generally any smooth infinite arc within the quadrants that admits an asymptotic direction. Each of these solutions of Painlevé IV for each sector q, Ψ (q) (λ), has an asymptotic expansion near λ = 0 with the following behaviour where the matrices C (q) defined by this equation are called "connection matrices" and G is a diagonalizing matrix for the 1 λ coefficient of A(λ; s). From this equation and using an argument of analytic continuation, we obtain the following monodromy relation Figure 4. Jumps in the Ψ(λ) Riemann-Hilbert problem.
The remaining connection matrices C (q) are obtained from C (0) by multiplying it by an appropriate sequence of the Stokes' matrices; for example The general solution to the isomonodromic equations (2.2) for given Θ 0 , Θ ∞ is parametrized by the choices of parameters S −1 , . . . and connection matrix C (0) that satisfy (2.8).

Special Riemann-Hilbert problem Ψ(λ) and Painlevé IV
Several particular solutions of the Painlevé IV have been considered in the literature, see, e.g., [7,14,29]. In our work we only need special values for the Stokes' parameters and Θ 0,∞ : (2.10) The resulting Riemann-Hilbert problem is the one depicted in Fig. 4, where we have more appropriately renamed Γ 1 and Γ ∞ . It is formalized below: admits non-tangential boundary values. Moreover

1.
The boundary values are related by (see Fig. 4): where Ψ ± are the boundary values on the left and right of the oriented contour Σ Ψ .

2.
Near λ = ∞ the matrix has the following sectorial behaviour Figure 5. Jumps in the M (λ) Riemann-Hilbert problem.
3. Near λ = 0 in the region Ω ∞ , Ψ(λ; s) has the behaviour We remind the reader that subcripts ± denote the boundary values from the left (+) or right (−) of an oriented countour. We will henceforth omit explicit notation for the s dependence of Ψ.
We now proceed to do a further modification of this Ψ(λ)-Riemann-Hilbert problem, that will be used in the next section for the study of the Fredholm determinant: • move the jump contour Γ ∞ = iR in the following way iR → iR + ε, with > 0; • join together (collapse) the two parts of the contour Γ 1 and R − . The use of the verb "move" means that we define a new matrix-valued function in the region between the old and new position of the contours that analytically continues Ψ across the original contours. This is possible because the jump matrices are all constants in λ. We still use the symbol Ψ for the resulting function. Finally, we introduce M as follows The powers of λ are intended as the principal determination with arg(λ) ∈ (−π, π). The matrix M solves a Riemann-Hilbert problem with jumps indicated in Fig. 5 and formalized below.
2. Near λ = ∞ the matrix has the following sectorial behaviour:

Fredholm determinant and the τ -function of Painlevé IV
The Riemann-Hilbert Problem 2.2 is associated to an integral operator on L 2 (R − ∪ Γ ∞ ), where Γ ∞ = iR+ , and will be used in the next section in order to compute the Fredholm determinant. We first introduce some general definitions. Let (X, dν) be a measurable space and K(x, y) : X × X → C a function. We define the operatorK : If X 2 |K(x, y)| 2 dxdy < +∞, the operatorK is a Hilbert-Schmidt operator and its spectrum is purely discrete, the eigenvalues can accumulate only at 0 and all the multiplicities of the nonzero eigenvalues are finite (see, e.g., [39]).
Definition 2.3. The τ -function is defined as the following Fredholm determinant where ρ is a nonzero complex parameter.
From the properties of the Fredholm determinant we know that τ (ρ) = 0 if and only if 1 ρ is an eigenvalue ofK. With this being said, since in our case we are interested in a different parametric dependence ofK and not its spectral properties, we shall set ρ = 1 (see (2.13)).

Riemann-Hilbert problems and Fredholm determinants
Following [28] let us consider a set of contours Σ ⊂ C and f, g : Σ → C p matrix-valued (smooth) functions that satisfy the condition f T (λ)g(λ) = 0. Then we define the following scalar kernel We denote byK : Theorem 2.4 (see [28]). The kernel, R(λ, µ), of the resolvent operator R is given by where the p × p matrix A(λ) solves the Riemann-Hilbert problem Moreover, the above Riemann-Hilbert problem admits a unique solution if and only if the Fredholm determinant det Id −K is non-zero.
The case of interest to us, below, is with p = 2.

Fredholm determinant and the Riemann-Hilbert problem M
The Riemann-Hilbert Problem 2.2 is a special example of the previous setup, as we now see. We first define the following vectors where χ ∞ and χ R − are the characteristic functions of Γ ∞ and R − , respectively. This leads to the following jump matrices for Theorem 2.4 which match those in the Riemann-Hilbert Problem 2.2 (1 denotes the 2 × 2 identity matrix). Let Σ = Γ ∞ ∪ R − , as a consequence of this, we can now write the Riemann-Hilbert Problem (2.2) for the matrix M in the following way It should also be noted that, due to our construction, we have f T (z)g(z) ≡ 0. Let us consider the operatorK : L 2 (Σ) → L 2 (Σ) associated to the above Riemann-Hilbert problem. Due to the decomposition of L 2 (Σ) = L 2 (R − ∪Γ ∞ ) = L 2 (R − )⊕L 2 (Γ ∞ ) and ϕ = ϕ 0 ⊕ϕ 1 with ϕ 0 ∈ L 2 (R − ) and ϕ 1 ∈ L 2 (Γ ∞ ), we can write the operatorK in the form Due to this fact, the operatorK can now be written in block formK ij witĥ In our case, due to the form of the vectors f and g in (2.11), (2.12), the diagonal blocks are null operators:K 00 =K 11 = 0, whilê Following the Definition 2.3 of the tau function, we can now compute it for the case whenK is given by this operator We thus see that the determinant ofK can be written as the Fredholm determinant of the operator on L 2 (R − ) given byK 01K10 and we are thus left with the equality (2.14) We defineF :=K 01K10 : The kernel of this operator is is finally written as where we recall that θ(x) = x 2 /4 + sx/2 so that the operatorF =F(γ, s) depends on the parameters s and γ. We observe that det Id −K is non zero as long as the sup-norm of the operatorF is less than one. The values of s which satisfy this condition correspond to a regular solution of the Painlevé IV equation. We also observe that the kernel F (z, w) is real-valued (and symmetric) and hence the tau function is real-valued.

Norm estimate for small s
For s negative and sufficiently large in absolute value, we can estimate the norm of the operator F =K 01K10 and guarantee that the determinant is non-zero. Introducing the quantities we can write the operatorF in the form π which is a bounded operator from L 2 (R − ) to itself with norm one: |||P||| = 1. Here and below, we denote by Given an essentially bounded function f ∈ L ∞ (R − ), we denote by M f the corresponding multiplication operator and recall that |||M f ||| = f ∞ . This leads to As a result of this procedure, we obtain the following estimate for the norm ofF A calculus exercise shows that the function a(z) in (2.16), satisfies a ∞ ≤ 1, for z ∈ R − , γ ∈ (0, 1), s < 0. We now estimate the values of h ∞ . By setting x = ε + it, t ∈ R in the formula for h(z) in (2.16) and s < 0 and ε = max(0, −s) we have The estimate for F becomes At this point, we need to see for which values of s < 0 the norm ofF is smaller than 1, which guarantees that the determinant (2.14) is not zero and hence that our Riemann-Hilbert Problem 2.2 admits a solution. From (2.17) we have (recall γ ∈ (0, 1)) and we can easily see that the normF is less than one for s < s 0 where s 0 −0.7701449782.
In summary, we have proven It would be desirable to show that the solution of the Painlevé equation (2.3) for our choice of monodromy data (2.9), (2.10) is pole-free on the whole real line. However, the numerical analysis which will be performed in the following section, shows that there is a discrete (in principle infinite) number of values of s for which the Fredholm determinant (2.14) vanishes and hence the solution of (2.3) has poles.

Gaussian quadrature and numerics
The following is an exposition, for the benefit of the reader, of the main idea of numerical evaluation of Fredholm determinants based on Nystrom's method, as explained by Bornemann [10]. The simple idea is to suitably "discretize" the integral operator using appropriate Gaussian quadrature (see [40]).
More specifically, the type of quadrature formulae that we will use are the so-called "Gauss-Hermite" quadratures, which are defined in the following way Definition 2.6. Considering an integral of the form the Gauss-Hermite quadrature is an approximation of the value of this integral and is defined as where m is the number of points used in this approximation: the points x We will have to adapt this definition so that it is suitable to be applied to our case. In order to use this technique for the computation of F (x, y) in (2.15) we rewrite it as H(x, y), Then Nystrom's method states that the Fredholm determinant is approximated by where the nodes x , where H jk = e ε 2 2 +sε sin (πγ) In order to compute the above integral, we once again apply the Gauss-Hermite quadrature Note that, in the above, we use all the 2n nodes because the integral is over the whole R. In principle, we could use a different number of nodes for this last quadrature, but this is the way the code was actually implemented for simplicity.
where A is the n × 2n matrix given by and 1 ≤ j ≤ n and 1 ≤ ≤ 2n. The result of this approximation is expressed in terms of the parameters s ∈ R, γ ∈ [0, 1) and n ∈ N, which is the number of nodes in the Gaussian quadrature. Therefore, the final approximation that we have for τ is We take ε = max [0, −s] so as to optimize the distance from the saddle-point of the function θ(x) = x 2 /4 + sx/2 while respecting the condition that the relative interior of the contours R − and Γ ∞ do not intersect. It is worth pointing out that the value of the Fredholm determinant (numerical considerations aside) is independent of ε thanks to Cauchy's theorem. Various numerical results are pictures in Figs. 6, 7, 8.

Orthogonal polynomials and their Riemann-Hilbert problem
The orthogonality relations for the polynomials p n can be rewritten C p n (λ)λ jd+ e −N W (λ) dA(λ) = 0, j = 0, . . . , k − 1, As a consequence, the n-th monic orthogonal polynomial has a discrete symmetry    Using this definition, the initial sequence of orthogonal polynomials {p n (λ)} ∞ n=0 can be split in d sub-sequences, each of which labelled by the remainder ≡ n mod d. Through a change of coordinate λ d = u these sequences of monic orthogonal polynomials are seen to satisfy the orthogonality relations With the further change of coordinate u = −t(z − 1), z ∈ C, we introduce the transformed monic polynomial Next we transform the orthogonality on the plane to the orthogonality on a contour as in [4] and [3] Theorem 3.1 ([4, Theorem 2.1]). For any polynomial q(u) the following identity is satisfied where γ ∈ (0, 1), j is an arbitrary non-negative integer, andΣ is a positively oriented simple closed loop enclosing u = 0 and u = t. Consider now q ( ) k (u) to be the monic polynomial of degree k with the orthogonality relations (3.1). Then according to Theorem 3.1 the transformed monic polynomial π k (z) of degree k given by (3.2), is characterised by the non-hermitian orthogonality relations where Σ is a simple, positively oriented contour encircling z = 0 and z = 1, as it can be seen in Fig. 9, and the function z z−1 γ is analytic in C \ [0, 1] and tends to one for |z| → ∞.

The Riemann-Hilbert problem
We now consider the polynomials π k (z) in the limit k → ∞ and N → ∞ in such a way that, for n = kd + , one has T = n − N > 0.
We set the notation
When the limit k → ∞ is taken, three different regimes arise • pre-critical case: 0 < t < t c , leading to z 0 > 1, • critical case: t = t c , leading to z 0 = 1, • post-critical case: t > t c , leading to z 0 < 1.
The pre-and post-critical case were already analyzed in [4] and we are now interested in analyzing the critical case. We begin by defining the so-called complex moments as where the dependency on k is omitted in order to simplify notation, and use this to introduce the auxiliary polynomial It can be seen that this polynomial is not necessarily monic and its degree may be less than k −1.
In order to guarantee the existence of such a polynomial, the determinant in the denominator must not vanish.
We can now reformulate the condition of orthogonality for the polynomials π k (z) as a Riemann-Hilbert boundary value problem. To do this, we define the matrix which is the unique solution of following Riemann-Hilbert problem for orthogonal polynomials [22].
Riemann-Hilbert Problem 3.3. The matrix Y (z) is analytic in C\Σ, where Σ is the oriented curve in Fig. 10. The limits Y ± (z) exist and are continuous along Σ. Moreover

Jump on Σ:
The continuous boundary values Y ± (z) are such that This jump can be seen in Fig. 10.
2. Behaviour at infinity: Y (z) has the following behaviour as z → ∞ 4 Riemann-Hilbert analysis

Transforming the Riemann-Hilbert problem
We will now begin by doing a simple transformation through the first undressing step and, after that, use the steepest descent method of Deift-Zhou [16] to proceed with further transformations and simplifications of our problem.

First undressing step
The first undressing step consists in a simplification of the Riemann-Hilbert Problem 3.3 that is done by defining a new matrix-valued functionỸ (z) as Clearly the matrixỸ now satisfies a new Riemann-Hilbert problem whose simple formulation we leave to the reader. From the solution of Riemann-Hilbert problem we can recover the orthogonal polynomials through

TransformationỸ → R
We now define the so-called g-function; this is a scalar function that is analytic off a contour Γ r (see Fig. 13), homotopically equivalent to Σ in C \ [0, 1] and suitably chosen. Following [4] let where the logarithm is the principal determination and Γ r = z ∈ C : Re(ϕ(z; r)) = 0, |z| ≤ z 0 ,  Figure 11. Jumps for the Riemann-Hilbert problem of the matrix R(z). and ϕ(z; r) = or equivalently, Note that g is such that In order to establish the transformationỸ → R we will now deform the contour Σ to the contour Γ r . The choice of r will be explained later. For the time being we define the matrix R as The matrix R(z) solves the following Riemann-Hilbert problem Riemann-Hilbert Problem 4.1. The matrix R(z) is analytic in C \ (Γ r ∪ [0, 1]) and admits non-tangential boundary values that satisfy: 1. Jumps on Γ r and [0, 1] (see Fig. 11): e −γπiσ 3 , z ∈ (0, 1). (4.5) 2. Large z boundary behaviour: 3. Endpoint behaviour: Figure 12. Jumps in the T (z) Riemann-Hilbert problem.
The orhtogonal polynomials can be recovered from this Riemann-Hilbert problem through (4.6)

Transformation R → T
The jump matrix of R(z) on z ∈ Γ given by (4.5) can be factorized in the following way where, in order to write the second line, we have used the definition of the ϕ(z) function defined in (4.3) and the function g(z) defined in (4.1) as well as the relation (4.4). We will now consider three different loops Γ i , Γ r and Γ e , so that the space is split into four different domains Ω 0 , Ω 1 , Ω 2 and Ω ∞ , as it can be seen in Fig. 12. Γ i is in the interior of Γ r and Γ e is in the exterior. Using this, we will now define a new matrix-valued function T (z) in the following way 1 Figure 13. The contours Γ r for z 0 = 1 and different values of 0 < r ≤ 1.

Choice of contour
Since our analysis will be performed in the critical case when z 0 → 1, the choice of the contour Γ r with 0 < r ≤ z 0 = 1 is forced (Fig. 13). The corresponding contours Γ i and Γ e have to be chosen in such a way that the jump matrix v T defined in (4.8) converges (as k → ∞ and hence also as n → ∞) exponentially fast to a constant. We make the choice of the deformed contours Γ i and Γ e as in Fig. 14. This is the correct choice because the sign of Re ϕ(z) on the contours Γ i and Γ e is negative (see Fig. 15 where the region Re ϕ(z) < 0 is plotted).
Therefore, the jump matrices on Γ i and Γ e converge exponentially fast to the identity matrix in any L p norm except p = ∞ because ϕ(z = 1) = 0.
With this choice of contours we arrive to the following theorem (whose proof is a simple inspection and hence omitted)  Figure 14. Contours in the critical case z 0 = 1. Figure 15. Contours Γ chosen for the critical case z 0 = r = 1. The orange region corresponds to the part of the complex plane where Re (ϕ(z)) < 0 and the green contour is the region where Re (ϕ(z)) > 0.
where D is a small disk surrounding the point z = 1.

Approximate solutions to T (z)
We are now ready to approximate the matrix T (z) with two solutions, one outside a neighbourhood of z = 1 and one inside. We call exterior parametrix the Riemann-Hilbert problem solved by the matrix M (z) with jump v ∞ (z). We call local parametrix the solution P (z) of the Riemann-Hilbert problem obtained within a neighbourhood of z = 1. These two solutions are approximations of the exact solution, T (z), in the limit k → ∞. In order to obtain the asymptotics of the orthogonal polynomials π k (z), we need sub-leading corrections to the matrices N (z) and P (z) and this will be accomplished by evaluating perturbatively the error matrix E(z), which is defined as We will first construct the matrix N (z) and then the matrix P (z).

Exterior parametrix
The exterior parametrix, N (z), is the piecewise analytic 2 × 2 matrix where χ L is the characteristic function that is one on the left of the contourΓ r=1 and is zero otherwise. By design, it satisfies N (∞) = 1 and where v ∞ has been defined in (4.9).

Local parametrix and double scaling limit
The local parametrix near the point z = 1 is the solution, P (z), of a matrix Riemann-Hilbert problem that has the same jumps as T (z) and that matches N (z) on the boundary of a disc centered at z = 1.

Riemann-Hilbert Problem 4.4. The matrix P (z) is analytic in
and admits non-tangential boundary values. Moreover 1. Jumps on Σ P : which are the jumps that can be seen in Fig. 16. 2. Behaviour at the boundary ∂D: as k → ∞ and z ∈ ∂D. (4.10)

The double scaling limit
Let φ(z; z 0 ) be defined by We note that φ(z; z 0 ) is the analytic continuation of ϕ(z; z 0 ) (4.3) (setting r = 1) from the exterior of Γ = Γ 1 . For z 0 ∼ 1 we see that it has a single critical point in the neighbourhood of z = 1; then Figure 16. Jumps of the Riemann-Hilbert problem for the matrix P (z).
Remark 4.6. A simple computation yields the following first few terms for the expansion (δ = z 0 − 1)

From this it follows
which manifests the fact that √ kζ homothetically expands the image of the disk D.
Definition 4.7 (double-scaling limit). The double-scaling limit is defined by taking k → ∞ and A → 0 (or equivalently z 0 → 1) so that lim k→∞, z 0 →1 √ kA = S with S in compact subsets of the complex plane. Equivalently, it can be defined as (4.12)

Model problem and the Painlevé IV equation
We will now recall the Riemann-Hilbert Problem 2.1 for the function Ψ(λ; s) associated to the Painlevé IV equation. For convenience, we replace the contour Γ 1 withΓ i and the contour Γ ∞ withΓ e in order to be consistent with the notation we are using in this section. Given that Θ ∞ = Θ 0 = γ 2 , the equations (2.2) and (2.7) reduce to a somewhat simpler form where Y = Y (s) solves the Painlevé IV equation (2.3) in the form For our purposes, we need to modify the Riemann-Hilbert problem for Ψ(λ) to be able to match the Riemann-Hilbert problem for P (z). Let us introduce the contourΓ r=1 on the λ plane, which is a contour betweenΓ i andΓ e and passing through λ = 0. Then, let us define where θ is given by (2.6) and χ L is the characteristic function that is one on the left of the con-tourΓ r=1 and is zero otherwise. It is now straightforward to check that the Riemann-Hilbert problem for Ψ(λ) is given by Riemann-Hilbert Problem 4.8. The matrix Ψ(λ)is analytic in C\Σ Ψ , Σ Ψ =Γ e ∪Γ i ∪Γ r=1 ∪ R − and admits non-tangential boundary values. Moreover

1.
Jumps on Σ Ψ : (4.14) 2. Behaviour for λ → ∞ Comparing the jump matrices for P (z) and Ψ(λ), we are now ready to obtain the local parametrix P (z) defined by the Riemann-Hilbert Problem 4.4, which is given by where ζ(z) = ζ(z; z 0 ) as defined in (4.5). We observe that the product of the first three terms of (4.16) is holomorphic in the neighbourhood of z = 1 and therefore does not change the Riemann-Hilbert problem. Furthermore, the first three terms have been inserted in order to have the behaviour (4.10) for z ∈ ∂D in the limit k → ∞ and A → 0 in such a way that In the analysis above and below, we assume that S belongs to compact sets where the solution of the Painlevé IV equation does not have poles. Using (4.15), (2.4) and (2.5), we obtain in the limit k → ∞ (4.17) From the above expansion we can see that the subleading terms are not uniformly small with respect to γ as k → ∞ since γ ∈ [0, 1). Note that the source of the non-uniformity in the error analysis is the element (2, 1) in (4.17). For this reason we need to introduce an improved parametrix in the next section.

Improved parametrix
To construct an improved parametrix we define where B is to be determined. In the same way we define where U is the (2, 1) entry of the subleading term of the expansion of Ψ(λ) for λ → ∞. Now we have in the limit k → ∞ P (z) = N (z)   1 + 0 −1 1 0 Painlevé IV Critical Asymptotics for Orthogonal Polynomials in the Complex Plane

29
The improved parametrics N (z) and P (z) have the same jump discontinuities as before, but P (z) might have poles at z = 1. The constant B in (4.18) is then determined by the requirement that P (z) is bounded at z = 1. This gives the constant B as We are now ready to compute the error matrix E(z).

Error matrix
The error matrix E(z) is defined as where the boundary of D is oriented clockwise. The matrix has also jumps on Γ i ∪ Γ e \ D which are exponentially close to the identity matrix and can be ignored for the purposes of the analysis. Then, the matrix E(z) satisfies the Riemann-Hilbert problem In the double scaling limit k → ∞, after using (4.14), (4.15), (2.4), (2.5) and the last equation in (4.10), the jump matrix v E (z) takes the form where the expansion is for z ∈ ∂D (while v E is exponentially close to the identity on Γ E \ D), and v (1) where σ + = ( 0 1 0 0 ). By the standard theory of small norm Riemann-Hilbert problems (see for example [17,Chapter 7]), one has a similar expansion for E(z), namely so that E (j) E (z), z ∈ ∂D, j = 1, 2, 3.
By solving the corresponding Riemann-Hilbert problem, we obtain, using the Plemelj-Sokhtski formula which gives   where we recall from (4.12) that z 0 = √ k √ k+S and that the function Z = Z(S) and U = U (S) are related to the Painlevé IV equation via (4.13). The curves in (4.26) approach Γ r=1 at the rate O(log k/k) and lie in Int(Γ r=1 ). The normalized counting measure of the zeros of π k (z) converges to the probability measure ν defined in (1.9).
The proof of Proposition 1.4 follows immediately from the proof of Proposition 4.9. The rest of the proof of Theorem 1.2 follows the steps obtained in [4] and for this reason we omit it.