111The research of J. Wei is partially supported by NSERC of Canada. Lei Zhang
is partially supported by a Simons Foundation Collaboration Grant

# Vanishing Estimates for Liouville equation with quantized singularities

Juncheng Wei  Department of Mathematics
University of British Columbia
Vancouver, BC V6T1Z2, Canada<EMAIL_ADDRESS>and  Lei Zhang  Department of
Mathematics
University of Florida
1400 Stadium Rd
Gainesville FL 32611<EMAIL_ADDRESS>
###### Abstract.

In this article we continue with the research initiated in our previous work
on singular Liouville equations with quantized singularity. The main goal of
this article is to prove that as long as the bubbling solutions violate the
spherical Harnack inequality near a singular source, the first derivatives of
coefficient functions must tend to zero.

###### 1991 Mathematics Subject Classification:

35J75,35J61

##  1\. Introduction

In this article we study bubbling solutions of

(1.1) $\Delta u+\mathrm{H}(x)e^{u}=4\pi\alpha\delta_{0}\quad\mbox{ in
}\quad\Omega\subset\mathbb{R}^{2}$

where $\Omega$ is an open, bounded subset of $\mathbb{R}^{2}$ that contains
the origin, $\alpha>-1$ is a constant and $\delta_{0}$ is the Dirac mass at
$0$, $\mathrm{H}$ is a positive and smooth function. Since

$\Delta(\frac{1}{2\pi}\log|x|)=\delta_{0},$

one can use the logarithmic function $2\pi\alpha\log|x|$ to remove the
singular source from the equation: let $u_{1}(x)=u(x)-2\alpha\log|x|$, then
$u_{1}$ satisfies

(1.2) $\Delta
u_{1}+|x|^{2\alpha}\mathrm{H}(x)e^{u_{1}}=0,\quad\mbox{in}\quad\Omega.$

If a sequence of solutions $\\{u^{k}\\}_{k=1}^{\infty}$ of (1.2) satisfies

$\lim_{k\to\infty}u^{k}(x_{k})=\infty,\quad\mbox{ for some $\bar{x}\in
B_{\tau}$ and $x_{k}\to\bar{x}$.}$

we say $u^{k}$ is a sequence of bubbling solutions or blowup solutions,
$\bar{x}$ is called a blowup point. For many reasons in applications it is
most interesting to consider $\alpha=N\in\mathbb{N}$ (the set of natural
numbers) and when $0$ is the only blowup point of $u^{k}$. Our set-up of
bubbling solutions is as follows: Let $\mathfrak{u}_{k}$ be a sequence of
solutions of

(1.3)
$\Delta\mathfrak{u}_{k}(x)+|x|^{2N}\mathrm{H}_{k}(x)e^{\mathfrak{u}_{k}}=0,\quad\mbox{in}\quad
B_{\tau}$

for some $\tau>0$ independent of $k$. $B_{\tau}$ is the ball centered at the
origin with radius $\tau$. In addition we postulate the usual assumptions on
$\mathfrak{u}_{k}$ and $\mathrm{H}_{k}$: For a positive constant $C$
independent of $k$, the following holds:

(1.4)
$\left\\{\begin{array}[]{ll}\|\mathrm{H}_{k}\|_{C^{3}(\bar{B}_{\tau})}\leq
C,\quad\frac{1}{C}\leq\mathrm{H}_{k}(x)\leq C,\quad x\in\bar{B}_{\tau},\\\ \\\
\int_{B_{\tau}}\mathrm{H}_{k}e^{\mathfrak{u}_{k}}\leq C,\\\ \\\
|\mathfrak{u}_{k}(x)-\mathfrak{u}_{k}(y)|\leq C,\quad\forall x,y\in\partial
B_{\tau},\end{array}\right.$

and since we study the asymptotic behavior of blowup solutions around the
singular source, we assume that there is no blowup point except at the origin:

(1.5) $\max_{K\subset\subset B_{\tau}\setminus\\{0\\}}\mathfrak{u}_{k}\leq
C(K).$

Also, we use the value of $\mathfrak{u}_{k}$ on $\partial B_{\tau}$ to define
a harmonic function $\phi_{k}(x)$:

(1.6) $\left\\{\begin{array}[]{ll}\Delta\phi_{k}(x)=0,\quad\mbox{in}\quad
B_{\tau},\\\ \\\
\phi_{k}(x)=\mathfrak{u}_{k}(x)-\frac{1}{2\pi\tau}\int_{\partial
B_{\tau}}\mathfrak{u}_{k}dS,\quad x\in\partial B_{\tau}.\end{array}\right.$

Clearly the mean value property of harmonic functions implies $\phi_{k}(0)=0$
and the finite oscillation of $\mathfrak{u}_{k}$ on $\partial B_{\tau}$ means
that all derivatives of $\phi_{k}$ are uniformly bounded in $B_{\tau/2}$. In
this article we consider the case that:

(1.7) $\max_{x\in B_{\tau}}\mathfrak{u}_{k}(x)+2(1+N)\log|x|\to\infty,$

which is equivalent to saying that the spherical Harnack inequality does not
hold for $\mathfrak{u}_{k}$. It is also mentioned in literature ( see [14, 19]
) that $0$ is called an non-simple blowup point. The main result of this
article is

######  Theorem 1.1.

Let $\\{\mathfrak{u}_{k}\\}$ be a sequence of solutions of (1.3) such that
(1.4),(1.5) and (1.7) hold. Then

$\nabla(\log\mathrm{H}_{k}+\phi_{k})(0)=o(1),\quad\mbox{as }k\to\infty,$

where $\phi_{k}$ is defined in (1.6).

When bubbling solutions satisfy (1.7), they are called non-simple blowup
solutions ( see [14]). Theorem 1.1 is a complement of Theorem 1.1 of [19],
which asserts that under certain conditions (see Theorem A below)
$\nabla(\log\mathrm{H}_{k}+\phi_{i}^{k})(0)$ tends to zero. Theorem 1.1
removes the restrictions in [19]. In other words, the combination of Theorem A
and Theorem 1.1 proves that $\nabla(\log\mathrm{H}_{k}+\phi_{k})(0)\to 0$ as
long as the non-simple blowup situation occurs. Besides the advancement of
analytical understanding, this conclusion is particularly important in
application. Theorem 1.1 can be applied to situations beyond single equations.
For certain systems of equations such as Toda systems, the bubble
accumulations can be described by a sequence of bubbling solutions with
quantized singular source. Theorem 1.1 is very useful to rule out complicated
bubbling accumulation pictures in Toda systems.

It remains an open question whether or not
$\nabla(\log\mathrm{H}_{k}+\phi_{k})(0)$ tends to $0$ if the bubbling
solutions satisfy the spherical Harnack inequality around the origin. We tend
to believe one can construct a sequence of bubbling solutions that satisfy
spherical Harnack inequality with nonzero first derivatives of the
coefficients functions. In particular the works of Del Pino-Esposito-Musso [8,
9] on two dimensional Euler flows seem to suggest that for simple blowup
solutions with quantized singular sources, the first derivatives of the
coefficient functions may not tend to zero at singular sources. Instead the
$(N+1)$-th derivatives should vanish. Our result, together with [8, 9],
demonstrates a striking contrast between simple and non-simple bubblings.

The non-simple bubbling situation and vanishing theorems have profound impact
to problems in geometry and physics. For example for the following mean field
equation defined on a Riemann surface $(M,g)$:

(1.8)
$\Delta_{g}u+\rho(\frac{h(x)e^{u(x)}}{\int_{M}he^{u}}-\frac{1}{Vol_{g}(M)})=4\pi\sum_{j}\alpha_{j}(\delta_{p_{j}}-\frac{1}{Vol_{g}(M)}),$

the solution $u$ represents a conformal metric with prescribed conic
singularities (see [10, 17, 18]) . In particular if the singular source is
quantized, the Liouville equation has close ties with Algebraic geometry,
integrable system, number theory and complex Monge-Ampere equations (see [7]).
In Physics the understanding of non-simple blowup phenomenon would be
extremely useful for the study of mean field limits of point vortices in the
Euler flow [4, 5] and models in the Chern-Simons-Higgs theory [13] and in the
electroweak theory [1], etc. It is also remarkable that non-simple bubbling
solutions also occur in systems. In [12], the non-simple blowup solutions are
studied for singular Liouville systems. Finally we remark that when the blowup
point is a location of a singular source, whether or not this point has to be
a critical point of coefficient functions has intrigued people for years. Our
previous result [19] is the first result for singular Liouville equation, the
second author proved a surprising vanishing theorem for singular Toda systems
in [22].

The organization of this paper is as follows. In section two we review a few
fundamental tools for the proof of the main theorem and invoke several key
estimates established in our previous work [19]. Then in section three we use
a sequence of global solutions to approximate our blowup solutions. The point-
wise estimates proved in this section are more precise than what is
established in [19] and are important for our argument. In section four we
prove a crucial estimate on the difference between blowup solution and the
global solutions as the first term in the approximation. As a consequence of
section four, we move to section five to complete the proof of the main
theorem. The proof in section four is similar to the proof of uniqueness
theorems for bubbling solutions in [20], [2], [15], etc.

Notation: We will use $B(x_{0},r)$ to denote a ball centered at $x_{0}$ with
radius $r$. If $x_{0}$ is the origin we use $B_{r}$. $C$ represents a positive
constant that may change from place to place.

##  2\. Preliminary discussions

For simple notation we set

(2.1) $u_{k}(x)=\mathfrak{u}_{k}(x)-\phi_{k}(x),\quad\mbox{and}$ (2.2)
$h_{k}(x)=\mathrm{H}_{k}(x)e^{\phi_{k}(x)}$

to write the equation of $u_{k}$ as

(2.3) $\Delta u_{k}(x)+|x|^{2N}h_{k}(x)e^{u_{k}}=0,\quad\mbox{ in }\quad
B_{\tau}$

Without loss of generality we assume

(2.4) $\lim_{k\to\infty}h_{k}(0)=1.$

Obviously (1.7) is equivalent to

(2.5) $\max_{x\in B_{\tau}}u_{k}(x)+2(1+N)\log|x|\to\infty.$

It is well known [14, 3] that $u_{k}$ exhibits a non-simple blowup profile. It
is established in [14, 3] that there are $N+1$ local maximum points of
$u_{k}$: $p_{0}^{k}$,….,$p_{N}^{k}$ and they are evenly distributed on
$\mathbb{S}^{1}$ after scaling according to their magnitude: Suppose along a
subsequence

$\lim_{k\to\infty}p_{0}^{k}/|p_{0}^{k}|=e^{i\theta_{0}},$

then

$\lim_{k\to\infty}\frac{p_{l}^{k}}{|p_{0}^{k}|}=e^{i(\theta_{0}+\frac{2\pi
l}{N+1})},\quad l=1,...,N.$

For many reasons it is convenient to denote $|p_{0}^{k}|$ as $\delta_{k}$ and
define $\mu_{k}$ as follows:

(2.6) $\delta_{k}=|p_{0}^{k}|\quad\mbox{and
}\quad\mu_{k}=u_{k}(p_{0}^{k})+2(1+N)\log\delta_{k}.$

Since $p_{l}^{k}$’s are evenly distributed around $\partial B_{\delta_{k}}$,
standard results for Liouville equations around a regular blowup point can be
applied to have $u_{k}(p_{l}^{k})=u_{k}(p_{0}^{k})+o(1)$. Also, (1.7) gives
$\mu_{k}\to\infty$. The interested readers may look into [14, 3] for more
detailed information.

In our previous work [19] we prove the following vanishing type estimates for
the first derivatives of the coefficient function $\log h_{k}$:

_Theorem A: Let $u_{k}$, $\phi_{k}$, $h_{k}$, $\delta_{k}$, $\mu_{k}$ be
defined by (2.3), (1.6), (2.2), (2.6) respectively. Then_

(2.7) $|\nabla\log
h_{k}(0)|=O(\delta_{k})+O(\delta_{k}^{-1}e^{-\mu_{k}}\mu_{k}).$

Here we observe that if $\mu_{k}e^{-\mu_{k}}=o(\delta_{k})$, we already have
$\nabla h_{k}(0)=o(1)$, which is equivalent to
$\nabla(\mathrm{H}_{k}e^{\phi_{k}})(0)=o(1)$. Thus throughout the paper we
assume

(2.8) $\delta_{k}\leq C\mu_{k}e^{-\mu_{k}}.$

Finally we shall use $E$ to denote a frequently appearing error term of the
size $O(\delta_{k}^{2})+O(\mu_{k}e^{-\mu_{k}})$. Because of (2.8),

$E=O(\mu_{k}e^{-\mu_{k}}).$

##  3\. Approximating bubbling solutions by global solutions

First we recall that $|p_{0}^{k}|=\delta_{k}$, so we write $p_{0}^{k}$ as
$p_{0}^{k}=\delta_{k}e^{i\theta_{k}}$ and define $v_{k}$ as

(3.1)
$v_{k}(y)=u_{k}(\delta_{k}ye^{i\theta_{k}})+2(N+1)\log\delta_{k},\quad|y|<\tau\delta_{k}^{-1}.$

If we write out each component, (3.1) is

$v_{k}(y_{1},y_{2})=u_{k}(\delta_{k}(y_{1}\cos\theta_{k}-y_{2}\sin\theta_{k}),\delta_{k}(y_{1}\sin\theta_{k}+y_{2}\cos\theta_{k}))+2(1+N)\log\delta_{k}.$

Then it is standard to verify that $v_{k}$ solves

(3.2) $\Delta
v_{k}(y)+|y|^{2N}\mathfrak{h}_{k}(\delta_{k}y)e^{v_{k}(y)}=0,\quad|y|<\tau/\delta_{k},$

where

(3.3) $\mathfrak{h}_{k}(x)=h_{k}(xe^{i\theta_{k}}),\quad|x|<\tau.$

Thus the image of $p_{0}^{k}$ after scaling is $Q_{1}^{k}=e_{1}=(1,0)$. Let
$Q_{1}^{k}$, $Q_{2}^{k}$,…,$Q_{N}^{k}$ be the images of $p_{i}^{k}$
$(i=1,...,N)$ after the scaling:

$Q_{l}^{k}=\frac{p_{l}^{k}}{\delta_{k}}e^{-i\theta_{k}},\quad l=1,...,N.$

It is established by Kuo-Lin in [14] and independently by Bartolucci-
Tarantello in [3] that

(3.4)
$\lim_{k\to\infty}Q_{l}^{k}=\lim_{k\to\infty}p_{l}^{k}/\delta_{k}=e^{\frac{2l\pi
i}{N+1}},\quad l=0,....,N.$

Then in our previous work [19] we obtained ( see (3.13) in [19])

(3.5) $Q_{l}^{k}-e^{\frac{2\pi li}{N+1}}=E.$

Choosing $3\epsilon>0$ small and independent of $k$, we can make disks
centered at $Q_{l}^{k}$ with radius $3\epsilon$ (denoted as
$B(Q_{l}^{k},3\epsilon)$) mutually disjoint. Let

(3.6) $\mu_{k}=\max_{B(Q_{0}^{k},\epsilon)}v_{k}.$

Since $Q_{l}^{k}$ are evenly distributed around $\partial B_{1}$, it is easy
to use standard estimates for single Liouville equations ([21, 11, 6]) to
obtain

$\max_{B(Q_{l}^{k},\epsilon)}v_{k}=\mu_{k}+o(1),\quad l=1,...,N.$

Let

(3.7)
$V_{k}(x)=\log\frac{e^{\mu_{k}}}{(1+\frac{e^{\mu_{k}}\mathfrak{h}_{k}(\delta_{k}e_{1})}{8(1+N)^{2}}|y^{N+1}-e_{1}|^{2})^{2}}.$

Clearly $V_{k}$ is a solution of

(3.8) $\Delta
V_{k}+\mathfrak{h}_{k}(\delta_{k}e_{1})|y|^{2N}e^{V_{k}}=0,\quad\mbox{in}\quad\mathbb{R}^{2},\quad
V_{k}(e_{1})=\mu_{k}.$

This expression is based on the classification theorem of Prajapat-Tarantello
[16].

The estimate of $v_{k}(x)-V_{k}(x)$ is important for the main theorem of this
article. For convenience we use

$\beta_{l}=\frac{2\pi
l}{N+1},\quad\mbox{so}\,\,e_{1}=e^{i\beta_{0}}=Q_{0}^{k},\quad
e^{i\beta_{l}}=Q_{l}^{k}+E,\,\,\mbox{ for }\,\,l=1,...,N.$

######  Proposition 3.1.

Let $l=0,...,N$ and $\delta$ be small so that $B(e^{i\beta_{l}},\delta)\cap
B(e^{i\beta_{s}},\delta)=\emptyset$ for $l\neq s$. In each
$B(e^{i\beta_{l}},\delta)$

(3.9)
$|v_{k}(x)-V_{k}(x)|\leq\left\\{\begin{array}[]{ll}C\mu_{k}e^{-\mu_{k}/2},\quad|x-e^{i\beta_{l}}|\leq
Ce^{-\mu_{k}/2},\\\ \\\
C\frac{\mu_{k}e^{-\mu_{k}}}{|x-e^{i\beta_{l}}|}+O(\mu_{k}^{2}e^{-\mu_{k}}),\quad
Ce^{-\mu_{k}/2}\leq|x-e^{i\beta_{l}}|\leq\delta.\end{array}\right.$

######  Remark 3.1.

Once (3.9) is established. We shall use a re-scaled version of Proposition
3.1. Let $\epsilon_{k}=e^{-\frac{1}{2}\mu_{k}}$, we have

(3.10)
$|v_{k}(e^{i\beta_{l}}+\epsilon_{k}y)-V_{k}(e^{i\beta_{l}}+\epsilon_{k}y)|\leq
C\mu_{k}^{2}\epsilon_{k}(1+|y|)^{-1},\quad 0<|y|<\delta_{0}\epsilon_{k}^{-1}.$

Proof of Proposition 3.1: The main idea of the proof is as follows. First from
the Green’s representation of $v_{k}$ we obtain a rather precise estimate of
$v_{k}$ in $B_{3}$ away from bubbling disks. On the other hand around each
$Q_{m}^{k}$ we invoke a standard pointwise estimate in [6, 21, 11] for
Liouville equation around a blowup point, which provides a precise description
of $v_{k}$ in a neighborhood of $Q_{m}^{k}$. The comparison of these two
estimates gives an accurate estimate of the maximum of $v_{k}$ around each
local maximum point.

Fixing the neighborhood of one $Q_{m}^{k}$, we first cite a result of Gluck
[11] (Appendix B of [19]) to write $v_{k}$ as

(3.11)
$v_{k}(y)=\log\frac{e^{\mu_{k,m}}}{(1+e^{\mu_{k,m}}\frac{|\tilde{Q}_{m}^{k}|^{2N}\mathfrak{h}_{k}(\delta_{k}\tilde{Q}_{m}^{k})}{8}|y-\tilde{Q}_{m}^{k}|^{2})^{2}}+\phi_{m}^{k}(y)+O(\mu_{k}^{2}e^{-\mu_{k}})$

where $\mu_{k,m}=v_{k}(\tilde{Q}_{m}^{k})$, $\phi_{m}^{k}$ is the harmonic
function taking $0$ at $Q_{m}^{k}$ that makes $v_{k}-\phi_{m}^{k}=constant$ on
$\partial B(Q_{m}^{k},\delta)$. $\tilde{Q}_{m}^{k}$ is where
$v_{k}-\phi_{m}^{k}$ takes its local maximum in a neighborhood of $Q_{m}^{k}$.
The difference between $\tilde{Q}_{m}^{k}$ and $Q_{m}^{k}$ is
$O(e^{-\mu_{k}})$. First we claim that

(3.12) $\mu_{k,m}-\mu_{k}=E.$

Let $\Omega_{k}=B(0,\tau\delta_{k}^{-1})$ and $G_{k}(y,\eta)$ be the Green’s
function on $\Omega_{k}$ with respect to Dirichlet boundary condition:

(3.13)
$G_{k}(y,\eta)=-\frac{1}{2\pi}\log|y-\eta|+\frac{1}{2\pi}\log(\frac{|\eta|}{\tau\delta_{k}^{-1}}|\frac{\tau^{2}\delta_{k}^{-2}\eta}{|\eta|^{2}}-y|),$

then for $y$ away from bubbling areas and $|y|\sim 1$, we have

$\displaystyle v_{k}(y)$
$\displaystyle=v_{k}|_{\partial\Omega_{k}}+\int_{\Omega_{k}}G_{k}(y,\eta)\mathfrak{h}_{k}(\eta)|\eta|^{2N}e^{v_{k}}d\eta,$
$\displaystyle=v_{k}|_{\partial\Omega_{k}}+\sum_{l=0}^{N}G(y,Q_{l}^{k})\int_{B(Q_{l}^{k},\epsilon)}|\eta|^{2N}\mathfrak{h}_{k}(\delta_{k}\eta)e^{v_{k}}d\eta$
$\displaystyle+\sum_{l}\int_{B(Q_{l}^{k},\epsilon)}(G_{k}(y,\eta)-G_{k}(y,Q_{l}^{k}))|\eta|^{2N}\mathfrak{h}_{k}(\delta_{k}\eta)e^{v_{k}}d\eta+E,$
$\displaystyle=v_{k}|_{\partial\Omega_{k}}+8\pi\sum_{l}G_{k}(y,Q_{l}^{k})+E.$

Later we shall use $H_{k}$ to denote the second term in (3.13). Note that we
use two standard estimates. First the integration outside bubbling disk is $E$
because

$v_{k}(x)\leq-\mu_{k}-(4(N+1)-o(1))\log|x|+C,\quad 3<|x|<\tau\delta_{k}^{-1}.$

Second, in the evaluation of the integral terms above we use standard bubble
expansion formula (see Gluck [11], for example) and symmetry properties. This
part is mentioned in Lemma 2.1 and Appendix B of [19]. It is important to
point out that the second estimate does not depend on $m$. In particular if we
consider $y$ located at $|y-Q_{m}^{k}|=\epsilon$, the expression of $v_{k}$
can be written as

(3.14) $\displaystyle v_{k}(y)$
$\displaystyle=v_{k}|_{\partial\Omega_{k}}-4\log|y-Q_{m}^{k}|+\phi_{m}^{k}$
$\displaystyle-4\sum_{l=0,l\neq
m}^{N}\log|Q_{m}^{k}-Q_{l}^{k}|+8\pi\sum_{l=0}^{N}H_{k}(Q_{m}^{k},Q_{l}^{k})+E,$

where

(3.15) $\phi_{m}^{k}=\sum_{l=0,l\neq
m}^{N}(-4)\log\frac{|y-Q_{l}^{k}|}{|Q_{m}^{k}-Q_{l}^{k}|}+8\pi\sum_{l=0}^{N}(H_{k}(y,Q_{l}^{k})-H_{k}(Q_{m}^{k},Q_{l}^{k}))$

is the harmonic function that takes $0$ at $Q_{m}^{k}$ and eliminates the
oscillation of $v_{k}$ on $\partial B(Q_{m}^{k},\epsilon)$. On the other hand
from (3.11) we have

(3.16)
$v_{k}(y)=-\mu_{k,m}-2\log\frac{|\tilde{Q}_{m}^{k}|^{2N}\mathfrak{h}_{k}(\delta_{k}\tilde{Q}_{m}^{k})}{8}-4\log|y-Q_{m}^{k}|+\phi_{m}^{k}+O(\mu_{k}^{2}e^{-\mu_{k}}).$

Comparing (3.16) and (3.14) on $|y-Q_{m}^{k}|=\epsilon$ we have

(3.17)
$\displaystyle-\mu_{m,k}-2\log\frac{|\tilde{Q}_{m}^{k}|^{2N}\mathfrak{h}_{k}(\delta_{k}\tilde{Q}_{m}^{k})}{8}$
$\displaystyle=$ $\displaystyle-4\sum_{l=0,l\neq
m}^{N}\log|Q_{m}^{k}-Q_{l}^{k}|+8\pi\sum_{l=0}^{N}H_{k}(Q_{m}^{k},Q_{l}^{k})+v_{k}|_{\partial\Omega_{k}}+O(\mu_{k}^{2}e^{-\mu_{k}}).$

To evaluate terms in (3.17) we observe that (see (3.5))

$|\tilde{Q}_{m}^{k}|^{2N}=1+E,\qquad\mathfrak{h}_{k}(\delta_{k}\tilde{Q}_{m}^{k})=1+E,$
$Q_{m}^{k}=e^{i\beta_{m}}+E,\qquad\tilde{Q}_{m}^{k}=Q_{m}^{k}+O(e^{-\mu_{k}}),$

and by the expression of $H_{k}(y,\eta)$:

$\displaystyle
H_{k}(y,\eta)=\frac{1}{2\pi}\log(\frac{|\eta|}{\tau\delta_{k}^{-1}}|\frac{\tau^{2}\delta_{k}^{-2}\eta}{|\eta|^{2}}-y|)$
$\displaystyle=\frac{1}{2\pi}\log(\tau\delta_{k}^{-1})+\frac{1}{2\pi}\log|\frac{\eta}{|\eta|}-\frac{|\eta|}{\tau^{2}}\delta_{k}^{2}y|$

we have

$H_{k}(Q_{m}^{k},Q_{l}^{k})=\frac{1}{2\pi}\log(\tau\delta_{k}^{-1})+E.$

Thus two terms in (3.17) are

(3.18)
$8\pi\sum_{l=0}^{N}H_{k}(Q_{m}^{k},Q_{l}^{k})=4(N+1)\log(\tau\delta_{k}^{-1})+E$
(3.19) $\sum_{l=0,l\neq m}^{N}\log|Q_{m}^{k}-Q_{l}^{k}|=\sum_{l=0,l\neq
m}^{N}\log|e^{i\beta_{m}}-e^{i\beta_{l}}|+E=\log(N+1)+E.$

The last equality can be verified by direct computation. Using (3.18) and
(3.19) in (3.17) we have

(3.20) $\displaystyle v_{k}|_{\partial\Omega_{k}}=$
$\displaystyle-\mu_{m,k}-2\log\frac{\mathfrak{h}_{k}(\delta_{k}e_{1})}{8}+4\log(1+N)-4(1+N)\log(\tau\delta_{k}^{-1})$
$\displaystyle+O(\mu_{k}^{2}e^{-\mu_{k}}),\quad m=0,1,...,N.$

The value $v_{k}|_{\partial\Omega_{k}}$ is independent of $m$. In particular
$\mu_{0,k}=\mu_{k}$. Thus the comparison of $\mu_{m,k}$ in (3.20) proves
(3.12). Next we observe that around $Q_{l}^{k}$

(3.21)
$V_{k}(y)=\log\frac{e^{\mu_{k}}}{(1+\frac{|\tilde{Q}_{l}^{k}|^{2N}\mathfrak{h}_{k}(\delta_{k}e_{1})e^{\mu_{k}}}{8}|y-\tilde{Q}_{l}^{k}|^{2})^{2}}+\tilde{\phi}_{l}^{k}(y)+O(\mu_{k}^{2}e^{-\mu_{k}}),$

where $y\in B(e^{i\beta_{l}},\delta_{0})$,
$\tilde{Q}_{l}^{k}=e^{i\beta_{l}}+O(e^{-\mu_{k}})$,

(3.22) $\tilde{\phi}_{l}^{k}(x)=\sum_{m=0,m\neq
l}^{N}(-4)\log\frac{|y-e^{i\beta_{m}|}}{|e^{i\beta_{m}}-e^{i\beta_{l}}|},$

$\delta_{0}$ is a small positive number independent of $k$. The way to prove
(3.21), by direct computation from the expression of $V_{k}$, is as follows:
It is easy to see that

$V_{k}(y)=-\frac{1}{2\pi}\int_{\mathbb{R}^{2}}\log|y-\eta|\mathfrak{h}_{k}(\delta_{k}e_{1})|\eta|^{2N}e^{V_{k}(\eta)}d\eta+C,\quad
y\in\mathbb{R}^{2}.$

Then $\mathfrak{h}_{k}(\delta_{k}e_{1})e^{V_{k}}|y|^{2N}$ weakly converges to
$8\pi\delta_{e^{i\beta_{l}}}$ in a small neighborhood of $e^{i\beta_{l}}$. For
$y\in\partial B(e^{i\beta_{l}},\delta)$ we have

$V_{k}(y)=-\sum_{l=0}^{N}4\log|y-e^{i\beta_{l}}|+C_{k}+O(\mu_{k}e^{-\mu_{k}}).$

From there we know that harmonic function around $e^{i\beta_{l}}$ that equals
$0$ at $e^{i\beta_{l}}$ is $\phi_{l}^{k}$ in (3.22). On the other hand the
equation of $v_{k}$ around $e^{i\beta_{l}}$ is (3.8). The standard expansion
(see [11]) for blowup solution leads to (3.21).

Since $Q_{m}^{k}-e^{i\beta_{m}}=E$, we can replace $\tilde{\phi}_{l}^{k}$ by
$\phi_{l}^{k}$ and have

$V_{k}(x)=\log\frac{e^{\mu_{k}}}{(1+\frac{\mathfrak{h}_{k}(\delta_{k}e_{1})e^{\mu_{k}}}{8}|y-e^{i\beta_{l}}|^{2})^{2}}+\phi_{l}^{k}+O(\mu_{k}^{2}e^{-\mu_{k}})$

in $B(Q_{l}^{k},e^{-\mu_{k}/2})$. Thus in the region
$B(Q_{l}^{k},e^{-\mu_{k}/2})$, the comparison between $v_{k}$ and $V_{k}$
boils down to the evaluation of:

(3.23)
$\log\frac{e^{\mu_{l,k}}}{(1+\frac{\mathfrak{h}_{k}(\delta_{k}e_{1})e^{\mu_{l,k}}}{8}|y-e^{i\beta_{l}}-p_{k}|^{2})^{2}}-\log\frac{e^{\mu_{k}}}{(1+\frac{\mathfrak{h}_{k}(\delta_{k}e_{1})e^{\mu_{k}}}{8}|y-e^{i\beta_{l}}|^{2})^{2}},$

for $|p_{k}|=E$. By elementary computation we see that the difference between
the two terms in (3.23) is $O(\mu_{k}e^{-\mu_{k}/2})$ if
$|y-e^{i\beta_{l}}|\leq Ce^{-\mu_{k}/2}$. On the other hand, for
$Ce^{-\mu_{k}/2}<|y-e^{i\beta_{l}}|<\epsilon/2$, the comparison of expressions
of $v_{k}$ and $U_{k}$ gives the difference upper bound as

$O(e^{-\mu_{k}})|y-e^{i\beta_{l}}|^{-1}+O(\mu_{k}^{2}e^{-\mu_{k}}),\quad
Ce^{-\mu_{k}/2}\leq|y-e^{i\beta_{l}}|<\epsilon.$

Moreover

(3.24) $v_{k}-V_{k}=O(\mu_{k}^{2}e^{-\mu_{k}})\quad\mbox{on}\quad\partial
B(Q_{l}^{k},\epsilon),\quad l=0,...,N.$

Also we observe from the expression of $V_{k}$ that

(3.25)
$V_{k}(x)=-\mu_{k}-2\log\frac{\mathfrak{h}_{k}(\delta_{k}e_{1})}{8}+4\log(N+1)-4(1+N)\log(\tau\delta_{k}^{-1})+E,$

for $x\in\partial\Omega_{k}$, thus

$v_{k}-V_{k}=O(\mu_{k}^{2}e^{-\mu_{k}})\quad\mbox{on}\quad\partial\Omega_{k}.$

Then the closeness of $v_{k}$ and $V_{k}$ on
$\Omega_{k}\setminus(\cup_{l}B(Q_{l}^{k},\epsilon))$ can be obtained by a
standard maximum principle argument: If we use $w_{k}$ to denote
$v_{k}-V_{k}$:

$w_{k}(z)=(v_{k}-V_{k})(z),$

then it is easy to see $w_{k}$ satisfies

$|\Delta w_{k}(z)|\leq Ce^{-\mu_{k}}|z|^{-4-2N},\quad\Omega_{k}\setminus
B_{2},\quad|w_{k}|\leq C\mu_{k}^{2}e^{-\mu_{k}}\mbox{ on }\partial
B_{2}\cup\partial\Omega_{k},$

Then $|w_{k}|$ can be majorized by
$Q(\mu_{k}^{2}e^{-\mu_{k}}-e^{-\mu_{k}}r^{-1-2N})$ for a large $Q>1$, which
yields the smallness of $v_{k}-V_{k}$ on $\Omega_{k}\setminus B_{2}$ as a
consequence. Proposition 3.1 is established. $\Box$

##  4\. First crucial bound for $v_{k}-V_{k}$

In this section we establish the first major estimate of $v_{k}-V_{k}$. The
main result in this section is

######  Proposition 4.1.

Let $w_{k}=v_{k}-V_{k}$, then

$|w_{k}(y)|\leq C\delta_{k},\quad y\in\Omega_{k}:=B(0,\tau\delta_{k}^{-1}).$

Proof of Proposition 4.1:

First we recall the equation for $v_{k}$ is (3.2), $v_{k}=$ constant on
$\partial B(0,\tau\delta_{k}^{-1})$. Moreover $v_{k}(e_{1})=\mu_{k}$. Recall
that $V_{k}$ defined in (3.7) satisfies

$\Delta
V_{k}+\mathfrak{h}_{k}(\delta_{k}e_{1})|y|^{2N}e^{V_{k}}=0,\quad\mbox{in}\quad\mathbb{R}^{2},\quad\int_{\mathbb{R}^{2}}|y|^{2N}e^{V_{k}}<\infty,$

$V_{k}$ has its local maximums at $e^{i\beta_{l}}$ for $l=0,...,N$ and
$V_{k}(e_{1})=\mu_{k}$. For $|y|\sim\delta_{k}^{-1}$,

$V_{k}(y)=-\mu_{k}-4(N+1)\log\delta_{k}^{-1}+C+O(\delta_{k}^{N+1})+O(e^{-\mu_{k}}).$

Let $\Omega_{k}=B(0,\tau\delta_{k}^{-1})$, we shall derive a precise, point-
wise estimate of $w_{k}$ in $B_{3}\setminus\cup_{l=1}^{N}B(Q_{l}^{k},\lambda)$
where $\lambda>0$ is a small number independent of $k$. Here we note that
among $N+1$ local maximum points, we already have $e_{1}$ as a common local
maximum point for both $v_{k}$ and $V_{k}$ and we shall prove that $w_{k}$ is
very small in $B_{3}$ if we exclude all bubbling disks except the one around
$e_{1}$. Before we carry out more specific computation we emphasize the
importance of

(4.1) $w_{k}(e_{1})=|\nabla w_{k}(e_{1})|=0.$

Now we write the equation of $w_{k}$ as

(4.2) $\Delta
w_{k}+\mathfrak{h}_{k}(\delta_{k}y)|y|^{2N}e^{\xi_{k}}w_{k}=(\mathfrak{h}_{k}(\delta_{k}e_{1})-\mathfrak{h}_{k}(\delta_{k}y))|y|^{2N}e^{V_{k}}$

in $\Omega_{k}$, where $\xi_{k}$ is obtained from the mean value theorem:

$e^{\xi_{k}(x)}=\left\\{\begin{array}[]{ll}\frac{e^{v_{k}(x)}-e^{V_{k}(x)}}{v_{k}(x)-V_{k}(x)},\quad\mbox{if}\quad
v_{k}(x)\neq V_{k}(x),\\\ \\\ e^{V_{k}(x)},\quad\mbox{if}\quad
v_{k}(x)=V_{k}(x).\end{array}\right.$

An equivalent form is

(4.3)
$e^{\xi_{k}(x)}=\int_{0}^{1}\frac{d}{dt}e^{tv_{k}(x)+(1-t)V_{k}(x)}dt=e^{V_{k}(x)}\big{(}1+\frac{1}{2}w_{k}(x)+O(w_{k}(x)^{2})\big{)}.$

For convenience we write the equation for $w_{k}$ as

(4.4) $\Delta
w_{k}+\mathfrak{h}_{k}(\delta_{k}y)|y|^{2N}e^{\xi_{k}}w_{k}=\delta_{k}\nabla\mathfrak{h}_{k}(\delta_{k}e_{1})\cdot(e_{1}-y)|y|^{2N}e^{V_{k}}+E_{1}$

where

$E_{1}=O(\delta_{k}^{2})|y-e_{1}|^{2}|y|^{2N}e^{V_{k}},\quad y\in\Omega_{k}.$

Let $M_{k}=\max_{x\in\bar{\Omega}_{k}}|w_{k}(x)|$. We shall get a
contradiction by assuming $M_{k}/\delta_{k}\to\infty$ at this moment. Set

$\tilde{w}_{k}(y)=w_{k}(y)/M_{k},\quad x\in\Omega_{k}.$

Clearly $\max_{x\in\Omega_{k}}|\tilde{w}_{k}(x)|=1$. The equation for
$\tilde{w}_{k}$ is

(4.5)
$\Delta\tilde{w}_{k}(y)+|y|^{2N}\mathfrak{h}_{k}(\delta_{k}e_{1})e^{\xi_{k}}\tilde{w}_{k}(y)=\frac{\delta_{k}}{M_{k}}\nabla\mathfrak{h}_{k}(\delta_{k}e_{1})\cdot(e_{1}-y)|y|^{2N}e^{V_{k}}+\tilde{E}_{1},$

in $\Omega_{k}$, where

(4.6) $\tilde{E}_{1}=o(\delta_{k})|y-e_{1}|^{2}|y|^{2N}e^{V_{k}},\quad
y\in\Omega_{k}.$

Now we give a more precise estimate of $e^{\xi_{k}}$. By Proposition 3.1

(4.7)
$\xi_{k}(y)=V_{k}(y)+\left\\{\begin{array}[]{ll}O(\mu_{k}e^{-\mu_{k}/2}),\quad|y-e_{1}|\leq
e^{-\mu_{k}/2},\\\ \\\
O(\mu_{k}^{2}e^{-\mu_{k}})|y-e_{1}|^{-1},\,\,e^{-\mu_{k}/2}\leq|y-e_{1}|\leq\delta_{0}.\end{array}\right.$

Since $V_{k}$ is not exactly symmetric around $e_{1}$, we shall replace the
re-scaled version of $V_{k}$ around $e_{1}$ by a radial function. Let $U_{k}$
be solutions of

(4.8) $\Delta
U_{k}+\mathfrak{h}_{k}(\delta_{k}e_{1})e^{U_{k}}=0,\quad\mbox{in}\quad\mathbb{R}^{2},\quad
U_{k}(0)=\max_{\mathbb{R}^{2}}U_{k}=0.$

Then we have

$U_{k}(z)=\log\frac{1}{(1+\frac{\mathfrak{h}_{k}(\delta_{k}e_{1})}{8}|z|^{2})^{2}}$

and

(4.9)
$V_{k}(e_{1}+\epsilon_{k}z)+2\log\epsilon_{k}=U_{k}(z)+O(\epsilon_{k})|z|+O(\mu_{k}^{2}\epsilon_{k}^{2}).$

Also we observe that

(4.10) $\log|e_{1}+\epsilon_{k}y|=O(\epsilon_{k})|y|.$

Thus, the combination of (4.7), (4.9) and (4.10) gives

(4.11) $\displaystyle
2N\log|e_{1}+\epsilon_{k}z|+\xi_{k}(e_{1}+\epsilon_{k}z)+2\log\epsilon_{k}-U_{k}(z)$
$\displaystyle=$ $\displaystyle O(\mu_{k}^{2}\epsilon_{k})(1+|z|)\quad
0\leq|z|<\delta_{0}\epsilon_{k}^{-1}.$

Since we shall use the re-scaled version, based on (4.11) we have

(4.12)
$\epsilon_{k}^{2}|e_{1}+\epsilon_{k}z|^{2N}e^{\xi_{k}(e_{1}+\epsilon_{k}z)}=e^{U_{k}(z)}+O(\mu_{k}^{2}\epsilon_{k})(1+|z|)^{-3}$

Here we note that the estimate in (4.11) is not optimal. In the following we
shall put the proof of Proposition 4.1 into a few estimates. In the first
estimate we prove

######  Lemma 4.1.

(4.13) $\tilde{w}_{k}(y)=o(1),\quad\nabla\tilde{w}_{k}=o(1)\quad\mbox{in}\quad
B(e_{1},\delta)\setminus B(e_{1},\delta/8)$

where $B(e_{1},3\delta)$ does not include other blowup points.

Proof of Lemma 4.1:

If (4.13) is not true, we have, without loss of generality that
$\tilde{w}_{k}\to c>0$. Note that $\tilde{w}_{k}$ tends to a global harmonic
function with removable singularity. So $\tilde{w}_{k}$ tends to constant.
Here we assume $c>0$ but the argument for $c<0$ is the same. Let

(4.14)
$W_{k}(z)=\tilde{w}_{k}(e_{1}+\epsilon_{k}z),\quad\epsilon_{k}=e^{-\frac{1}{2}\mu_{k}},$

then if we use $W$ to denote the limit of $W_{k}$, we have

$\Delta W+e^{U}W=0,\quad\mathbb{R}^{2},\quad|W|\leq 1,$

and $U$ is a solution of $\Delta U+e^{U}=0$ in $\mathbb{R}^{2}$ with
$\int_{\mathbb{R}^{2}}e^{U}<\infty$. Since $0$ is the local maximum of $U$,

$U(z)=\log\frac{1}{(1+\frac{1}{8}|z|^{2})^{2}}.$

Here we further claim that $W\equiv 0$ in $\mathbb{R}^{2}$ because
$W(0)=|\nabla W(0)|=0$, a fact well known based on the classification of the
kernel of the linearized operator. Going back to $W_{k}$, we have

$W_{k}(z)=o(1),\quad|z|\leq R_{k}\mbox{ for some }\quad R_{k}\to\infty.$

Based on the expression of $\tilde{w}_{k}$, (4.9) and (4.12) we write the
equation of $W_{k}$ as

(4.15) $\Delta
W_{k}(z)+\mathfrak{h}_{k}(\delta_{k}e_{1})e^{U_{k}(z)}W_{k}(z)=-\frac{\delta_{k}}{M_{k}}\nabla\mathfrak{h}_{k}(\delta_{k}e_{1})\cdot
z\epsilon_{k}e^{U_{k}(z)}+E_{2}^{k},$

for $|z|<\delta_{0}\epsilon_{k}^{-1}$ where

$E_{2}^{k}(z)=o(1)\mu_{k}^{2}\epsilon_{k}(1+|z|)^{-3}.$

Let

(4.16) $g_{0}^{k}(r)=\frac{1}{2\pi}\int_{0}^{2\pi}W_{k}(r,\theta)d\theta.$

Then clearly $g_{0}^{k}(r)\to c>0$ for $r\sim\epsilon_{k}^{-1}$. The equation
for $g_{0}^{k}$ is

$\displaystyle\frac{d^{2}}{dr^{2}}g_{0}^{k}(r)+\frac{1}{r}\frac{d}{dr}g_{0}^{k}(r)+\mathfrak{h}_{k}(\delta_{k}e_{1})e^{U_{k}(r)}g_{0}^{k}(r)=\tilde{E}_{0}^{k}(r)$
$\displaystyle g_{0}^{k}(0)=\frac{d}{dr}g_{0}^{k}(0)=0.$

where $\tilde{E}_{0}^{k}(r)$ has the same upper bound as that of
$E_{2}^{k}(r)$:

$|\tilde{E}_{0}^{k}(r)|\leq C\mu_{k}^{2}\epsilon_{k}(1+r)^{-3}.$

For the homogeneous equation, the two fundamental solutions are known:
$g_{01}$, $g_{02}$, where

$g_{01}=\frac{1-c_{1}r^{2}}{1+c_{1}r^{2}},\quad
c_{1}=\frac{\mathfrak{h}_{k}(\delta_{k}e_{1})}{8}.$

By the standard reduction of order process, $g_{02}(r)=O(\log r)$ for $r>1$.
Then it is easy to obtain, assuming $|W_{k}(z)|\leq 1$, that

$\displaystyle|g_{0}(r)|\leq
C|g_{01}(r)|\int_{0}^{r}s|\tilde{E}_{0}^{k}(s)g_{02}(s)|ds+C|g_{02}(r)|\int_{0}^{r}s|g_{01}(s)\tilde{E}_{0}^{k}(s)|ds$
$\displaystyle\leq C\mu_{k}^{2}\epsilon_{k}\log(2+r).\quad
0<r<\delta_{0}\epsilon_{k}^{-1}.$

Clearly this is a contradiction to (4.16). We have proved $c=0$, which means
$\tilde{w}_{k}=o(1)$ in $B(e_{1},\delta_{0})\setminus B(e_{1},\delta_{0}/8)$.
Then it is easy to use the equation for $\tilde{w}_{k}$ and standard Harnack
inequality to prove $\nabla\tilde{w}_{k}=o(1)$ in the same region. Lemma 4.1
is established. $\Box$

The second estimate is a more precise description of $\tilde{w}_{k}$ around
$e_{1}$:

######  Lemma 4.2.

For any given $\sigma\in(0,1)$ there exists $C>0$ such that

(4.17) $|\tilde{w}_{k}(e_{1}+\epsilon_{k}z)|\leq
C\epsilon_{k}^{\sigma}(1+|z|)^{\sigma},\quad 0<|z|<\tau\epsilon_{k}^{-1}.$

for some $\tau>0$.

Proof of Lemma 4.2: Let $W_{k}$ be defined as in (4.14) we recall the equation
for $W_{k}$ as (4.15). Here we let $\Omega_{wk}=B(0,\tau\epsilon_{k}^{-1})$
where $\tau\in(0,1)$ is chosen such that $W_{k}=o(1)$ on
$\partial\Omega_{wk}$. and

$\Lambda_{k}=\max_{z\in\Omega_{wk}}\frac{|W_{k}(z)|}{\epsilon_{k}^{\sigma}(1+|z|)^{\sigma}}.$

If (4.17) does not hold, $\Lambda_{k}\to\infty$ and we use $z_{k}$ to denote
where $\Lambda_{k}$ is attained. Note that because of the smallness of $W_{k}$
on $\partial\Omega_{wk}$, $z_{k}$ is an interior point. Let

$g_{k}(z)=\frac{W_{k}(z)}{\Lambda_{k}(1+|z_{k}|)^{\sigma}\epsilon_{k}^{\sigma}},\quad
z\in\Omega_{k,1},$

we see immediately that

(4.18)
$|g_{k}(z)|=\frac{|W_{k}(z)|}{\epsilon_{k}^{\sigma}\Lambda_{k}(1+|z|)^{\sigma}}\cdot\frac{(1+|z|)^{\sigma}}{(1+|z_{k}|)^{\sigma}}\leq\frac{(1+|z|)^{\sigma}}{(1+|z_{k}|)^{\sigma}}.$

Note that $\sigma$ can be as close to $1$ as needed. The equation of $g_{k}$
is

$\Delta
g_{k}(z)+\mathfrak{h}_{k}(\delta_{k}e_{1})e^{\xi_{k}}g_{k}=o(\epsilon_{k}^{1-\sigma})\frac{(1+|z|)^{-3}}{(1+|z_{k}|)^{\sigma}},\quad\mbox{in}\quad\Omega_{k,1}.$

Then we can obtain a contradiction to $|g_{k}(z_{k})|=1$ as follows: If
$\lim_{k\to\infty}z_{k}=P\in\mathbb{R}^{2}$, this is not possible because that
fact that $g_{k}(0)=|\nabla g_{k}(0)|=0$ and the sub-linear growth of $g_{k}$
in (4.18) implies that $g_{k}\to 0$ over any compact subset of
$\mathbb{R}^{2}$. So we have $|z_{k}|\to\infty$. But this would lead to a
contradiction again by using the Green’s representation of $g_{k}$:

(4.19) $\displaystyle\pm 1=g_{k}(z_{k})=g_{k}(z_{k})-g_{k}(0)$
$\displaystyle=\int_{\Omega_{k,1}}(G_{k}(z_{k},\eta)-G_{k}(0,\eta))(\mathfrak{h}_{k}(\delta_{k}e_{1})e^{\xi_{k}}g_{k}(\eta)+o(\epsilon_{k}^{1-\sigma})\frac{(1+|\eta|)^{-3}}{(1+|z_{k}|)^{\sigma}})d\eta+o(1).$

where $G_{k}(y,\eta)$ is the Green’s function on $\Omega_{k,1}$ and $o(1)$ in
the equation above comes from the smallness of $W_{k}$ on
$\partial\Omega_{k,1}$. Let $L_{k}=\tau\epsilon_{k}^{-1}$. From (3.13) we have

$G_{k}(z_{k},\eta)-G_{k}(0,\eta)=-\frac{1}{2\pi}\log|z_{k}-\eta|+\frac{1}{2\pi}\log|\frac{z_{k}}{|z_{k}|}-\frac{\eta
z_{k}}{L_{k}}|+\frac{1}{2\pi}\log|\eta|.$

Using this expression in (4.19) we obtain from elementary computation that the
right hand side of (4.19) is $o(1)$, a contradiction to $|g_{k}(z_{k})|=1$.
Lemma 4.2 is established. $\Box$

The smallness of $\tilde{w}_{k}$ around $e_{1}$ can be used to obtain the
following third key estimate:

######  Lemma 4.3.

(4.20) $\tilde{w}_{k}=o(1)\quad\mbox{in}\quad B(e^{i\beta_{l}},\tau)\quad
l=1,..,N.$

Proof of Lemma 4.3: We abuse the notation $W_{k}$ by defining it as

$W_{k}(z)=\tilde{w}_{k}(e^{i\beta_{l}}+\epsilon_{k}z),\quad
z\in\Omega_{k,l}:=B(0,\tau\epsilon_{k}^{-1}).$

Here we point out that the smallness of $\delta_{k}$ (see 2.8)implies that
$\epsilon_{k}^{-1}|Q_{l}^{k}-e^{i\beta_{l}}|\to 0$. So the scaling around
$e^{i\beta_{l}}$ or $Q_{l}^{k}$ does not affect the limit function.

$|e^{i\beta_{l}}+\epsilon_{k}z|^{2N}\mathfrak{h}_{k}(\delta_{k}e_{1})e^{\xi_{k}(e^{i\beta_{l}}+\epsilon_{k}z)}\to
e^{U(z)}$

where $U(z)$ is a solution of

$\Delta
U+e^{U}=0,\quad\mbox{in}\quad\mathbb{R}^{2},\quad\int_{\mathbb{R}^{2}}e^{U}<\infty.$

Here we recall that $\lim_{k\to\infty}\mathfrak{h}_{k}(\delta_{k}e_{1})=1$.
Since $W_{k}$ converges to a solution of the linearized equation:

$\Delta W+e^{U}W=0,\quad\mbox{in}\quad\mathbb{R}^{2}.$

$W$ can be written as a linear combination of three functions:

$W(x)=c_{0}\phi_{0}+c_{1}\phi_{1}+c_{2}\phi_{2},$

where

$\phi_{0}=\frac{1-\frac{1}{8}|x|^{2}}{1+\frac{1}{8}|x|^{2}}$
$\phi_{1}=\frac{x_{1}}{1+\frac{1}{8}|x|^{2}},\quad\phi_{2}=\frac{x_{2}}{1+\frac{1}{8}|x|^{2}}.$

The remaining part of the proof consisting of proving $c_{0}=0$ and
$c_{1}=c_{2}=0$. First we prove $c_{0}=0$.

Step one: $c_{0}=0$. First we write the equation for $W_{k}$ in a convenient
form. Since

$|e^{i\beta_{l}}+\epsilon_{k}z|^{2N}\mathfrak{h}_{k}(\delta_{k}e_{1})=\mathfrak{h}_{k}(\delta_{k}e_{1})+O(\epsilon_{k}z),$

and

$\epsilon_{k}^{2}e^{\xi_{k}(e^{i\beta_{l}}+\epsilon_{k}z)}=e^{U_{k}(z)}+O(\epsilon_{k})(1+|z|)^{-3}.$

Based on (4.5) we write the equation for $W_{k}$ as

(4.21) $\Delta
W_{k}(z)+\mathfrak{h}_{k}(\delta_{k}e_{1})e^{U_{k}}W_{k}=E_{l}^{k}(z)$

where

$E_{l}^{k}(z)=O(\epsilon_{k}^{1/2})(1+|z|)^{-3}+\sigma_{k}\mathfrak{h}_{k}(\delta_{k}e_{1})e^{U_{k}}\quad\mbox{in}\quad\Omega_{k,l}$

where $\sigma_{k}=\frac{\delta_{k}}{M_{k}}\to 0$. In order to prove $c_{0}=0$,
the key is to control the derivative of $W_{0}^{k}(r)$ where

$W_{0}^{k}(r)=\frac{1}{2\pi r}\int_{\partial B_{r}}W_{k}(re^{i\theta})dS,\quad
0<r<\tau\epsilon_{k}^{-1}.$

To obtain a control of $\frac{d}{dr}W_{0}^{k}(r)$ we use $\phi_{0}^{k}(r)$ as
the radial solution of

$\Delta\phi_{0}^{k}+\mathfrak{h}_{k}(\delta_{k}e_{1})e^{U_{k}}\phi_{0}^{k}=0,\quad\mbox{in
}\quad\mathbb{R}^{2}.$

When $k\to\infty$, $\phi_{0}^{k}\to c_{0}\phi_{0}$. Thus using the equation
for $\phi_{0}^{k}$ and $W_{k}$, we have

(4.22) $\int_{\partial
B_{r}}(\partial_{\nu}W_{k}\phi_{0}^{k}-\partial_{\nu}\phi_{0}^{k}W_{k})=o(\epsilon_{k}^{1/2})+\sigma_{k}\int_{B_{r}}\mathfrak{h}_{k}(\delta_{k}e_{1})e^{U_{k}}\phi_{k}.$

The last term above is

$\sigma_{k}\int_{B_{r}}\mathfrak{h}_{k}(\delta_{k}e_{1})e^{U_{k}}\phi_{k}=\sigma_{k}\int_{\partial
B_{r}}\partial_{\nu}\phi_{0}^{k}=o(1)r^{-3},\quad r>1.$

Thus from (4.22) we have

(4.23) $\frac{d}{dr}W_{0}^{k}(r)=o(\epsilon_{k}^{1/2})/r+O(1/r^{3}),\quad
1<r<\tau\epsilon_{k}^{-1}.$

Since we have known that

$W_{0}^{k}(\tau\epsilon_{k}^{-1})=o(1).$

By the fundamental theorem of calculus we have

$W_{0}^{k}(r)=W_{0}^{k}(\tau\epsilon_{k}^{-1})+\int_{\tau\epsilon_{k}^{-1}}^{r}(\frac{o(\epsilon_{k}^{\frac{1}{2}})}{s}+O(s^{-3}))ds=O(1/r^{2}),\quad
r\sim 1,\quad r\mbox{ is large}.$

$c_{0}=0$ is established because $W_{0}^{k}(r)\to c_{0}\phi_{0}$.

Step two $c_{1}=c_{2}=0$. We first observe that Lemma 4.3 follows from this.
Indeed, once we have proved $c_{1}=c_{2}=c_{0}=0$ around each
$e^{i\beta_{l}}$, it is easy to use maximum principle to prove
$\tilde{w}_{k}=o(1)$ in $B_{3}$ using $\tilde{w}_{k}=o(1)$ on $\partial B_{3}$
and the Green’s representation of $\tilde{w}_{k}$. The smallness of
$\tilde{w}_{k}$ immediately implies $\tilde{w}_{k}=o(1)$ in $B_{R}$ for any
fixed $R>>1$. Outside $B_{R}$, a crude estimate of $v_{k}$ is

$v_{k}(y)\leq-\mu_{k}-4(N+1)\log|y|+C,\quad 3<|y|<\tau\delta_{k}^{-1}.$

Using this and the Green’s representation of $w_{k}$ we can first observe that
the oscillation on each $\partial B_{r}$ is $o(1)$
($R<r<\tau\delta_{k}^{-1}/2$) and then by the Green’s representation of
$\tilde{w}_{k}$ and fast decay rate of $e^{V_{k}}$ we obtain
$\tilde{w}_{k}=o(1)$ in $B(0,\tau\delta_{k}^{-1})$. A contradiction to
$\max|\tilde{w}_{k}|=1$.

There are $N+1$ local maximums with one of them being $e_{1}$. Correspondingly
there are $N+1$ global solutions $V_{l,k}$ that approximate $v_{k}$ accurately
near $Q_{l}^{k}$ for $l=0,...,N$. Note that $Q_{0}^{k}=e_{1}$. For $V_{l,k}$
the expression is

$V_{l,k}=\log\frac{e^{\mu_{l}^{k}}}{(1+\frac{e^{\mu_{l}^{k}}}{D_{l}^{k}}|y^{N+1}-(e_{1}+p_{l}^{k})|^{2})^{2}},\quad
l=0,...,N,$

where $p_{l}^{k}=E$ and
$D_{l}^{k}=8(N+1)^{2}/\mathfrak{h}_{k}(\delta_{k}Q_{l}^{k})$. Note that a
simple computation gives

$Q_{l}^{k}=e^{i\beta_{l}}+\frac{p_{l}^{k}e^{i\beta_{l}}}{N+1}+O(|p_{l}^{k}|^{2}),\quad\beta_{l}=\frac{2l\pi}{N+2}.$

Let $M_{l,k}$ be the maximum of $|v_{k}-V_{l,k}|$ and we claim that all these
$M_{l,k}$ are comparable:

(4.24) $M_{l,k}\sim M_{s,k},\quad\forall s\neq l.$

The proof of (4.24) is as follows: We use $L_{s,l}$ to denote the limit of
$v_{k}-V_{l,k}$ around $Q_{s}^{k}$:

$\frac{(v_{k}-V_{l,k})(Q_{s}^{k}+\epsilon_{k}z)}{M_{l,k}}=L_{s,l}+o(1),\quad|z|\leq\tau\epsilon_{k}^{-1}$

where

$L_{s,l}=c_{1,s,l}\frac{z_{1}}{1+\frac{1}{8}|z|^{2}}+c_{2,s,l}\frac{z_{2}}{1+\frac{1}{8}|z|^{2}},\quad\mbox{and}\quad
L_{l,l}=0,\quad s=0,...,N.$

If all $c_{1,s,l}$ and $c_{2,s,l}$ are zero there is nothing to prove. So at
least one of them is not zero. For each $s\neq l$, by Lemma 4.2 we have

(4.25)
$v_{k}(Q_{s}^{k}+\epsilon_{k}z)-V_{s,k}(Q_{s}^{k}+\epsilon_{k}z)=O(\epsilon_{k}^{\sigma})(1+|z|)^{\sigma}M_{s,k},\quad|z|<\tau\epsilon_{k}^{-1}.$

Let $M_{k}=\max_{i}M_{i,k}$ ($i=0,...,N$) and we suppose $M_{k}=M_{l,k}$. Then
to determine $L_{s,l}$ we see that

$\frac{v_{k}(Q_{s}^{k}+\epsilon_{k}z)-V_{l,k}(Q_{s}^{k}+\epsilon_{k}z)}{M_{k}}=o(\epsilon_{k}^{\sigma})(1+|z|)^{\sigma}+\frac{V_{s,k}(Q_{s}^{k}+\epsilon_{k}z)-V_{l,k}(Q_{s}^{k}+\epsilon_{k}z)}{M_{k}}.$

We write $V_{s,k}(y)-V_{l,k}(y)$ as

$V_{s,k}(y)-V_{l,k}(y)=\mu_{s,k}-\mu_{l,k}+2A-A^{2}+O(|A|^{3})$

where

$A(y)=\frac{\frac{e^{\mu_{l}^{k}}}{D_{l}^{k}}|y^{N+1}-e_{1}-p_{l}^{k}|^{2}-\frac{e^{\mu_{s}^{k}}}{D_{s}^{k}}|y^{N+1}-e_{1}-p_{s}^{k}|^{2}}{1+\frac{e^{\mu_{s}^{k}}}{D_{s}^{k}}|y^{N+1}-e_{1}-p_{s}^{k}|^{2}}.$

From $A$ we have

(4.26) $\displaystyle
V_{s,k}(Q_{s}^{k}+\epsilon_{k}z)-V_{l,k}(Q_{s}^{k}+\epsilon_{k}z)$
$\displaystyle=$
$\displaystyle\phi_{1}+\phi_{2}+\phi_{3}+\phi_{4}+\mathfrak{R},$

where

$\displaystyle\phi_{1}=(\mu_{l,k}-\mu_{s,k})(1-\frac{(N+1)^{2}}{D_{s}^{k}}|z+\frac{N}{2}z^{2}e^{-i\beta_{s}}\epsilon_{k}|^{2})/B,$
$\displaystyle\phi_{2}=\frac{2(N+1)^{2}}{D_{s}^{k}}\delta_{k}\nabla\mathfrak{h}_{k}(\delta_{k}Q_{s}^{k})(Q_{l}^{k}-Q_{s}^{k})|z|^{2}/B$
$\displaystyle\phi_{3}=\frac{4e^{\mu_{l,k}-\mu_{s,k}}(N+1)^{2}}{D_{s}^{k}B}Re((z+\frac{N}{2}\epsilon_{k}e^{-i\beta_{s}}z^{2})(\frac{\bar{p}_{s}^{k}-\bar{p}_{l}^{k}}{\epsilon_{k}}e^{-i\beta_{s}}))$
$\displaystyle\phi_{4}=\frac{|p_{s}^{k}-p_{l}^{k}|^{2}}{\epsilon_{k}^{2}}\bigg{(}\frac{2}{D_{s}^{k}B}-\frac{2(N+1)^{2}|z|^{2}}{D_{s}^{2}B^{2}}-\frac{2(N+1)^{2}}{D_{s}^{2}B}|z|^{2}\cos(2\theta-2\theta_{st}-2\beta_{s})\bigg{)},$
$\displaystyle
B=1+\frac{(N+1)^{2}}{D_{s}^{k}}|z+\frac{N}{2}z^{2}e^{-i\beta_{s}}\epsilon_{k}|^{2},$

and $\mathfrak{R}_{k}$ is the collections of other insignificant terms. Here
we note that if we set $\epsilon_{l,k}=e^{-\mu_{l}^{k}/2}$, there is no
essential difference between $\epsilon_{l,k}$ and
$\epsilon_{k}=e^{-\frac{1}{2}\mu_{1,k}}$ because
$\epsilon_{l,k}=\epsilon_{k}+O(\epsilon_{k}E)$. If
$|\mu_{s,k}-\mu_{l,k}|/M_{k}\geq C$ there is no way to obtain a limit in the
form of $L_{s,l}$ mentioned before. Thus we must have
$|\mu_{s,k}-\mu_{l,k}|/M_{k}\to 0$. After simplification (see $\phi_{3}$ of
(4.26)) we have

(4.27)
$c_{1,l,s}=\lim_{k\to\infty}\frac{|p_{s}^{k}-p_{l}^{k}|}{2(N+1)M_{k}\epsilon_{k}}\cos(\beta_{s}+\theta_{sl}),\quad
c_{2,l,s}=\lim_{k\to\infty}\frac{|p_{s}^{k}-p_{l}^{k}|}{2(N+1)\epsilon_{k}M_{k}}sin(\beta_{s}+\theta_{sl})$

and $\theta_{sl}$ comes from
$p_{s}^{k}-p_{l}^{k}=|p_{s}^{k}-p_{l}^{k}|e^{i\theta_{sl}}$. We omit $k$ for
convenience. It is also important to observe that even if
$M_{k}=o(\epsilon_{k})$ we still have
$M_{k}\sim\max_{s}|p_{s}^{k}-p_{l}^{k}|/\epsilon_{k}$. Since each
$|p_{l}^{k}|=E$, an upper bound for $M_{k}$ is $M_{k}\leq
C\mu_{k}\epsilon_{k}$.

Equation (4.27) gives us the important observation:
$|c_{1,l,s}|+|c_{2,l,s}|\sim|p_{s}^{k}-p_{l}^{k}|/M_{k}$. So whenever
$|c_{1,l,s}|+|c_{2,l,s}|\neq 0$ we have
$|p_{s}^{k}-p_{l}^{k}|\sim\epsilon_{k}M_{k}$. In other words for each $l$,
$\epsilon_{k}M_{l,k}\sim\max_{t\neq l}|p_{t}^{k}-p_{l}^{k}|$. Hence for any
$t$, if $|p_{t}^{k}-p_{l}^{k}|\sim\epsilon_{k}M_{k}$, let $M_{t,k}$ be the
maximum of $|v_{k}-V_{t,k}|$, we have $M_{t,k}\sim M_{k}$. If all
$|p_{t}^{k}-p_{l}^{k}|\sim\epsilon_{k}M_{k}$ (4.24) is proved. So we assume
that at least one of them is too close to $p_{l}^{k}$, say
$|p_{s}^{k}-p_{l}^{k}|=o(1)\epsilon_{k}M_{k}$. We choose $t\neq l$ such that
$|p_{t}^{k}-p_{l}^{k}|\sim\epsilon_{k}M_{k}$. Then

$|p_{s}^{k}-p_{t}^{k}|\geq|p_{t}^{k}-p_{l}^{k}|-|p_{s}^{k}-p_{l}^{k}|\geq\frac{1}{2}|p_{t}^{k}-p_{l}^{k}|\sim\epsilon_{k}M_{k}.$

This means $M_{s,k}\sim M_{k}$. Thus (4.24) is established. From now on for
convenience we shall just use $M_{k}$.

Set $w_{l,k}=(v_{k}-V_{l,k})/M_{k}$, then we have $w_{l,k}(Q_{l}^{k})=|\nabla
w_{l,k}(Q_{l}^{k})|=0$,

(4.28)
$\lim_{k\to\infty}w_{l,k}(Q_{s}^{k}+\epsilon_{k}z)/M_{k}=\frac{c_{1,l,s}z_{1}+c_{2,l,s}z_{2}}{1+\frac{1}{8}|z|^{2}}$

and around each $Q_{s}^{k}$ (4.25) holds with $M_{s,k}$ replaced by $M_{k}$.

Now we use the estimates above to evaluate $\nabla w_{l,k}(Q_{l}^{k})=0$. By
elementary computation,

(4.29) $\displaystyle\nabla w_{l,k}(Q_{l}^{k})$ $\displaystyle=$
$\displaystyle-\frac{1}{2\pi}\int_{\Omega_{k}}\frac{Q_{l}^{k}-\eta}{|Q_{l}^{k}-\eta|^{2}}(w_{l,k}(\eta)\mathfrak{h}_{k}(\delta_{k}Q_{l}^{k})+\sigma_{k}\nabla\mathfrak{h}_{k}(\delta_{k}Q_{l}^{k})(\eta-
Q_{l}^{k}))|\eta|^{2N}e^{V_{k}}d\eta+o(\epsilon_{k}).$

Since the left hand side is $0$, the right hand side gives

(4.30) $\displaystyle\sum_{s\neq
l}\bigg{(}\int_{B(Q_{s}^{k},\tau)}\frac{Q_{l}^{k}-\eta}{|Q_{l}^{k}-\eta|^{2}}w_{l,k}(\eta)\mathfrak{h}_{k}(\delta_{k}Q_{l}^{k})|\eta|^{2N}e^{V_{k}}d\eta$
$\displaystyle+8\pi\sigma_{k}\nabla\log\mathfrak{h}_{k}(\delta_{k}Q_{s}^{k})(Q_{s}^{k}-Q_{l}^{k})\bigg{)}=o(\epsilon_{k}),$

where $\sigma_{k}=\delta_{k}/M_{k}\to 0$. To evaluate the first term in
(4.30), around each $Q_{s}^{k}$ we write

$\frac{Q_{l}^{k}-\eta}{|Q_{l}^{k}-\eta|^{2}}=\frac{Q_{l}^{k}-Q_{s}^{k}}{|Q_{l}^{k}-Q_{s}^{k}|^{2}}+F_{ls}(\eta)$

where

$F_{ls}(\eta)=(F_{ls}^{1}(\eta),F_{ls}^{2}(\eta))=\frac{Q_{l}^{k}-\eta}{|Q_{l}^{k}-\eta|^{2}}-\frac{Q_{l}^{k}-Q_{s}^{k}}{|Q_{l}^{k}-Q_{s}^{k}|^{2}}.$

We use (4.26) to obtain

(4.31)
$\displaystyle\int_{B(Q_{s}^{k},\tau)}(w_{l,k}(\eta)\mathfrak{h}_{k}(\delta_{k}Q_{l}^{k})|\eta|^{2N}e^{V_{k}}d\eta$
$\displaystyle+8\pi\sigma_{k}\nabla\log\mathfrak{h}_{k}(\delta_{k}Q_{s}^{k})(Q_{s}^{k}-Q_{l}^{k})=O(\epsilon_{k}^{\sigma})$

Note that in the computation above, the terms of $\phi_{1}$ and $\phi_{3}$
lead to $o(\epsilon_{k})$, the integration involving $\phi_{2}$ cancels with
the second term of (4.31). The integration involving $\phi_{4}$ provides the
leading term. More detailed information is the following: First for a global
solution

$V_{\mu,p}=\log\frac{e^{\mu}}{(1+\frac{e^{\mu}}{\lambda}|z^{N+1}-p|^{2})^{2}}$

of

$\Delta
V_{\mu,p}+\frac{8(N+1)^{2}}{\lambda}|z|^{2N}e^{V_{\mu,p}}=0,\quad\mbox{in
}\quad\mathbb{R}^{2},$

we have

(4.32)
$\int_{\mathbb{R}^{2}}\partial_{\mu}V_{\mu,p}|z|^{2N}e^{V_{\mu,p}}=\int_{\mathbb{R}^{2}}\frac{(1-\frac{e^{\mu}}{\lambda}|z^{N+1}-P|^{2})|z|^{2N}}{(1+\frac{e^{\mu}}{\lambda}|z^{N+1}-P|^{2})^{3}}dz=0.$

From this we use standard estimate that includes some cancellation to have

$\int_{B(0,\tau\epsilon_{k}^{-1})}\frac{\phi_{1}}{M_{k}}B^{-3}=o(\epsilon_{k}).$

From $V_{\mu,p}$ we also have

$\int_{\mathbb{R}^{2}}\partial_{P}V_{\mu,p}|y|^{2N}e^{V_{\mu,p}}=\int_{\mathbb{R}^{2}}\partial_{\bar{P}}V_{\mu,p}|y|^{2N}e^{V_{\mu,p}}=0,$

which gives

(4.33)
$\int_{\mathbb{R}^{2}}\frac{\frac{e^{\mu}}{\lambda}(\bar{z}^{N+1}-\bar{P})|z|^{2N}}{(1+\frac{e^{\mu}}{\lambda}|z^{N+1}-P|^{2})^{3}}=\int_{\mathbb{R}^{2}}\frac{\frac{e^{\mu}}{\lambda}(z^{N+1}-P)|z|^{2N}}{(1+\frac{e^{\mu}}{\lambda}|z^{N+1}-P|^{2})^{3}}=0.$

From here we use scaling and cancellation to have

$\int_{B(0,\tau\epsilon_{k}^{-1})}\frac{\phi_{3}}{M_{k}}B^{-3}=o(\epsilon_{k}).$

Thus (4.31) holds. Equation (4.31) also leads to a more accurate estimate of
$w_{l,k}$ in regions between bubbling disks. By the Green’s representation
formula of $w_{l,k}$ it is easy to have

$\displaystyle w_{l,k}(y)$ $\displaystyle=$
$\displaystyle-\frac{1}{2\pi}\int_{\Omega_{k}}\log\frac{|y-\eta|}{|Q_{l}^{k}-\eta|}(w_{l,k}(\eta)\mathfrak{h}_{k}(\delta_{k}Q_{l}^{k})+\nabla\mathfrak{h}_{k}(\delta_{k}Q_{l}^{k})(\eta-
Q_{l}^{k})|\eta|^{2N}e^{V_{k}}d\eta+o(\epsilon_{k})$

Writing

$\log\frac{|y-\eta|}{|Q_{l}^{k}-\eta|}=\log\frac{|y-Q_{s}^{k}|}{|Q_{l}^{k}-Q_{s}^{k}|}+(\log\frac{|y-\eta|}{|Q_{l}^{k}-\eta|}-\log\frac{|y-Q_{s}^{k}|}{|Q_{l}^{k}-Q_{s}^{k}|}),$

the integration related to the second term is $O(\epsilon_{k})$. The
integration involving the first term is $O(\epsilon_{k}^{\sigma})$ by (4.31).
Therefore

$|w_{l,k}(y)|=o(\epsilon_{k}^{\sigma}),\quad y\in
B_{3}\setminus\cup_{s=0}^{N}B(Q_{s}^{k},\tau)$

for $\sigma\in(0,1)$. Thus this extra control of $w_{l,k}$ away from bubbling
disks gives a better estimate than (4.25) around each $Q_{s}^{k}$:

(4.34) $|w_{s,k}(Q_{s}^{k}+\epsilon_{k}z)|/M_{k}\leq
o(\epsilon_{k})(1+|z|)^{\sigma},\quad|z|<\tau\epsilon_{k}^{-1}.$

From the decomposition in (4.26) we can now compute (4.29) in more detail:

(4.35)
$\displaystyle\int_{B(Q_{s}^{k},\tau)}w_{l,k}(\eta)\mathfrak{h}_{k}(\delta_{k}Q_{l}^{k})|\eta|^{2N}e^{V_{k}}d\eta+8\pi\sigma_{k}\nabla\log\mathfrak{h}_{k}(\delta_{k}Q_{s}^{k})(Q_{s}^{k}-Q_{l}^{k})$
$\displaystyle=\frac{\pi}{(N+1)^{2}}\frac{|p_{s}^{k}-p_{l}^{k}|^{2}}{\epsilon_{k}^{2}M_{k}^{2}}M_{k}+O(\epsilon_{k}^{\sigma}).$

(4.36)
$\displaystyle\int_{B(Q_{s}^{k},\tau)}F_{ls}^{1}(\eta)(w_{l,k}(\eta)\mathfrak{h}_{k}(\delta_{k}Q_{l}^{k})|\eta|^{2N}e^{V_{k}}d\eta$
$\displaystyle=\epsilon_{k}\int_{B(0,\tau\epsilon_{k}^{-1})}(\partial_{1}F_{ls}^{1}(Q_{s}^{k})z_{1}+\partial_{2}F^{1}_{ls}(Q_{s}^{k})z_{2})\frac{w_{l,k}(Q_{s}^{k}+\epsilon_{k}z)}{(1+\frac{1}{8}|z|^{2})^{2}}dz+o(\epsilon_{k}).$
$\displaystyle=4\pi\epsilon_{k}\frac{-c_{1,l,s}\cos(\beta_{s}+\beta_{l})-c_{2,l,s}\sin(\beta_{l}+\beta_{s})}{\sin^{2}(\frac{\beta_{l}-\beta_{s}}{2})}+o(\epsilon_{k}).$

(4.37)
$\displaystyle\int_{B(Q_{s}^{k},\tau)}F_{ls}^{2}(\eta)(w_{l,k}(\eta)\mathfrak{h}_{k}(\delta_{k}Q_{l}^{k})|\eta|^{2N}e^{V_{k}}d\eta$
$\displaystyle=\epsilon_{k}\int_{B(0,\tau\epsilon_{k}^{-1})}(\partial_{1}F_{ls}^{2}(Q_{s}^{k})z_{1}+\partial_{2}F^{2}_{ls}(Q_{s}^{k})z_{2})\frac{w_{l,k}(Q_{s}^{k}+\epsilon_{k}z)}{(1+\frac{1}{8}|z|^{2})^{2}}dz+o(\epsilon_{k}).$
$\displaystyle=4\pi\epsilon_{k}\frac{-c_{1,l,s}\sin(\beta_{s}+\beta_{l})+c_{2,l,s}\cos(\beta_{l}+\beta_{s})}{\sin^{2}(\frac{\beta_{l}-\beta_{s}}{2})}+o(\epsilon_{k}).$

Using (4.35),(4.36),(4.37) in (4.30), and by the expressions of
$c_{1,l,s},c_{2,l,s}$ in (4.27) we have

$\sum_{s\neq
l}\bigg{(}\frac{Q_{l}^{k}-Q_{s}^{k}}{|Q_{l}^{k}-Q_{s}^{k}|^{2}}\frac{\pi}{(N+1)^{2}}\frac{|p_{s}^{k}-p_{l}^{k}|^{2}}{\epsilon_{k}^{2}M_{k}^{2}}M_{k}-\frac{2\pi\epsilon_{k}}{\sin^{2}(\frac{\beta_{l}-\beta_{s}}{2})}e^{i(\beta_{l}-\theta_{sl})}\frac{|p_{s}^{k}-p_{l}^{k}|}{\epsilon_{k}M_{k}}\bigg{)}=o(\epsilon_{k}).$

Since $Q_{s}^{k}=e^{i\beta_{s}}+E$, the equation above can be simplified to

(4.38) $\displaystyle\sum_{s\neq
l}\bigg{(}(\frac{1}{2}-\frac{i\sin(\beta_{l}-\beta_{s})}{4\sin^{2}(\frac{\beta_{s}-\beta_{l}}{2})})\frac{\pi}{(N+1)^{2}}(\frac{|p_{l}^{k}-p_{s}^{k}|}{\epsilon_{k}M_{k}})^{2}M_{k}$
$\displaystyle-\frac{2\pi\epsilon_{k}e^{-i\theta_{sl}}}{\sin^{2}(\frac{p_{s}^{k}-p_{l}^{k}}{2})}\frac{|p_{s}^{k}-p_{l}^{k}|}{\epsilon_{k}M_{k}}\bigg{)}=o(\epsilon_{k}).$

Note that
$e^{-i\theta_{sl}}|p_{s}^{k}-p_{l}^{k}|=\bar{p}_{s}^{k}-\bar{p}_{l}^{k}$.
Taking the sum of all $l$ in (4.38) we have

$CM_{k}\sum_{s\neq
l}(|p_{s}^{k}-p_{l}^{k}|^{2}/(\epsilon_{k}M_{k})^{2})=o(\epsilon_{k}),\quad
C>0.$

If $M_{k}\geq C\epsilon_{k}$ for some $C>0$, we have
$|c_{1,s,l}|+|c_{2,s,l}|=0$ for all $s,l$. Finally if $M_{k}=o(\epsilon_{k})$,
the first term in (4.38) is ignored and (4.38) becomes

$\sum_{s\neq
l}\frac{\bar{p}_{s}^{k}-\bar{p}_{l}^{k}}{\sin^{2}(\frac{\beta_{s}-\beta_{l}}{2})\epsilon_{k}M_{k}}=o(1),\quad
l=0,1,...,N.$

From here it is also easy to obtain all
$|p_{l}^{k}|/(\epsilon_{k}M_{k})=o(1)$, which implies $c_{1,s,l}=c_{2,s,l}=0$
for all $s,l$. Thus Lemma 4.3 is established. $\Box$

Proposition 4.1 is an immediate consequence of Lemma 4.3. $\Box$.

##  5\. Proof of Theorem 1.1

Let $\hat{w}_{k}=w_{k}/\delta_{k}$. Then the equation for $\hat{w}_{k}$ is

(5.1)
$\Delta\hat{w}_{k}+|y|^{2N}e^{\xi_{k}}\hat{w}_{k}=\nabla\mathfrak{h}_{k}(0)\cdot(e_{1}-y)|y|^{2N}e^{V_{k}}+O(\delta_{k})e^{V_{k}}|y-e_{1}|^{2},$

in $\Omega_{k}$. By Proposition 4.1, $|\hat{w}_{k}(y)|\leq C$. Before we carry
out the remaining part of the proof we observe that $\hat{w}_{k}$ converges to
a harmonic function in $\mathbb{R}^{2}$ minus finite singular points. Since
$\hat{w}_{k}$ is bounded, all these singularities are removable. Thus
$\hat{w}_{k}$ converges to a constant. Based on the information around
$e_{1}$, we shall prove that this constant is $0$. However, looking at the
right hand side the equation,

$\nabla\mathfrak{h}_{k}(0)\cdot(e_{1}-y)|y|^{2N}e^{V_{k}}\rightharpoonup\sum_{l=1}^{N}8\pi\nabla\mathfrak{h}_{k}(0)\cdot(e_{1}-e^{i\beta_{l}})\delta_{e^{i\beta_{l}}}.$

If $\nabla\mathfrak{h}_{k}(0)\neq 0$ we would get a contradiction by comparing
the Pohozaev identities of $v_{k}$ and $V_{k}$.

Now we use the notation $W_{k}$ again and use Proposition 4.1 to rewrite the
equation for $W_{k}$. Let

$W_{k}(z)=\hat{w}_{k}(e_{1}+\epsilon_{k}z),\quad|z|<\delta_{0}\epsilon_{k}^{-1}$

for $\delta_{0}>0$ small. Then from Proposition 4.1 we have

(5.2)
$\mathfrak{h}_{k}(\delta_{k}y)=\mathfrak{h}_{k}(\delta_{k}e_{1})+\delta_{k}\nabla\mathfrak{h}_{k}(\delta_{k}e_{1})(y-e_{1})+O(\delta_{k}^{2})|y-e_{1}|^{2},$
(5.3) $|y|^{2N}=|e_{1}+\epsilon_{k}z|^{2N}=1+O(\epsilon_{k})|z|,$ (5.4)
$V_{k}(e_{1}+\epsilon_{k}z)+2\log\epsilon_{k}=U_{k}(z)+O(\epsilon_{k})|z|+O(\epsilon_{k}^{2})(\log(1+|z|))^{2}$

and

(5.5)
$\xi_{k}(e_{1}+\epsilon_{k}z)+2\log\epsilon_{k}=U_{k}(z)+O(\epsilon_{k})(1+|z|).$

Using (5.2),(5.3),(5.4) and (5.5) in (5.1) we write the equation of $W_{k}$ as

(5.6) $\Delta
W_{k}+\mathfrak{h}_{k}(\delta_{k}e_{1})e^{U_{k}(z)}W_{k}=-\epsilon_{k}\nabla\mathfrak{h}_{k}(0)\cdot
ze^{U_{k}(z)}+E_{w},\quad 0<|z|<\delta_{0}\epsilon_{k}^{-1}$

where

(5.7)
$E_{w}=O(\delta_{k})(1+|z|)^{-4}+O(\epsilon_{k})(1+|z|)^{-3}W_{k}(z),\quad|z|<\delta_{0}\epsilon_{k}^{-1}.$

At this moment we use $|W_{k}(z)|\leq C$ and a rough estimate of $E_{w}$ is

(5.8)
$E_{w}(z)=O(\epsilon_{k})(1+|z|)^{-3},\quad|z|<\delta_{0}\epsilon_{k}^{-1}.$

Since $\hat{w}_{k}$ obviously converges to a global harmonic function with
removable singularity, we have $\hat{w}_{k}\to\bar{c}$ for some
$\bar{c}\in\mathbb{R}$. Then we claim that

######  Lemma 5.1.

$\bar{c}=0$.

Proof of Lemma 5.1:

If $\bar{c}_{1}\neq 0$, we use $W_{k}(z)=\bar{c}+o(1)$ on
$B(0,\delta_{0}\epsilon_{k}^{-1})\setminus
B(0,\frac{1}{2}\delta_{0}\epsilon_{k}^{-1})$ and consider the projection of
$W_{k}$ on $1$:

$g_{0}(r)=\frac{1}{2\pi}\int_{0}^{2\pi}W_{k}(re^{i\theta})d\theta.$

If we use $F_{0}$ to denote the projection to $1$ of the right hand side we
have, using the rough estimate of $E_{w}$ in (5.8)

$g_{0}^{\prime\prime}(r)+\frac{1}{r}g_{0}^{\prime}(r)+\mathfrak{h}_{k}(\delta_{k}e_{1})e^{U_{k}(r)}g_{0}(r)=F_{0},\quad
0<r<\delta_{0}\epsilon_{k}^{-1}$

where

$F_{0}(r)=O(\epsilon_{k})(1+|z|)^{-3}.$

In addition we also have

$\lim_{k\to\infty}g_{0}(\delta_{0}\epsilon_{k}^{-1})=\bar{c}_{1}+o(1).$

Here the $O(\delta_{k})(1+|z|)^{-4}$ is absorbed because of the smallness of
$\delta$. For simplicity we omit $k$ in some notations. By the same argument
as in Lemma 4.1, we have

$g_{0}(r)=O(\epsilon_{k})(\log(2+r))^{2},\quad
0<r<\delta_{0}\epsilon_{k}^{-1}.$

Thus $\bar{c}_{1}=0$. Lemma 5.1 is established. $\Box$

Based on Lemma 5.1 and standard Harnack inequality for elliptic equations we
have

(5.9) $\tilde{w}_{k}(x)=o(1),\,\,\nabla\tilde{w}_{k}(x)=o(1),\,\,x\in
B_{3}\setminus(\cup_{l=1}^{N}(B(e^{i\beta_{l}},\delta_{0})\setminus
B(e^{i\beta_{l}},\delta_{0}/8))).$

Equation (5.9) is equivalent to $w_{k}=o(\delta_{k})$ and $\nabla
w_{k}=o(\delta_{k})$ in the same region.

Proof of Theorem 1.1:

For $s=1,...,N$ we consider the Pohozaev identity around $Q_{s}^{k}$. Let
$\Omega_{s,k}=B(Q_{s}^{k},r)$ for small $r>0$. For $v_{k}$ we have

(5.10)
$\displaystyle\int_{\Omega_{s,k}}\partial_{\xi}(|y|^{2N}\mathfrak{h}_{k}(\delta_{k}y))e^{v_{k}}-\int_{\partial\Omega_{s,k}}e^{v_{k}}|y|^{2N}\mathfrak{h}_{k}(\delta_{k}y)(\xi\cdot\nu)$
$\displaystyle=\int_{\partial\Omega_{s,k}}(\partial_{\nu}v_{k}\partial_{\xi}v_{k}-\frac{1}{2}|\nabla
v_{k}|^{2}(\xi\cdot\nu))dS.$

where $\xi$ is an arbitrary unit vector. Correspondingly the Pohozaev identity
for $V_{k}$ is

(5.11)
$\displaystyle\int_{\Omega_{s,k}}\partial_{\xi}(|y|^{2N}\mathfrak{h}_{k}(\delta_{k}e_{1}))e^{V_{k}}-\int_{\partial\Omega_{s,k}}e^{V_{k}}|y|^{2N}\mathfrak{h}_{k}(\delta_{k}e_{1})(\xi\cdot\nu)$
$\displaystyle=\int_{\partial\Omega_{s,k}}(\partial_{\nu}V_{k}\partial_{\xi}V_{k}-\frac{1}{2}|\nabla
V_{k}|^{2}(\xi\cdot\nu))dS.$

Using $w_{k}=v_{k}-V_{k}$ and $|w_{k}(y)|\leq C\delta_{k}$ we have

$\displaystyle\int_{\partial\Omega_{s,k}}(\partial_{\nu}v_{k}\partial_{\xi}v_{k}-\frac{1}{2}|\nabla
v_{k}|^{2}(\xi\cdot\nu))dS$ $\displaystyle=$
$\displaystyle\int_{\partial\Omega_{s,k}}(\partial_{\nu}V_{k}\partial_{\xi}V_{k}-\frac{1}{2}|\nabla
V_{k}|^{2}(\xi\cdot\nu))dS$
$\displaystyle+\int_{\partial\Omega_{s,k}}(\partial_{\nu}V_{k}\partial_{\xi}w_{k}+\partial_{\nu}w_{k}\partial_{\xi}V_{k}-(\nabla
V_{k}\cdot\nabla w_{k})(\xi\cdot\nu))dS+O(\delta_{k}^{2}).$

If we just use crude estimate: $\nabla w_{k}=o(\delta_{k})$, then

$\int_{\partial\Omega_{s,k}}(\partial_{\nu}v_{k}\partial_{\xi}v_{k}-\frac{1}{2}|\nabla
v_{k}|^{2}(\xi\cdot\nu))dS-\int_{\partial\Omega_{s,k}}(\partial_{\nu}V_{k}\partial_{\xi}V_{k}-\frac{1}{2}|\nabla
V_{k}|^{2}(\xi\cdot\nu))dS=o(\delta_{k}).$

The difference on the second terms is minor:

$\int_{\partial\Omega_{s,k}}e^{v_{k}}|y|^{2N}\mathfrak{h}_{k}(\delta_{k}y)(\xi\cdot\nu)-\int_{\partial\Omega_{s,k}}e^{V_{k}}|y|^{2N}\mathfrak{h}_{k}(\delta_{k}e_{1})(\xi\cdot\nu)=O(\delta_{k}\epsilon_{k}^{2}).$

To evaluate the first term, we use

(5.12)
$\displaystyle\partial_{\xi}(|y|^{2N}\mathfrak{h}_{k}(\delta_{k}y))e^{v_{k}}$
$\displaystyle=$
$\displaystyle\partial_{\xi}(|y|^{2N}\mathfrak{h}_{k}(\delta_{k}e_{1})+|y|^{2N}\delta_{k}\nabla\mathfrak{h}_{k}(\delta_{k}e_{1})(y-e_{1})+O(\delta_{k}^{2}))e^{V_{k}}(1+w_{k}+O(\delta_{k}^{2}))$
$\displaystyle=$
$\displaystyle\partial_{\xi}(|y|^{2N})\mathfrak{h}_{k}(\delta_{k}e_{1})e^{V_{k}}+\delta_{k}\partial_{\xi}(|y|^{2N}\nabla\mathfrak{h}_{k}(\delta_{k}e_{1})(y-e_{1}))e^{V_{k}}$
$\displaystyle+\partial_{\xi}(|y|^{2N}\mathfrak{h}_{k}(\delta_{k}e_{1}))e^{V_{k}}w_{k}+O(\delta_{k}^{2})e^{V_{k}}.$

For the third term on the right hand side of (5.12) we use the equation for
$w_{k}$:

$\Delta
w_{k}+\mathfrak{h}_{k}(\delta_{k}e_{1})e^{V_{k}}|y|^{2N}w_{k}=-\delta_{k}\nabla\mathfrak{h}_{k}(\delta_{k}e_{1})\cdot(y-e_{1})|y|^{2N}e^{V_{k}}+O(\delta_{k}^{2})e^{V_{k}}|y|^{2N}.$

From integration by parts we have

$\displaystyle\int_{\Omega_{s,k}}\partial_{\xi}(|y|^{2N})\mathfrak{h}_{k}(\delta_{k}e_{1})e^{V_{k}}w_{k}$
$\displaystyle=$ $\displaystyle
2N\int_{\Omega_{s,k}}|y|^{2N-2}y_{\xi}\mathfrak{h}_{k}(\delta_{k}e_{1})e^{V_{k}}w_{k}$
$\displaystyle=$ $\displaystyle
2N\int_{\Omega_{s,k}}\frac{y_{\xi}}{|y|^{2}}(-\Delta
w_{k}-\delta_{k}\nabla\mathfrak{h}_{k}(\delta_{k}e_{1})(y-e_{1})|y|^{2N}e^{V_{k}}+O(\delta_{k}^{2})e^{V_{k}})$
$\displaystyle=$
$\displaystyle-2N\delta_{k}\int_{\Omega_{s,k}}\frac{y_{\xi}}{|y|^{2}}\nabla\mathfrak{h}_{k}(\delta_{k}e_{1})(y-e_{1})|y|^{2N}e^{V_{k}}$
$\displaystyle+2N\int_{\partial\Omega_{s,k}}(\partial_{\nu}(\frac{y_{\xi}}{|y|^{2}})w_{k}-\partial_{\nu}w_{k}\frac{y_{\xi}}{|y|^{2}})+O(\delta_{k}^{2})$
(5.13) $\displaystyle=$
$\displaystyle-16N\delta_{k}\pi(e^{i\beta_{s}}\cdot\xi)\nabla\mathfrak{h}_{k}(\delta_{k}e_{1})(e^{i\beta_{s}}-e_{1})+o(\delta_{k}),$

where we have used $\nabla w_{k},w_{k}=o(\delta_{k})$ on
$\partial\Omega_{s,k}$. For the second term on the right hand side of (5.12),
we have

(5.14)
$\displaystyle\int_{\Omega_{s,k}}\delta_{k}\partial_{\xi}(|y|^{2N}\nabla\mathfrak{h}_{k}(\delta_{k}e_{1})(y-e_{1}))e^{V_{k}}$
$\displaystyle=$ $\displaystyle
2N\delta_{k}\int_{\Omega_{s,k}}y_{\xi}|y|^{2N-2}\nabla\mathfrak{h}_{k}(\delta_{k}e_{1})(y-e_{1})e^{V_{k}}+\delta_{k}\int_{\Omega_{s,k}}|y|^{2N}\partial_{\xi}\mathfrak{h}_{k}(\delta_{k}e_{1})e^{V_{k}}$
$\displaystyle=$ $\displaystyle
16N\pi\delta_{k}(e^{i\beta_{s}}\cdot\xi)\nabla\mathfrak{h}_{k}(\delta_{k}e_{1})(e^{i\beta_{s}}-e_{1})+8\pi\delta_{k}\partial_{\xi}\mathfrak{h}_{k}(\delta_{k}e_{1})+o(\delta_{k}).$

Using (5) and (5.14) in the difference between (5.10) and (5.11), we have

$\delta_{k}\partial_{\xi}\mathfrak{h}_{k}(\delta_{k}e_{1})=o(\delta_{k}).$

Thus $\nabla\mathfrak{h}_{k}(\delta_{k}e_{1})=o(1)$. Theorem 1.1 is
established. $\Box$

## References

  * [1] J. Ambjorn and P. Olesen, A condensate solution of the electroweak theory which interpolates between the broken and the symmetry phase, Nucl. Phys. B 330 (1990), 193-204.
  * [2] Bartolucci, Daniele; Jevnikar, Aleks; Lee, Youngae; Yang, Wen Uniqueness of bubbling solutions of mean field equations. J. Math. Pures Appl.  (9) 123 (2019), 78–126.
  * [3] D. Bartolucci, G. Tarantello, Asymptotic Blow-up Analysis for Singular Liouville type Equations with Applications, J. Diff. Eq.  262, n. 7 (2017), 3887-3931
  * [4] E. Caglioti, P.L. Lions, C. Marchioro , M. Pulvirenti, A special class of stationary flows for two-dimensional Euler equations: A statistical mechanics description, Comm. Math. Phys. 143 (1992), 501-525.
  * [5] E. Caglioti, P.L. Lions, C. Marchioro , M. Pulvirenti, A special class of stationary flows for two-dimensional Euler equations: A statistical mechanics description, part II, Comm. Math. Phys. 174 (1995), 229-260.
  * [6] C. C. Chen, C. S. Lin, Sharp estimates for solutions of multi-bubbles in compact Riemann surfaces. Comm. Pure Appl. Math.  55 (2002), no. 6, 728-771.
  * [7] Chen, Chiun-Chuan; Lin, Chang-Shou Mean field equation of Liouville type with singular data: topological degree. Comm. Pure Appl. Math. 68 (2015), no. 6, 887-947.
  * [8] M. Del Pino, P. Esposito , M. Musso, Nondegeneracy of entire solutions of a singular Liouville equation, Proc. of AMS  140 (2012), no. 2 581-588.
  * [9] M. del Pino, P. Esposito, M. Musso, Two-dimensional Euler flows with concentrated vorticities, Trans. Amer. Math. Soc. 362 (2010), no. 12, 6381-6395.
  * [10] A. Eremenko, A. Gabrielov, V. Tarasov, Metrics with conic singularities and spherical polygons, Illinois J. Math, vol 58, no3, 2014, 739-755.
  * [11] M. Gluck, Asymptotic behavior of blow up solutions to a class of prescribing Gauss curvature equations. Nonlinear Anal. 75 (2012), no. 15, 5787-5796.
  * [12] Y. Gu, L. Zhang, Degree counting theorems for singular Liouville systems, to appear on Annali della Scuola Normale Superiore di Pisa, Classe di Scienze (5) Vol. XXI (2020), 1103-1135.
  * [13] R. Jackiw , E.J. Weinberg, Self-Dual Chern Simons vortices, Phys. Rev. Lett. 64 (1990), 2234-2237.
  * [14] Kuo, Ting-Jung; Lin, Chang-Shou Estimates of the mean field equations with integer singular sources: non-simple blowup. J. Differential Geom. 103 (2016), no. 3, 377-424.
  * [15] C. S. Lin, S. Yan, On the mean field type bubbling solutions for Chern–Simons–Higgs equation, Adv. Math. 338 (2018) 1141–1188.
  * [16] Prajapat, J.; Tarantello, G. On a class of elliptic problems in R2: symmetry and uniqueness results. Proc. Roy. Soc. Edinburgh Sect. A 131 (2001), no. 4, 967-985.
  * [17] M. Troyanov, Prescribing curvature on compact surfaces with conical singularities, Trans. AMS, 324 (1991) 793-821.
  * [18] M. Troyanov, Metrics of constant curvature on a sphere with two conical singularities, Lect. Notes Math., 1410, Springer, NY, 1989, 296-308.
  * [19] Juncheng Wei, Lei Zhang, Estimates for Liouville equation with quantized singularities. Adv. Math.  380 (2021), 107606.
  * [20] Lina Wu, Lei Zhang, Uniqueness of bubbling solutions of mean field equations with non-quantized singularities, Communications in Contemporary Mathematics,(2021), no. 4, 2050038, 39 pp.
  * [21] L. Zhang, Blowup solutions of some nonlinear elliptic equations involving exponential nonlinearities. Comm. Math. Phys. 268 (2006), no. 1, 105–133.
  * [22] L. Zhang, Vanishing estimates for fully bubbling solutions of $SU(n+1)$ Toda Systems at a singular source, Int. Math. Res. Not.  2020, no. 18, 5774-–5795.