# Anisotropy fluctuation and correlation in central $\alpha$-clustered
12C+197Au collisions

L. Ma<EMAIL_ADDRESS>Key Laboratory of Nuclear Physics and Ion-beam
Application (MOE), Institute of Modern Physics, Fudan University, Shanghai
200433, China  Y. G. Ma<EMAIL_ADDRESS>Key Laboratory of Nuclear
Physics and Ion-beam Application (MOE), Institute of Modern Physics, Fudan
University, Shanghai 200433, China  Shanghai Institute of Applied Physics,
Chinese Academy of Sciences, Shanghai 201800, China  S. Zhang
<EMAIL_ADDRESS>Key Laboratory of Nuclear Physics and Ion-beam
Application (MOE), Institute of Modern Physics, Fudan University, Shanghai
200433, China

###### Abstract

In the framework of a multi-phase transport model, fluctuation and correlation
of the azimuthal anisotropies in 12C+197Au collisions at $\sqrt{S_{NN}}$ = 200
GeV are explored. Properties of the initial eccentricity and final harmonic
flow fluctuation resulted from 12C+197Au collisions with Woods-Saxon
configuration and triangular $\alpha$-clustering configurations of 12C are
investigated via scaled variance, skewness and kurtosis. Comparisons are made
between results from $\alpha$-clustered configurations and Woods-Saxon
configuration. The triangular flow fluctuation is found to have particular
sensitivity in distinguishing triangular $\alpha$-clustering structure of 12C.
Furthermore, correlations between initial eccentricities and final flow
harmonics are studied with strong multiplicity dependence observed for the
correlation functions.

###### pacs:

25.75.-q

##  I Introduction

Anisotropic flow as a typical collective behavior of particles emitted was
proved to be a good observable to study the equation-of-state of the hot and
dense matter created in the early stage of ultra-relativistic heavy-ion
collisions Ollitrault (1992); Voloshin et al. (2008a); Kolb et al. (1999);
Ackermann et al. (2001); Teaney et al. (2001); Romatschke and Romatschke
(2007); Luo and Xu (2017). The flow harmonics defined as $v_{n}$ ($n$ = 2,3..)
which are the Fourier expansion coefficients of the azimuthal anisotropy of
the final-state particles are suggested to be sensitive to not only the early
partonic dynamics but also the transportation properties of the source
evolution Adams et al. (2004); Adamczyk et al. (2013); Adare et al. (2011).
Theoretically, the hydrodynamic model and the multi-phase transport model have
been widely used to make predictions and was successful in giving comparable
descriptions on experimental measurements Heinz (2005); Gale et al. (2013);
Qiu and Heinz (2012); Song et al. (2011, 2017); Chen et al. (2004); Huang et
al. (2020).

The essential role of the initial collision geometry was realized when looking
into the flow harmonics of different collision systems scaled by initial
eccentricities Alver et al. (2007); Voloshin (2006). Important information
about the fluctuating initial conditions is believed to be transferred to the
final flow anisotropy during source expansion. In the past two decades,
significant efforts have been made on the studies of initial geometry
fluctuation effect on the final flow observables Voloshin et al. (2008b);
Lacey et al. (2011, 2014); Alver and Roland (2010); Derradi de Souza et al.
(2012); Ma and Wang (2011); Ma et al. (2014); Han et al. (2011); Wang et al.
(2013). In addition, experimental measurements of the event-by-event
anisotropic flow fluctuation suggest the close correlation between flow
fluctuation and the fluctuations of the initial geometry Alver et al. (2010);
Agakishiev et al. (2012); Sorensen (2007); Margutti et al. (2019). It was
realized that the flow fluctuation on the event-by-event basis provides good
access to the initial eccentricity fluctuation, elucidating both source
evolution and initial fluctuation properties Andrade et al. (2006); Petersen
et al. (2010); Ma et al. (2016).

The nuclear cluster is one of the essential features of a nuclear system.
Systematic studies on the $\alpha$ cluster in light nuclei have been performed
for more than 40 years Brink et al. (1970); von Oertzen et al. (2006); Freer
et al. (2018); He et al. (2014, 2016); Kanada-En’yo et al. (2015); Liu and Ye
(2018); Huang et al. (2017); Huang and Ma (2020). Recent studies suggest that
the intrinsic structure of light nuclei can be captured by the “snapshots”
made in relativistic nuclear collisions by colliding light nucleus on the
heavy nuclei target Broniowski and Ruiz Arriola (2014). As the final-state
anisotropy inherits geometric information from the initial state, it was
proposed that the nuclear clustering configurations can be explored by
harmonic flow in the final stage of the nuclear collisions. Extensive
theoretical studies have been performed on harmonic flow as a probe for the
$\alpha$-clustering structure of light nuclei in nuclear collisions involving
12C or 16O as the projectile and it is found that the intrinsic nuclear
structure, predetermined by the arrangement of the $\alpha$ clusters, leads to
quantitative difference of the flow measurements Zhang et al. (2017, 2018);
Guo et al. (2019); Xu et al. (2018); Guo et al. (2017). In addition, flow
fluctuation characterized by the ratio of the cumulant flow is also proposed
to be an approach for investigating clustering in light nuclei Bożek et al.
(2014); Rybczyński et al. (2018). Studying the effects of initial
$\alpha$-clustering on flow fluctuation offers valuable information about the
sensitivity of experimental measurement of flow fluctuation in probing the
signatures of $\alpha$-clustering configurations in light nuclei.

In this paper, we present systematic simulation study of the fluctuation and
correlation properties of initial eccentricity and final flow harmonics for
12C+197Au collision at 200 GeV with a multi-phase transport model. Influence
of $\alpha$-clustering on flow fluctuation and initial-final anisotropy
correlation is investigated in particular. This paper is organized as follows:
In Sec. II, model and research methods are described. In Sec. III, results and
discussion are presented. In Sec. IV, a brief summary is given.

##  II Model and method

A multi-phase transport model (AMPT) is employed for studying 12C+197Au
collisions. The model consists of four main components: the initial condition,
partonic interactions, hadronization, and hadronic interactions Zhang et al.
(2000); Lin et al. (2005). Starting from Monte Carlo Glauber initial
conditions, phase space distributions of minijet partons and soft string
excitations are generated from the Heavy Ion Jet Interaction Generator
(HIJING) model Wang and Gyulassy (1991). In the string melting scenario, both
excited strings and minijet partons are decomposed into partons followed by
elastic partonic scatterings. Scatterings among partons are treated according
to a parton cascade model - Zhang’s parton cascade (ZPC) model which includes
parton-parton elastic scattering with cross sections obtained from the
theoretical calculations Zhang (1998). After partons stop interacting with
each other, a simple quark coalescence model is used to combine partons into
hadrons. Partonic matter is then turned into hadronic matter and the
subsequential hadronic interactions are simulated using a relativistic
transport model (ART) including both elastic and inelastic scattering
descriptions for baryon-baryon, baryon-meson and meson-meson interactions Li
and Ko (1995). With properly choosing partonic scattering cross section, the
AMPT model was successful in describing many experimental observations in
heavy-ion collisions at RHIC and LHC energies Xu et al. (2016); Zhou et al.
(2016); Jin et al. (2018); Wang et al. (2019).

Figure 1: (Color online) Nucleon radial density distributions of 12C for the
Woods-Saxon (W-S) and Triangle $\alpha$-clustering configurations. Results are
normalized to the number of nucleons of 12C.

Simulation events are generated for 12C+197Au collisions at center-of-mass
energy of 200 GeV. By default, the initial nucleon distribution in 12C and
197Au are initialized event-by-event according to Woods-Saxon distribution
from the HIJING model. For studying the $\alpha$-clustering effect, 12C is
configured with three $\alpha$ clusters in triangular structure based on the
information given by an extended quantum molecular dynamics model (EQMD)
Maruyama et al. (1996). The distribution of the radial centers of the $\alpha$
clusters are assumed to be Gaussian and initialized according to
$e^{-0.5(\frac{r-r_{c}(\alpha)}{\sigma_{r_{c}}(\alpha)})^{2}}$, where
$r_{c}(\alpha)$ is the averaged radial center coordinate of an $\alpha$
cluster and $\sigma_{r_{c}(\alpha)}$ is the width of the distribution.
Nucleons inside each $\alpha$ cluster are initialized according to Woods-Saxon
distribution. In our simulation study, the parameters for configuring the 12C
are taken from the EQMD calculation with $r_{c}(\alpha)$ = 1.8 fm and
$\sigma_{r_{c}}(\alpha)$ = 0.1 fm for the triangular $\alpha$-clustering
configuration. Details about the methodology applied for the initialization of
the collision system can be found in Refs. He et al. (2014, 2016). Fig. 1
shows the comparison of the normalized nucleon radial density distributions
initialized with Woods-Saxon and triangular $\alpha$-clustering configurations
of 12C in this study. One can find that though the root mean square radius of
the two configurations are comparable, the difference in the shape of the
density distribution can be clearly seen.

##  III Results and Discussion

###  III.1 Eccentricity fluctuation in 12C+197Au collisions

For nuclear-nuclear collision, the initial geometric anisotropy of the
collision zone in the transverse plane (perpendicular to the beam direction)
can be characterized with eccentricity. It was argued that the magnitude and
trend of the eccentricity and its fluctuation imply specifically testable
predictions for final stage harmonic flow Lacey et al. (2011); Drescher and
Nara (2007); Broniowski et al. (2007). Definition of the eccentricity for the
$n$th-order harmonic in the coordinate space of the participant nucleons or
partons for a single collision event is in the form:

$\varepsilon_{n}\left\\{P\right\\}=\frac{{\sqrt{\left\langle{r^{n}\cos(n\varphi)}\right\rangle^{2}+\left\langle{r^{n}\sin(n\varphi)}\right\rangle^{2}}}}{{\left\langle{r^{n}}\right\rangle}},$
(1)

where $r$ and $\varphi$ are position and azimuthal angles of each nucleon or
parton in the transverse plane. $\varepsilon_{n}\left\\{P\right\\}$
characterizes the eccentricity through the distribution of participant
nucleons or partons which naturally contains event-by-event fluctuation.
$\varepsilon_{n}\left\\{P\right\\}$ defined in this way is usually named as a
“participant eccentricity”.

Initial eccentricity can also be quantified by cumulants of
$\varepsilon_{n}\left\\{P\right\\}$ following the same way as in Ref. Miller
and Snellings (2003). The definitions in terms of two-particle and four-
particle cumulant moments of the eccentricities are in the form:

$\begin{split}&\varepsilon_{n}\left\\{2\right\\}=\sqrt{\left\langle{\varepsilon^{2}_{n}\left\\{P\right\\}}\right\rangle},\\\
&\varepsilon_{n}\left\\{4\right\\}=(2\left\langle{\varepsilon^{2}_{n}\left\\{P\right\\}}\right\rangle^{2}-\left\langle{\varepsilon^{4}_{n}\left\\{P\right\\}}\right\rangle)^{1/4}.\end{split}$
(2)

Figure 2: (Color online) Initial eccentricity $\varepsilon_{n}$ (n=2,3,4) as a
function of $N_{track}$ for 12C+197Au collisions at $\sqrt{S_{NN}}$ = 200 GeV
in AMPT model. Results from the Triangle $\alpha$-clustering configuration are
compared with Woods-Saxon (W-S) configuration without $\alpha$ clusterization.

Figure 3: (Color online) Eccentricity fluctuation $R_{\varepsilon_{n}}$ as a
function of $N_{track}$ for 12C+197Au collisions at $\sqrt{S_{NN}}$ = 200 GeV
in AMPT model. Results from the Triangle clustering configuration are compared
with Woods-Saxon (W-S) configuration without $\alpha$ clusterization. Panels
(a-c) show scaled standard deviations $R_{\varepsilon_{n}}$ for different
order of harmonics.

In central nuclear collision involving clustered light nuclei colliding on a
heavy-ion target, the profile of the initial collision zone inherits intrinsic
geometry of the nucleon distribution of the light nuclei. The clustering
configuration of nuclei is expected to be reflected by the eccentricity
coefficients. In Fig. 2, simulation results of $\varepsilon_{n}$ are shown up
to fourth order of the harmonic as a function of $N_{track}$ for central
12C+197Au collisions with different $\alpha$-clustering configurations of 12C.
$N_{track}$ is defined to be the number of charged particles within rapidity
window, 1.0 $<$ y $<$ 1.0 and transverse momentum window 0.2 $<$ $p_{T}$ $<$
6.0 GeV/c in the collision event. $\varepsilon_{n}$ results from two- and
four-particle cumulant definitions by Eq. 2 are compared with participant
definition by Eq. 1. It is found $\varepsilon_{n}\left\\{P\right\\}$ (n=2,3,4)
for all the configurations are quantitatively between
$\varepsilon_{n}\left\\{2\right\\}$ and $\varepsilon_{n}\left\\{4\right\\}$
over the whole $N_{track}$ range. Non-monotonous dependence on $N_{track}$ can
be seen for $\varepsilon_{n}$ resulted from the Triangle and unclustered
Woods-Saxon (W-S) configurations. In comparison, Triangle and W-S
configurations exhibit clear decreasing trend for $N_{track}$ $<$ 100\.
$\varepsilon_{3}$ and $\varepsilon_{4}$ are similar in trend as a function of
$N_{track}$ but show different magnitude orderings for different
$\alpha$-clustered cases.

Event-by-event fluctuation of the eccentricity is studied by first looking
into the scaled standard deviations (scaled variance) of eccentricity
distribution defined as:

$R_{\varepsilon_{n}}=\sigma_{\varepsilon_{n}}/\langle\varepsilon_{n}\rangle=\sqrt{\frac{\langle\varepsilon^{2}_{n}\rangle-\langle\varepsilon_{n}\rangle^{2}}{\langle\varepsilon_{n}\rangle^{2}}},$
(3)

where $\sigma_{\varepsilon_{n}}$ is the variance which quantifies the absolute
fluctuation. The brackets denote event averaging. Fig. 3 comparatively shows
the scaled standard deviation of eccentricity fluctuation
$R_{\varepsilon_{n}}$ (n=2,3,4) for W-S and $\alpha$-clustering configurations
of 12C. In comparison, $R_{\varepsilon_{2}}$ and $R_{\varepsilon_{4}}$ show a
similar decreasing trend as a function of $N_{track}$. This is consistent with
the expectation that quadrangularity $\varepsilon_{4}$ is strongly correlated
with ellipticity $\varepsilon_{2}$ and their fluctuations should behave
similarly Qiu and Heinz (2011). It is interesting to see the difference
between W-S and Triangle configurations that $R_{\varepsilon_{3}}$ shows
different ordering compared with $R_{\varepsilon_{2}}$ and
$R_{\varepsilon_{4}}$. It could be understood in view of the geometrical
origin that triangular structure from the intrinsic geometric bias contributes
more to the ellipticity or quadrangularity fluctuation but less triangularity
fluctuation. It is also found that for both $\varepsilon_{n}$ and
$R_{\varepsilon_{n}}$ the discrepancy between different configurations tend to
converge when approaching low $N_{track}$ and become diverge at higher
$N_{track}$.

In addition to the standard deviation, skewness and kurtosis which
characterize the non-Gaussian fluctuation properties have been used to study
flow fluctuation and explain possible contributions to the splitting of
higher-order flow cumulants Giacalone et al. (2017); Bhalerao et al. (2019).
The standardized skewness and kurtosis of $\varepsilon$ fluctuations are
defined as:

$\begin{split}&S_{\varepsilon_{n}}=\frac{\langle(\varepsilon_{n}-\langle\varepsilon_{n}\rangle)^{3}\rangle}{\langle(\varepsilon_{n}-\langle\varepsilon_{n}\rangle)^{2}\rangle^{3/2}},\\\
&K_{\varepsilon_{n}}=\frac{\langle(\varepsilon_{n}-\langle\varepsilon_{n}\rangle)^{4}\rangle}{\langle(\varepsilon_{n}-\langle\varepsilon_{n}\rangle)^{2}\rangle^{2}}-3\end{split}$
(4)

where angular brackets denote an average over events. Both skewness and
kurtosis vanish if the distribution is Guassian.

Figure 4: (Color online) The AMPT results on skewness and kurtosis of
$\varepsilon_{n}$ (n=2,3,4) fluctuation for 12C+197Au collisions at 200 GeV
with $\alpha$-clustering and W-S configurations of 12C, as a function of
$N_{track}$. Upper panels: Skewness of $\varepsilon_{n}$ fluctuation. Lower
panels: Kurtosis of $\varepsilon_{n}$ fluctuation.

In addition to the scaled variance, we further studied skewness and kurtosis
of the eccentricity fluctuation. Fig. 4 displays skewness and kurtosis of
$\varepsilon_{n}$ for W-S and triangle $\alpha$-clustering configurations of
12C. The skewness of eccentricity fluctuation $S_{\varepsilon_{n}}$ defined by
Eq.( 4) takes a similar trend as $R_{\varepsilon_{n}}$ for a different order
of harmonics. In comparison, $S_{\varepsilon_{2,4}}$ and
$K_{\varepsilon_{2,4}}$ for triangle and W-S configurations show tiny
multiplicity dependence. In addition, a significant difference in both
magnitude and trend can be seen for $S_{\varepsilon_{3}}$ and
$K_{\varepsilon_{3}}$ as a function of $N_{track}$ in the cases of different
initial configurations. The pronounced dependencies of $N_{track}$ of skewness
and kurtosis of $\varepsilon_{2}$ and $\varepsilon_{3}$ fluctuations imply
possible approaches for studying $\alpha$ clusterization of the collision
system by looking into the high order fluctuation properties.

###  III.2 Flow fluctuation in 12C+197Au collisions

It is predicted that the azimuthal anisotropy characterized by anisotropic
flow coefficients $v_{n}$ (n=2,3..) of final-state particle distribution is
correlated in magnitude and phase with initial eccentricity. Considering the
event-by-event fluctuation, calculation of $v_{n}$ can be done referred to the
participant plane angle Voloshin et al. (2008b) defined with the coordinate
information of initial participant partons:

$\psi_{n}\left\\{PP\right\\}=\frac{1}{n}\left[\arctan\frac{\left\langle{r^{2}\sin(n\varphi_{PP})}\right\rangle}{\left\langle{r^{2}\cos(n\varphi_{PP})}\right\rangle}+\pi\right],$
(5)

where $n$ denotes the $n$th-order participant plane, $r$ and $\varphi_{PP}$
are the position and azimuthal angle of each parton in AMPT initial stage and
the bracket denotes density weighted average. Flow coefficients $v_{n}$
calculated with respect to the participant plane $\psi_{n}$ are defined as

$v_{n}\left\\{PP\right\\}=\left\langle
cos[n(\phi-\psi_{n}\left\\{PP\right\\})]\right\rangle,$ (6)

where $\phi$ is azimuthal angle of final-state charged hadrons in the momentum
space, and the average $\langle\cdots\rangle$ denotes event average. Similar
to the definition of the eccentricity, this method for calculation of $v_{n}$
is referred to as participant plane $v_{n}$ which was widely used for flow
calculations in different models.

One can also characterize the different orders of azimuthal anisotropies with
Event plane (EP) method Poskanzer and Voloshin (1998); Voloshin et al.
(2008a). Unlike participant plane which is not accessible experimentally,
event plane can be reconstructed using final-state charged particles.
Definition of the $n$th-order event plane angle is in the form:

$\psi_{n}\left\\{EP\right\\}=\frac{1}{n}\arctan\frac{\left\langle{\omega\sin(n\phi)}\right\rangle}{\left\langle{\omega\cos(n\phi)}\right\rangle},$
(7)

where $\phi$ and $\omega$ are azimuthal angle and weight for the final
particle, respectively. Flow coefficients $v_{n}$ w.r.t the n-th order event
plane $\psi_{n}\left\\{EP\right\\}$ is defined as:

$v_{n}\left\\{EP\right\\}=\frac{\langle
cos(n[\phi-\psi_{n}\left\\{EP\right\\}])\rangle}{Res\left\\{\psi_{n}\left\\{EP\right\\}\right\\}},$
(8)

where $Res\left\\{\psi_{n}\left\\{EP\right\\}\right\\}$ is the resolution of
event plane angle and the brackets indicate average over particles.

Figure 5: (Color online) The AMPT results on elliptic and triangular flow
fluctuation as a function of $N_{track}$ for W-S and Triangle
$\alpha$-clustering configurations of 12C in 12C+197Au collisions at
$\sqrt{s_{NN}}=200GeV$. Results are shown from both event plane method and
participant plane methods. Left panel: $v_{2}$ fluctuation. Right panel:
$v_{3}$ fluctuation.

Ideal hydrodynamics predicts linear response to the initial eccentricities of
the final flow harmonics Qiu and Heinz (2011). Approximate proportionality of
the flow coefficients $v_{n}$ to the eccentricities $\varepsilon_{n}$ is
suggested to hold for n=2 and 3, $v_{n}$=$\kappa$$\varepsilon_{n}$ (n=2,3),
where $\kappa$ is the linear response coefficient. Assuming linear flow
response, flow fluctuation quantified with scaled standard deviation should be
approximately equal to the eccentricity fluctuation quantified in the same
way.

$\frac{\sigma_{v_{n}}}{\langle
v_{n}\rangle}\approx\frac{\sigma_{\varepsilon_{n}}}{\langle\varepsilon_{n}\rangle}$
(9)

This relation can be applied only on the premise that initial-final
correlation is dominated by linear response. Though fluctuations of $v_{n}$
mainly stem from the fluctuations of $\varepsilon_{n}$, during the source
evolution, non-linear responses may play an important role in the development
of final flow fluctuation with the presence of nonzero higher order effects
Yan et al. (2015); Renk and Niemi (2014); Retinskaya et al. (2014); Noronha-
Hostler et al. (2016); Bhalerao et al. (2011). Nevertheless, we can still
examine the flow fluctuation to see how sensitive it is in distinguishing
initial nuclear clustering structure.

Defined in the same way as Eq.( 3), Fig. 5 shows the simulation results of the
scaled variance of the elliptic and triangular flow fluctuation of the charged
hadrons at mid-rapidity (-1.0 $<$ y $<$ 1.0) for 12C+197Au collisions at
$\sqrt{s_{NN}}$ = 200 GeV. Comparisons are made between $\alpha$-clustered
case and Woods-Saxon case where the nucleons are distributed without
clusterization. $R_{v_{n}}\left\\{EP\right\\}$ which quantifies
$v_{n}\left\\{EP\right\\}$ fluctuation are shown in addition to
$v_{n}\left\\{PP\right\\}$ fluctuation $R_{v_{n}}\left\\{PP\right\\}$ based on
the participant plane method. $R_{v_{n}}\left\\{EP\right\\}$ ($n$=2,3) show
similar behavior as $R_{v_{n}}\left\\{PP\right\\}$ but comparably much higher
for low-multiplicity collisions owing to the smearing effect brought by the
source evolution.

It is found that the scaled variances of both $v_{2}$ and $v_{3}$ present a
monotonous decreasing trend at low $N_{track}$ ($N_{track}$$<$100) and change
slowly towards higher $N_{track}$. $R_{v_{2}}$ is insensitive for
distinguishing Traingle and W-S configurations as they are close in magnitude.
A magnitude ordering can be seen for $R_{v_{3}}$ resulted from W-S and
clustering configurations. It is observed that $R_{v_{3}}$ for Triangle
configuration is lower than W-S configuration similar to the behavior of
$R\varepsilon_{3}$ which is consistent with the picture of the flow response.

Experimentally, one could take $R_{v_{n}}\left\\{EP\right\\}$ as a probe of
$\alpha$-clustering structure in 12C. Similar to probing nuclear clustering
structure with the ratio of flow harmonics, by choosing a isobar nucleus close
to 12C with non-exotic structure colliding against 197Au as a reference,
quantitative difference of $R_{v_{n}}\left\\{EP\right\\}$ between
$\alpha$-clustering collision system and unclustered reference system can
serve as a good probe to distinguish the clustering nuclear structure.

Figure 6: (Color online) Ratio of the relative flow fluctuation
$R_{v_{3}}$/$R_{v_{2}}$ as a function of $N_{track}$ for 12C+197Au collisions
with $\alpha$-clustered 12C in different structures. Flow fluctuations are
calculated based on the event-plane method.

We further study the ratio of flow fluctuations. Fig. 6 displays the ratio of
the relative fluctuation of $v_{3}$ and $v_{2}$. Significant difference in
magnitude can be seen for the $N_{track}$ dependence of
$R_{v_{3}}$/$R_{v_{2}}$ for W-S and Triangle $\alpha$-clustering
configurations. It is found that both Triangle and W-S configurations show
monotonic $N_{track}$ dependence. The ratio is more sensitive for
distinguishing Triangle and W-S configurations of 12C at large $N_{track}$.

Figure 7: (Color online) Skewness and kurtosis of $v_{2}$ and $v_{3}$
fluctuation as a function of $N_{track}$ for $\alpha$-clustered 12C+197Au
collisions at 200 GeV. Comparisons are made between results from event plane
method and participant plane method. Upper panels: Skewness of $v_{2}$ and
$v_{3}$ fluctuation. Lower panels: Kurtosis of $v_{2}$ and $v_{3}$
fluctuation.

Assuming linear flow response, the approximation relation shown by Eq.( 9)
also applies to the kurtosis and skewness of flow fluctuations. We study the
skewness ($S_{v_{n}}$) and kurtosis ($K_{v_{n}}$) of $v_{n}$ ($n$=2,3)
fluctuation in a similar way as corresponding $\varepsilon_{n}$ fluctuation as
approaches for investigating the initial $\alpha$-clustering effect. Results
from both the participant plane method and the event plane method are shown in
comparison in Fig. 7. Different from the eccentricity fluctuation, one can
find that both skewness and kurtosis of the flow harmonic are insensitive for
distinguishing the unclustered W-S configuration and the triangular
$\alpha$-clustering configuration. It could be from the source evolution that
in the final-state of AMPT the $v_{n}$ ($n$=2,3) distributions are nearly
Gaussian, thus the skewness and kurtosis of elliptic and triangular flow
fluctuations especially at large $N_{track}$ are almost consistent with zero
within statistical uncertainties.

###  III.3 Correlation of the initial eccentricity with final flow harmonic

Impressive progresses have been made in studying the sensitivity of flow
response to the initial geometry in relativistic heavy-ion collisions Petersen
et al. (2012); Niemi et al. (2013). We understand that elliptic flow $v_{2}$
and triangular flow $v_{3}$ are driven mainly by the linear response to the
initially produced ellipticity and triangularity of the source geometry. For
asymmetric colliding systems involving $\alpha$-clustered nucleus, the
conversion of the initial geometry to the final flow could become more
complicated. Quantitative study of the initial-final correlation in event-by-
event basis in model simulation is important for understanding the sensitivity
of flow and flow fluctuation in probing clustering configurations.

In particular, we study the correlations of the initial eccentricities with
final flow harmonics for 12C+197Au collisions with $\alpha$-clustered 12C
configurations. Pearson coefficient is used to quantify the strength of the
correlation. We define the coefficient in the following form which takes into
account both the magnitude and the angle:

$C_{v_{n},\varepsilon_{n}}=\frac{\langle
v_{n}\varepsilon_{n}cos(n[\psi-\psi_{PP}])\rangle}{\sqrt{\langle|\varepsilon_{n}|^{2}\rangle\langle|v_{n}|^{2}\rangle}},$
(10)

where $\psi_{PP}$ is the participant plane angle, $\psi$ is the phase of the
flow coefficient $v_{n}$ (e.g. $\psi$ is the event-plane angle if $v_{n}$ is
calculated based on event plane method). As $C_{v_{n},\varepsilon_{n}}$
approaches one, good linear correlation can be expected, whereas
$C_{v_{n},\varepsilon_{n}}$ approaches zero other contributions beyond linear
response could contribute. We remark here the Pearson coefficent can be used
to quantify correlations even non-linear effect contributes.

Figure 8: (Color online) 2-D plots illustrating the event-by-event correlation
of $v_{n}$ and $\varepsilon_{n}$ for 12C+197Au collisions with W-S and
Triangle $\alpha$-clustering configuration of 12C. The values of the
correlation coefficients $C_{v_{n},\varepsilon_{n}}$ are shown in the plots.
Upper panels: $v_{2}$ vs $\varepsilon_{2}$ for 80$<$$N_{track}$$<$100\. Lower
panels: $v_{3}$ vs $\varepsilon_{3}$ for 80$<$$N_{track}$$<$100\. Red and blue
data points present the $\varepsilon_{n}$-bin averaged profile.

The two-dimensional plots in Fig. 8 show the correlations between $v_{n}$ and
$\varepsilon_{n}$ calculated using participant plane method. As can be seen in
the figures, for collision events with 80$<$$N_{track}$$<$100, the $v_{2}$ and
$v_{3}$ coefficients display a strong linear correlation to their
corresponding initial-state coefficients for all the $\alpha$-clustering cases
considered.

Figure 9: (Color online) Correlation coefficients $C_{v_{n},\varepsilon_{n}}$
(n=2,3) as a function of $N_{track}$ for 12C+197Au collisions with different
12C $\alpha$-clustering configurations.

Fig. 9 compares the model predictions of the $N_{track}$ dependence of the
correlation function defined by Eq.( 10) for 12C+197Au collisions with W-S and
$\alpha$-clustered 12C configurations. $C_{v_{n},\varepsilon_{n}}$ ($n$=2,3)
are generally seen to follow a smooth increasing trend as a function of
$N_{track}$ indicating stronger linear $v_{n}$-$\varepsilon_{n}$ correlation
at larger $N_{track}$. In comparison, W-S results are comparable with triangle
configuration for the correlation between $v_{2}$ and $\varepsilon_{2}$
whereas triangle configuration presents stronger $v_{3}$-$\varepsilon_{3}$
correlation than W-S especially at large $N_{track}$.

##  IV Summary

In summary, we studied the $\alpha$-clustering effects on the initial
eccentricity and final flow fluctuations in $\alpha$-clustered 12C+197Au
collisions at center-of-mass energy of 200 GeV using a multi-phase transport
model.

Event-by-event fluctuations of the initial eccentricity $\varepsilon_{n}$ and
final anisotropic flow $v_{n}$ are characterized by scaled standard variance,
skewness and kurtosis. Differences in the multiplicity dependence of the flow
fluctuation are observed for collision systems with W-S and Triangle
$\alpha$-clustering configurations of 12C. The triangular flow fluctuation is
shown to have particular sensitivity in distinguishing the triangle
$\alpha$-clustering configuration which is consistent with the picture of the
flow response to the initial geometry. The ratio of the triangular flow
fluctuation over elliptic flow fluctuation decreases with increasing
$N_{track}$ for both W-S and triangle configurations but show difference in
magnitude. We also note that the differences in flow fluctuation between
clustered and unclustered W-S structures are more significant in high-
multiplicity collisions. Furthermore, correlations between flow harmonic and
initial eccentricity are investigated. The correlation functions for 12C+197Au
collisions generally show an increasing trend as a function of $N_{track}$,
indicating stronger linear correlations at larger $N_{track}$. In the picture
of the flow response to the initial geometry, it is expected that harmonic
flow or flow fluctuation can better probe the $\alpha$-clustering structure in
higher multiplicity collision events.

Experimentally, as the ground state of 12C could exist as a configuration
mixing state in reality, the actual significance of the fluctuation observable
in distinguishing clustered configurations could be lower than we obtained in
the ideal cases in our study. Nevertheless, because of the particular
sensitivity of the flow fluctuation in reflecting the clustering structure,
one could take flow fluctuation in addition to the flow measurement as an
effective probe for the $\alpha$-clustering structure.

As in high energy nuclear collisions, study of the event-by-event anisotropy
fluctuation is crucial for understanding not only system initial conditions
but also the afterburner evolution properties, future studies by looking into
the flow fluctuation in other light-heavy collision systems (e.g. 16O+197Au,
16O+208Pb) with $\alpha$-clustering configurations will be promising for
providing more important information about the nuclear clustering effect.

## Acknowledgements

This work is supported in part by the National Natural Science Foundation of
China under contract Nos. 11890714, 11421505, 11905034, 11925502, 11935001,
11961141003, the Key Research Program of Frontier Sciences of the CAS under
Grant No. QYZDJ-SSW-SLH002 and the Key Research Program of the CAS under Grant
No. XDPB09.

## References

  * (1)
  * Ollitrault (1992) J.-Y. Ollitrault, Phys. Rev. D 46, 229 (1992), URL http://link.aps.org/doi/10.1103/PhysRevD.46.229.
  * Voloshin et al. (2008a) S. A. Voloshin, A. M. Poskanzer, and R. Snellings, Relativistic Heavy Ion Physics 46 (2008a).
  * Kolb et al. (1999) P. F. Kolb, J. Sollfrank, and U. W. Heinz, Phys. Lett. B 459, 667 (1999).
  * Ackermann et al. (2001) K. H. Ackermann et al. (STAR Collaboration), Phys. Rev. Lett. 86, 402 (2001), URL http://link.aps.org/doi/10.1103/PhysRevLett.86.402.
  * Teaney et al. (2001) D. Teaney, J. Lauret, and E. V. Shuryak, Phys. Rev. Lett. 86, 4783 (2001), URL http://link.aps.org/doi/10.1103/PhysRevLett.86.4783.
  * Romatschke and Romatschke (2007) P. Romatschke and U. Romatschke, Phys. Rev. Lett. 99, 172301 (2007), URL http://link.aps.org/doi/10.1103/PhysRevLett.99.172301.
  * Luo and Xu (2017) X. Luo and N. Xu, Nucl. Sci. Tech. 28, 112 (2017).
  * Adams et al. (2004) J. Adams, C. Adler, et al. (STAR Collaboration), Phys. Rev. Lett. 92, 062301 (2004), URL http://link.aps.org/doi/10.1103/PhysRevLett.92.062301.
  * Adamczyk et al. (2013) L. Adamczyk et al. (STAR Collaboration), Phys. Rev. C88, 014904 (2013).
  * Adare et al. (2011) A. Adare et al. (PHENIX Collaboration), Phys. Rev. Lett. 107, 252301 (2011).
  * Heinz (2005) U. Heinz, Journal of Physics G: Nuclear and Particle Physics 31, S717 (2005), URL http://stacks.iop.org/0954-3899/31/i=6/a=012.
  * Gale et al. (2013) C. Gale, S. Jeon, and B. Schenke, Int. J. Mod. Phys. A28, 1340011 (2013).
  * Qiu and Heinz (2012) Z. Qiu and U. Heinz, AIP Conf. Proc. 1441, 774 (2012), eprint 1108.1714.
  * Song et al. (2011) H. Song, S. A. Bass, U. Heinz, T. Hirano, and C. Shen, Phys. Rev. Lett. 106, 192301 (2011), URL http://link.aps.org/doi/10.1103/PhysRevLett.106.192301.
  * Song et al. (2017) H. Song, Y. Zhou, and K. Gajdošová, Nuclear Science and Techniques 28, 99 (2017).
  * Chen et al. (2004) L.-W. Chen, C. M. Ko, and Z.-W. Lin, Phys. Rev. C 69, 031901(R) (2004), URL http://link.aps.org/doi/10.1103/PhysRevC.69.031901.
  * Huang et al. (2020) S. Huang, Z. Chen, W. Li, and J. Jia, Phys. Rev. C 101, 021901(R) (2020).
  * Alver et al. (2007) B. Alver, B. B. Back, et al., Phys. Rev. Lett. 98, 242302 (2007), URL http://link.aps.org/doi/10.1103/PhysRevLett.98.242302.
  * Voloshin (2006) S. A. Voloshin, arXiv preprint nucl-th/0606022 (2006).
  * Voloshin et al. (2008b) S. A. Voloshin, A. M. Poskanzer, A. Tang, and G. Wang, Phys. Lett. B659, 537 (2008b).
  * Lacey et al. (2011) R. A. Lacey, R. Wei, N. N. Ajitanand, and A. Taranenko, Phys. Rev. C83, 044902 (2011).
  * Lacey et al. (2014) R. A. Lacey, A. Taranenko, J. Jia, D. Reynolds, N. N. Ajitanand, J. M. Alexander, Y. Gu, and A. Mwai, Phys. Rev. Lett. 112, 082302 (2014), URL http://link.aps.org/doi/10.1103/PhysRevLett.112.082302.
  * Alver and Roland (2010) B. Alver and G. Roland, Phys. Rev. C 81, 054905 (2010), URL http://link.aps.org/doi/10.1103/PhysRevC.81.054905.
  * Derradi de Souza et al. (2012) R.D. deSouza, J. Takahashi, T. Kodama, and P. Sorensen, Phys. Rev. C85, 054909 (2012).
  * Ma and Wang (2011) G.-L. Ma and X.-N. Wang, Phys. Rev. Lett. 106, 162301 (2011), URL http://link.aps.org/doi/10.1103/PhysRevLett.106.162301.
  * Ma et al. (2014) L. Ma, G. L. Ma, and Y. G. Ma, Phys. Rev. C 89, 044907 (2014), URL http://link.aps.org/doi/10.1103/PhysRevC.89.044907.
  * Han et al. (2011) L. X. Han, G. L. Ma, Y. G. Ma, X. Z. Cai, J. H. Chen, S. Zhang, and C. Zhong, Phys. Rev. C 84, 064907 (2011), URL http://link.aps.org/doi/10.1103/PhysRevC.84.064907.
  * Wang et al. (2013) J. Wang, Y. G. Ma, G. Q. Zhang, D. Q Fang, L. X. Han, and W. Q. Shen, Nucl. Sci. Tech. 24, 30501 (2013), URL https://inspirehep.net/files/7aa9f962f34f205bbbccc43d41847509.
  * Alver et al. (2010) B. Alver et al., Phys. Rev. Lett. 104, 142301 (2010), URL http://link.aps.org/doi/10.1103/PhysRevLett.104.142301.
  * Agakishiev et al. (2012) G. Agakishiev et al. (STAR Collaboration), Phys. Rev. C86, 014904 (2012).
  * Sorensen (2007) P. Sorensen, Journal of Physics G: Nuclear and Particle Physics 34, S897 (2007), URL http://stacks.iop.org/0954-3899/34/i=8/a=S121.
  * Margutti et al. (2019) J. Margutti et al., Nuclear Physics A 982, 367 (2019).
  * Andrade et al. (2006) R. Andrade, F. Grassi, Y. Hama, T. Kodama, and O. Socolowski, Phys. Rev. Lett. 97, 202302 (2006), URL http://link.aps.org/doi/10.1103/PhysRevLett.97.202302.
  * Petersen et al. (2010) H. Petersen, G.-Y. Qin, S. A. Bass, and B. Müller, Phys. Rev. C 82, 041901(R) (2010), URL http://link.aps.org/doi/10.1103/PhysRevC.82.041901.
  * Ma et al. (2016) L. Ma, G. L. Ma, and Y. G. Ma, Phys. Rev. C 94, 044915 (2016).
  * Brink et al. (1970) D. Brink, H. Friedrich, A. Weiguny, and C. Wong, Phys. Lett. B 33, 143 (1970).
  * von Oertzen et al. (2006) W. von Oertzen, M. Freer, and Y. Kanada-En’yo, Physics Reports 432, 43 (2006), ISSN 0370-1573, URL http://www.sciencedirect.com/science/article/pii/S0370157306002626.
  * Freer et al. (2018) M. Freer, H. Horiuchi, Y. Kanada-En’yo, D. Lee, and Ulf-G. Meißner, Rev. Mod. Phys. 90, 035004 (2018).
  * He et al. (2014) W. B. He, Y. G. Ma, X. G. Cao, X. Z. Cai, G. Q. Zhang, et al., Phys. Rev. Lett. 113, 032506 (2014).
  * He et al. (2016) W. B. He, Y. G. Ma, X. G. Cao, X. Z. Cai, G. Q. Zhang, et al., Phys. Rev. C 94, 014301 (2016).
  * Kanada-En’yo et al. (2015) Y. Kanada-En’yo et al., Nucl. Sci. Tech. 26, 20501 (2015).
  * Liu and Ye (2018) Y. Liu and Y.-L. Ye, Nucl. Sci. Tech. 29, 184 (2018).
  * Huang et al. (2017) B. S. Huang, Y. G. Ma, and W. B. He, Phys. Rev. C 95, 034606 (2017).
  * Huang and Ma (2020) B. S. Huang and Y. G. Ma, Phys. Rev. C 101, 034615 (2020).
  * Broniowski and Ruiz Arriola (2014) W. Broniowski and E. Ruiz Arriola, Phys. Rev. Lett. 112, 112501 (2014), URL https://link.aps.org/doi/10.1103/PhysRevLett.112.112501.
  * Zhang et al. (2017) S. Zhang, Y. G. Ma, J. H. Chen, W. B. He, and C. Zhong, Phys. Rev. C 95, 064904 (2017), URL https://link.aps.org/doi/10.1103/PhysRevC.95.064904.
  * Zhang et al. (2018) S. Zhang, Y. G. Ma, J. H. Chen, W. B. He, and C. Zhong, The European Physical Journal A 54, 161 (2018).
  * Guo et al. (2019) C.-C. Guo, Y.-G. Ma, Z.-D. An, and B.-S. Huang, Phys. Rev. C 99, 044607 (2019), URL https://link.aps.org/doi/10.1103/PhysRevC.99.044607.
  * Xu et al. (2018) Z.-W. Xu, S. Zhang, Y.-G. Ma, J.-H. Chen, and C. Zhong, Nucl. Sci. Tech. 29, 186 (2018).
  * Guo et al. (2017) C.-C. Guo, W.-B. He, and Y.-G. Ma, Chinese Physics Letters 34, 092101 (2017).
  * Bożek et al. (2014) P. Bozek, W. Broniowski, E. R. Arriola, and M. Rybczyński, Phys. Rev. C 90, 064902 (2014), URL https://link.aps.org/doi/10.1103/PhysRevC.90.064902.
  * Rybczyński et al. (2018) M. Rybczyński, M. Piotrowska, and W. Broniowski, Phys. Rev. C 97, 034912 (2018).
  * Zhang et al. (2000) B. Zhang, C. M. Ko, B.-A. Li, and Z. Lin, Phys. Rev. C 61, 067901 (2000), URL http://link.aps.org/doi/10.1103/PhysRevC.61.067901.
  * Lin et al. (2005) Z.-W. Lin, C. M. Ko, B.-A. Li, B. Zhang, and S. Pal, Phys. Rev. C 72, 064901 (2005), URL http://link.aps.org/doi/10.1103/PhysRevC.72.064901.
  * Wang and Gyulassy (1991) X.-N. Wang and M. Gyulassy, Phys. Rev. D 44, 3501 (1991), URL http://link.aps.org/doi/10.1103/PhysRevD.44.3501.
  * Zhang (1998) B. Zhang, Comput. Phys. Commun. 109, 193 (1998).
  * Li and Ko (1995) B.-A. Li and C. M. Ko, Phys. Rev. C 52, 2037 (1995), URL https://link.aps.org/doi/10.1103/PhysRevC.52.2037.
  * Xu et al. (2016) Y.-F. Xu, Y.-J. Ye, J.-H. Chen, Y.-G. Ma, S. Zhang, and C. Zhong, Nucl. Sci. Tech. 27, 87 (2016).
  * Zhou et al. (2016) Y. Zhou, K. Xiao, Z. Feng, F. Liu, and R. Snellings, Phys. Rev. C93, 034909 (2016).
  * Jin et al. (2018) X.-H. Jin, J.-H. Chen, Y.-G. Ma, S. Zhang, C.-J. Zhang, and C. Zhong, Nucl. Sci. Tech. 29, 54 (2018).
  * Wang et al. (2019) H. Wang, J.-H. Chen, Y.-G. Ma, and S. Zhang, Nucl. Sci. Tech. 30, 185 (2019).
  * Maruyama et al. (1996) T. Maruyama, K. Niita, and A. Iwamoto, Phys. Rev. C 53, 297 (1996).
  * Drescher and Nara (2007) H.-J. Drescher and Y. Nara, Phys. Rev. C 76, 041903(R) (2007), URL http://link.aps.org/doi/10.1103/PhysRevC.76.041903.
  * Broniowski et al. (2007) W. Broniowski, P. Bozek, and M. Rybczyński, Phys. Rev. C 76, 054905 (2007), URL http://link.aps.org/doi/10.1103/PhysRevC.76.054905.
  * Miller and Snellings (2003) M. Miller and R. Snellings (2003), eprint nucl-ex/0312008.
  * Qiu and Heinz (2011) Z. Qiu and U. W. Heinz, Phys. Rev. C84, 024911 (2011).
  * Giacalone et al. (2017) G. Giacalone, L. Yan, J. Noronha-Hostler, and J.-Y. Ollitrault, Phys. Rev. C 95, 014913 (2017).
  * Bhalerao et al. (2019) R. S. Bhalerao, G. Giacalone, and J.-Y. Ollitrault, Phys. Rev. C 99, 014907 (2019).
  * Poskanzer and Voloshin (1998) A. M. Poskanzer and S. A. Voloshin, Phys. Rev. C 58, 1671 (1998).
  * Yan et al. (2015) L. Yan, J.-Y. Ollitrault, and A. M. Poskanzer, Phys. Lett. B 742, 290 (2015).
  * Renk and Niemi (2014) T. Renk and H. Niemi, Phys. Rev. C 89, 064907 (2014), URL https://link.aps.org/doi/10.1103/PhysRevC.89.064907.
  * Retinskaya et al. (2014) E. Retinskaya, M. Luzum, and J.-Y. Ollitrault, Phys. Rev. C 89, 014902 (2014), URL https://link.aps.org/doi/10.1103/PhysRevC.89.014902.
  * Noronha-Hostler et al. (2016) J. Noronha-Hostler, L. Yan, F. G. Gardim, and J.-Y. Ollitrault, Phys. Rev. C 93, 014909 (2016).
  * Bhalerao et al. (2011) R. S. Bhalerao, M. Luzum, and J.-Y. Ollitrault, Phys. Rev. C 84, 054901 (2011), URL http://link.aps.org/doi/10.1103/PhysRevC.84.054901.
  * Petersen et al. (2012) H. Petersen, R. La Placa, and S. A. Bass, J. Phys. G39, 055102 (2012).
  * Niemi et al. (2013) H. Niemi, G.S. Denicol, H. Holopainen, and P. Huovinen, Phys. Rev. C 87, 054901 (2013).


# The Brezis-Nirenberg problem in 4D

Angela Pistoia  Angela Pistoia, Dipartimento di Scienze di Base e Applicate
per l’Ingegneria, Sapienza Università di Roma, Via Antonio Scarpa 10, 00161
Roma (Italy<EMAIL_ADDRESS>and  Serena Rocci  Serena Rocci,
Dipartimento di Scienze di Base e Applicate per l’Ingegneria, Sapienza
Università di Roma, Via Antonio Scarpa 10, 00161 Roma (Italy)
<EMAIL_ADDRESS>
###### Abstract.

We address the existence of blowing-up solutions for the Brezis-Nirenberg
problem in 4D.

Dedicated to Yihong Du’s 60th birthday

###### Key words and phrases:

Brezis-Nirenberg problem, blow-up solutions, Ljapunov-Schmidt construction

###### 1991 Mathematics Subject Classification:

Primary: 35J25. Secondary: 35B09

The authors are partially supported by the group GNAMPA of the Istituto
Nazionale di Alta Matematica (INdAM)

##  1\. Introduction

The problem

$-\Delta u=|u|^{4\over n-2}u+\lambda Vu\ \hbox{in}\ \Omega,\ u=0\ \hbox{on}\
\partial\Omega$ (1.1)

where $\Omega$ is a bounded regular domain in $\mathbb{R}^{n}$,
$\lambda\in\mathbb{R}$ and $V\in C^{0}(\overline{\Omega})$ was introduced by
Brezis and Nirenberg in the pivotal paper [4], where they address the
existence of positive solutions in the autonomous case, i.e. the potential $V$
is constant. Since then, a huge aumont of work has been done. In the following
we will make a brief history highlighting the results which are much closer to
the problem we wish to study in the present paper.

It is well known that existence and multiplicity of positive and sign-changing
solutions to (1.1) is strictly affected by the geometry of the domain, the
dimension of the euclidean space where the domain lies and the values of the
parameter $\lambda$. For example, if $\lambda\leq 0,$ $V$ is constant and
$\Omega$ is a starshaped domain, problem (1.1) does not have any solutions.
A particular feature of problem (1.1), due to the critical behaviour of the
non-linearity which appears on the R.H.S., is the possible existence of
solutions which blow-up at one or more points in the domain as the parameter
$\lambda$ approaches $0$ in dimensione $n\geq 4$ or a positive number
$\lambda_{*}$ in dimension $n=3.$ The description of the profile of the
positive blowing-up solutions has been the subject of a wide literature.
Starting from the pioneering paper by Brezis and Peletier [5] where the
authors consider the radial case, a lot of results have been obtained
concerning the asymptotic profiles of solutions to (1.1). When $n\geq 4$, the
first result is due to Han [9] and Rey [19] who study the profile of one-peak
solutions in a general domain when the potential V is constant, while the most
recent one has been obtained by Konig and Laurin [10] who give an accurate
description of the blow-up profile of a solution with multiple blow-up points
in the non-autonomous case. The 3-dimensional case has been firstly faced by
Druet [7] in the case of one blow-up point and very recently by Konig and
Laurain [11] in the presence of multiple blow-up points.
In particular, the asymptotic analysys of blowing-up solutions ensures that
the blow-up points are nothing but the critical points of a suitable function
which involves the function V, the Green’s function of $-\Delta$ in $\Omega$
with Dirichlet boundary condition and the Robin’s function.
While the behaviour of blowing-up positive solutions is by now clear, the
profile of the sign-changing solutions is far being completely understood. As
far as we know the complete scenario of sign-changing solutions has been
obtained only in the radial case by Esposito, Ghossoub, Pistoia and Vaira in
[8] in dimensions $n\geq 7$ and by Amadori, Gladiali, Grossi, Pistoia and
Varia [1] in lower dimensions.

A parallel and prolific field of research consists in finding blowing-up
solutions whose profile is predicted by the asymptotic analysis developed in
the previous papers. The first result is due to Rey [19] who build solutions
blowing-up at a non-degenerate critical point of the Robin function when $V$
is constant and $n\geq 5.$ Successively, positive solutions blowing-up at
multiple points have been constructed by Musso and Pistoia [15]. If $n\geq 5$
in the non-autonomous case, solutions with one blow-up point have been found
by Pistoia and Molle [13], while Micheletti and Pistoia [12] build positive
and sign-changing solutions with multiple blow-up points.
While all the blow-up points of positive solutions are always isolated and
simple, in high dimensions $n\geq 7$ sign-changing solutions can have towering
blow-up points (see Premoselli [18] and Morabito-Pistoia-Vaira [14]) and
clustering blow-up points (see Pistoia Vaira [17]).
In the 3-dimensional case, when $V$ is constant, single and multiple blowing-
up solutions have been built by del Pino, Dolbeault and Musso [6] and Musso
and Salazar [16], respectively.

Even if the problem has been widely studied in the last decades and a huge
number of results has been produced, many problems are still open.

Here we focus on the problem (1.1) in dimension $n=4$:

$\begin{cases}-\Delta u=u^{3}+\varepsilon Vu&\mbox{in }\Omega\\\ u=0&\mbox{on
}\partial\Omega\end{cases}$ (1.2)

when $\Omega$ is a bounded regular domain in ${\mathbb{R}^{4}}$, $V\in
C^{1}(\Omega)\cap C^{0}(\overline{\Omega})$ and we prove the existence of
solutions which blows-up at a single point as the parameter $\epsilon$
approaches zero.
More precisely, if $\tau$ denotes the Robin function, our main result reads as
follows.

######  Theorem 1.1.

Let ${\xi}_{0}\in\Omega$ be a non-degenerate critical point of
$f(\xi):={\tau(\xi)\over V(\xi)}$ with $V(\xi_{0})>0.$
If $\epsilon$ is small enough there exists a solution of problem (1.2) which
blows-up at the point $\xi_{0}$ as $\varepsilon\to 0$.

The proof relies on a classical Ljapunov-Schmidt procedure. The setting of the
problem (see Section 2) and the reduction process (see Section 3) can be
carried out as usual. However, the last step in the procedure needs new ideas.
In fact, the rate of the error term is not small enough to argue as in the
higher dimensional case and the reduced problem is solved using some local
Pohozaev identities (see Section 4).

Remarkably, differently from solutions with one blow-up point, in the case of
multiple blow-up points it is harder to derive the concentration speed.
Indeed, it seems that they appear at the second order expansion and so to
catch them a more accurate description of the ansatz is needed. This will be
the topic of a forthcoming paper in collaboration with Monica Musso.

##  2\. Setting of the problem

###  2.1. The bubbles

All the positive solutions to the limit problem

$-\Delta U=U^{3}\mbox{ in }{\mathbb{R}^{4}}$

are the so called bubbles (see [2, 20])

$U_{\delta,\xi}(x)=\frac{1}{\delta}U\left(x-\xi\over\delta\right),\
x,\xi\in\mathbb{R}^{4},\ \delta>0$

where

$U(x):=\mathfrak{c}\frac{1}{1+|x|^{2}},\ \mathfrak{c}:=2\sqrt{2}.$

It is useful to introduce the projection of the bubble $U_{\delta,\xi}$ onto
$H^{1}_{0}(\Omega)$, namely the solution of the problem

$-\Delta PU_{\delta,\xi}=-\Delta U_{\delta,\xi}\ \mbox{ in }\Omega,\
PU_{\delta,\xi}=0\ \mbox{ on }\partial\Omega.$

Let $G$ be the Green function of $-\Delta$ on $H_{0}^{1}(\Omega)$ and $H$ be
its regular part, i.e.

$G(x,y)=\frac{1}{2\omega}\frac{1}{|x-y|^{2}}-H(x,y)$

where $\omega$ denotes the measure of the unit sphere
$S^{3}\subset{\mathbb{R}^{4}}$. The Robin function is defined as
$\tau(x)=H(x,x).$

It is well known that

$PU_{\delta,\xi}(x)=U_{\delta,\xi}(x)-\mathfrak{C}\delta
H(x,\xi)+\mathcal{O}(\delta^{3}),\ \mathfrak{C}:=2\mathfrak{c}\omega$

uniformly with respect to $x\in\Omega$ and $\xi$ in compact sets of $\Omega$
and

$PU_{\delta,\xi}(x)=\mathfrak{C}\delta G(x,\xi)+\mathcal{O}(\delta^{3})$

uniformly with respect to $x$ in compact sets of $\Omega\setminus\\{\xi\\}$
and $\xi$ in compact sets of $\Omega$

###  2.2. Some background material

Let $H^{1}_{0}(\Omega)$ be the Hilbert space equipped with the usual inner
product and the usual norm

$\langle u,v\rangle=\int_{\Omega}\nabla u\cdot\nabla v\ \hbox{and}\
\|u\|:=\|u\|_{H^{1}_{0}(\Omega)}=\left(\int_{\Omega}|\nabla
u|^{2}\right)^{1/2}.$

For $r\in[1,+\infty)$ the space $L^{r}(\Omega)$ is also equipped with the
standard norm

$\|u\|_{r}=\left(\int_{\Omega}|u|^{r}\right)^{\frac{1}{r}}.$

Now, let us introduce ${\mathtt{i}}^{*}:L^{4/3}(\Omega)\to H^{1}_{0}(\Omega)$
as the adjoint operator of the embedding
$\mathtt{i}:H^{1}_{0}(\Omega)\hookrightarrow L^{4}(\Omega)$, i.e.
$u={\mathtt{i}}^{*}(f)$ if and only if

$\langle u,\varphi\rangle=\int_{\Omega}f(x)\varphi(x)dx\mbox{ for all
}\varphi\in H^{1}_{0}(\Omega)$

or equivalently

$-\Delta u=f\mbox{ in }\Omega,\ u=0\mbox{ on }\partial\Omega$

The operator ${\mathtt{i}}^{*}:L^{4/3}(\Omega)\to H^{1}_{0}(\Omega)$ is
continuous as

$\|{\mathtt{i}}^{*}(f)\|_{H^{1}_{0}(\Omega)}\leqslant S^{-1}\|f\|_{4/3}$

where $S$ is the best constant for the Sobolev embedding.

Therefore, the problem (1.2) can be rewritten as

$u={\mathtt{i}}^{*}(u^{3}+\varepsilon Vu),\quad u\in H^{1}_{0}(\Omega).$ (2.1)

###  2.3. The ansatz

Let $k\geq 1$ be a fixed integer. We look for a solution to (1.2) of the form

$u=U_{\delta,\xi}+\phi,$ (2.2)

where $U_{\delta,\xi}=PU_{\delta,\xi}$ whose blow-up rate is
$\delta=\delta(\varepsilon)\to 0$ and blow-up point is
$\xi=\xi(\varepsilon)\in\Omega$. Moreover, the lower order term $\phi$
satisfies a set of orthogonal conditions. More precisely, let us consider the
linear problem

$-\Delta\psi=3U^{2}\psi\mbox{ in }{\mathbb{R}^{4}}.$

It is well known [3] that the set of solutions is a $5-$dimensional linear
space spanned by

$\psi^{0}(x):=U(x)+\frac{1}{2}\nabla U(x)\cdot
x=-\mathfrak{c}\frac{1-|x|^{2}}{(1+|x|^{2})^{2}}$

and

$\psi^{j}(x):=\frac{\partial U}{\partial
x_{j}}(x)=-2\mathfrak{c}\frac{x_{j}}{(1+|x|^{2})^{2}},\ j=1,\dots,4.$

For $j=0,1,2,3,4$ we set

$\psi_{\delta,\xi}^{j}(x)=\frac{1}{\delta}\psi^{j}\left(x-\xi\over\delta\right).$

We introduce their projections
$P\psi_{\delta,\xi}^{j}={\mathtt{i}}^{*}(3U_{\delta,\xi}^{2}\psi_{\delta,\xi}^{j})$
onto $H^{1}_{0}(\Omega),$ namely the solutions to the problem

$-\Delta
P\psi_{\delta,\xi}^{j}=-\Delta\psi_{\delta,\xi}^{j}=3U_{\delta,\xi}^{2}\psi_{\delta,\xi}^{j}\
\hbox{in}\ \Omega,\ P\psi_{\delta,\xi}^{j}=0\ \hbox{in}\ \partial\Omega.$

Finally, we introduce the linear space

$K_{{\delta},{\xi}}=\operatorname{span}\\{P\psi_{\delta,\xi}^{j}\ |\
j=0,\cdots,4\\}$

and its orthogonal space

$K_{\delta,\xi}^{\perp}=\\{\phi\in
H_{0}^{1}(\Omega):\langle\phi,P\psi_{\delta,\xi}^{j}\rangle=0,\mbox{ for all
}j=0,\cdots,4\\}.$

The function $\phi$ belongs to the space $K_{\delta,\xi}^{\perp}.$

It is useful to remind the well known properties

$P\psi^{j}_{\delta,\xi}=\psi^{j}_{\delta,\xi}-\delta^{2}\mathfrak{C}\partial_{\xi_{j}}H(x,\xi)+\mathcal{O}(\delta^{3}),\quad
j=1,\ldots,4$

uniformly for $x\in\Omega$ and $\xi$ in compact sets of $\Omega$ and

$P\psi^{0}_{\delta,\xi}=\psi^{0}_{\delta,\xi}-\delta\mathfrak{C}H(x,\xi)+\mathcal{O}(\delta^{2})$

uniformly for $x$ in compact sets of $\Omega\setminus\\{\xi\\}$ and $\xi$ in
compact sets of $\Omega.$

###  2.4. An equivalent system

Let us introduce the linear projection $\Pi_{{\delta},{\xi}}:K\to K$ and
$\Pi_{{\delta},{\xi}}^{\perp}:K^{\perp}\to K^{\perp}$ which are defined by

$\Pi_{{\delta},{\xi}}(\phi)=\sum_{j=0,\cdots,4}\langle\phi,P\psi_{\delta,\xi}^{j}\rangle
P\psi_{\delta,\xi}^{j}\qquad\hbox{and}\qquad\Pi_{{\delta},{\xi}}^{\perp}(\phi)=\phi-\Pi_{{\delta},{\xi}}(\phi).$

Equation (2.1) can be rewritten as the following system of two equations

$\Pi_{{\delta},{\xi}}^{\perp}\big{[}\mathcal{L}_{\delta,\xi}(\phi)-\mathcal{E}_{\delta,\xi}-\mathcal{N}_{\delta,\xi}(\phi)\big{]}=0$
(2.3)

and

$\Pi_{{\delta},{\xi}}\big{[}\mathcal{L}_{\delta,\xi}(\phi)-\mathcal{E}_{\delta,\xi}-\mathcal{N}_{\delta,\xi}(\phi)\big{]}=0,$
(2.4)

where the linear operator $\mathcal{L}_{\delta,\xi}$ is

$\mathcal{L}_{\delta,\xi}(\phi)=\phi-{\mathtt{i}}^{*}\left(3\phi
W_{{\delta},{\xi}}^{2}+\varepsilon\phi\right),$

the error term $\mathcal{E}_{\delta,\xi}$ is

$\mathcal{E}_{\delta,\xi}={\mathtt{i}}^{*}\left(W_{{\delta},{\xi}}^{3}+\varepsilon
W_{{\delta},{\xi}}\right)-W_{{\delta},{\xi}}$

and the nonlinear term $\mathcal{N}_{\delta,\xi}$ is

$\mathcal{N}_{\delta,\xi}(\phi)={\mathtt{i}}^{*}\left(\phi^{3}+3\phi^{2}W_{{\delta},{\xi}}\right).$

##  3\. Solving equation (2.3)

Given $\rho>0$ small, let $\mathcal{O}_{\rho}=\left\\{{\xi}\in\Omega\ |\
\operatorname{dist}(\xi,\partial\Omega)\geq\rho\right\\}.$

First of all, we estimate the error term $\mathcal{E}_{\delta,\xi}.$

######  Lemma 3.1.

For any $\rho>0$ small enough there exist $\varepsilon_{0}>0$ and $c>0$ such
that for any $\xi\in\mathcal{O}_{\rho}$ and for any
$\varepsilon,\delta\in(0,\varepsilon_{0})$ it holds

$\|\mathcal{E}_{\delta,\xi}\|\lesssim|{\delta}|^{2}+\varepsilon|{\delta}|.$

###### Proof.

First of all, we remark that
$W_{{\delta},{\xi}}={\mathtt{i}}^{*}\left(U_{\delta,\xi}^{3}\right).$ Then by
a straightforward computation

$\displaystyle\|\mathcal{E}_{\delta,\xi}\|$
$\displaystyle\lesssim\|(PU_{\delta,\xi})^{3}-U_{\delta,\xi}^{3}\|_{\frac{4}{3}}+\varepsilon\|V\|_{\infty}\|PU_{\delta,\xi}\|_{\frac{4}{3}}$

where

$\displaystyle\|PU_{\delta,\xi}^{3}-U_{\delta,\xi}^{3}\|$
$\displaystyle\lesssim\delta\|U_{\delta,\xi}^{2}\|_{\frac{4}{3}}+\delta^{2}\|U_{\delta,\xi}\|_{\frac{4}{3}}+\delta^{3}\lesssim\delta^{2}.$

and $\|PU_{\delta,\xi}\|_{\frac{4}{3}}=\mathcal{O}(\delta)$. ∎

Next, we state the invertibility of the linear operator
$\mathcal{L}_{\delta,\xi}.$

######  Lemma 3.2.

For any $\rho>0$ small enough there exist $\varepsilon_{0}>0$ and $c>0$ such
that for any $\varepsilon,\delta\in(0,\varepsilon_{0}),$ $i=1,\cdots,k$, and
for any ${\xi}\in\mathcal{O}_{\rho}$

$\|\left(\Pi^{\perp}_{\delta,\xi}\circ\mathcal{L}_{\delta,\xi}\right)(\phi)\|\geq
c\|\phi\|\mbox{ for all }\phi\in K_{\delta,\xi}^{\perp}.$

Furthermore, the operator
$\Pi^{\perp}_{\delta,\xi}\circ\mathcal{L}_{\delta,\xi}$ is invertible and its
inverse is continuous.

###### Proof.

We omit the proof because it is enough to apply the arguments used in the
proof of Lemma 1.7 in [15] to the $4-$dimensional case. ∎

Finally, using a classical fixed point argument, we can solve equation (2.3).

######  Proposition 3.3.

For any $\rho>0$ small enough there exists $\varepsilon_{0}>0$ such that for
any $\varepsilon,\delta\in(0,\varepsilon_{0})$ and $\xi\in\mathcal{O}_{\rho}$,
there exists a unique $\phi=\phi_{\delta,\xi}\in K_{\delta,\xi}^{\perp}$
solving (2.3), which also satisfies

$\|\phi\|\lesssim{\delta}^{2}+\varepsilon{\delta}.$ (3.1)

##  4\. Solving equation (2.4)

Let $u=U_{\delta,\xi}+\phi$ (see (2.2)). By (2.3) we deduce that

$-\Delta u-u^{3}-\epsilon
Vu=\sum\limits_{j=0,\dots,4}c^{j}U_{\delta,\xi}^{2}\psi_{\delta,\xi}^{j}$
(4.1)

for some real numbers $c^{j}=c^{j}(\delta,\xi)$. Then solving equation (2.4)
is equivalent to find $(\delta,\xi)$ such that the $c^{j}$’s are zero.

First of all, we prove a sufficient condition which ensures that all the
$c^{j}$’s in (4.1) are zero. Set $\partial_{j}u:={\partial u\over\partial
x_{j}}.$

######  Lemma 4.1.

If

$\int\limits_{\Omega}(-\Delta u-u^{3}-\epsilon Vu)\psi_{\delta,\xi}^{0}=0$
(4.2)

and for some $\eta>0$

$\int\limits_{B(\xi,\eta)}(-\Delta u-u^{3}-\epsilon
Vu)\partial_{j}u=0,\quad\mbox{ for all }j=1,2,3,4$ (4.3)

then $c^{\ell}=0$ for any $\ell=0,1,\dots,4.$

###### Proof.

By (4.2) and (4.3)

$\sum\limits_{\ell=0}^{4}c^{\ell}\int\limits_{\Omega}U_{\delta,\xi}^{2}\psi_{\delta,\xi}^{\ell}\psi_{\delta,\xi}^{0}=\sum\limits_{\ell=0}^{4}c^{\ell}\int\limits_{B(\xi,\eta)}U_{\delta,\xi}^{2}\psi_{\delta,\xi}^{\ell}\partial_{j}u=0\quad\mbox{
for all }j=0,\dots,4$

and the linear system in the $c^{\ell}$’s is diagonally dominant. Indeed

$\int\limits_{\Omega}U_{\delta,\xi}^{2}\psi_{\delta,\xi}^{\ell}\psi_{\delta,\xi}^{0}=\begin{cases}a+\mathcal{O}(\delta^{4})&\mbox{if
}\ell=0\\\ \mathcal{O}(\delta^{5})&\mbox{if }\ell=1,\cdots,4\\\ \end{cases}$

for some constant $a\not=0.$ Moreover,

$\int\limits_{B(\xi,\eta)}U_{\delta,\xi}^{2}\psi_{\delta,\xi}^{\ell}\partial_{j}u=\begin{cases}\frac{1}{\delta}b+\mathcal{O}(\varepsilon)&\mbox{if
}\ell=j\\\ \mathcal{O}(\varepsilon)&\mbox{if }\ell\neq j\\\ \end{cases}$

for some constant $b\not=0,$ because by (2.2) and (3.1),

$\int\limits_{B(\xi_{i},\eta_{i})}U_{\delta,\xi}^{2}\psi_{\delta,\xi}^{\ell}\partial_{j}\phi=\mathcal{O}(\|\phi\|)\left(\int_{B(\xi,\eta)}|U_{\delta,\xi}\psi_{\delta,\xi}^{\ell}|^{2}\right)^{1/2}=\mathcal{O}(\varepsilon)\quad\mbox{for
all }\ell=0,\cdots,4$

and

$\int\limits_{B(\xi,\eta)}U_{\delta,\xi}^{2}\psi_{\delta,\xi}^{\ell}\partial_{j}PU_{\delta,\xi}=\begin{cases}\frac{1}{\delta}b+\mathcal{O}(\delta^{2})&\mbox{if
}\ell=j\\\ \mathcal{O}(\delta^{2})&\mbox{if }\ell\neq j\\\ \end{cases}$

∎

Next, we write the first order term of (4.2).

######  Lemma 4.2.

For any $\rho>0$ small enough

$\int\limits_{\Omega}(-\Delta u-u^{3}-\epsilon
V(x)u)\psi_{\delta,\xi}^{0}=\mathfrak{C}^{2}\delta^{2}\tau(\xi)+\varepsilon\delta^{2}\ln\delta\mathfrak{c}^{2}\omega
V(\xi)+o(\delta^{2}),$ (4.4)

as $\varepsilon,\delta\to 0$, uniformly respect to $\xi\in\mathcal{O}_{\rho}$.

###### Proof.

We point out that

$\displaystyle\int\limits_{\Omega}(-\Delta u-u^{3}-\epsilon
Vu)\psi_{\delta,\xi}^{0}$
$\displaystyle=\underbrace{\int\limits_{\Omega}\left(-\Delta
U_{\delta,\xi}-U_{\delta,\xi}^{3}-\epsilon
U_{\delta,\xi}\right)\psi_{\delta,\xi}^{0}}_{=:I_{1}}$
$\displaystyle+\underbrace{\int\limits_{\Omega}\left(-\Delta\phi-3U_{\delta,\xi}^{2}\phi\right)\psi_{\delta,\xi}^{0}}_{=:I_{2}}$
$\displaystyle+\underbrace{\int\limits_{\Omega}\left(\phi^{3}+3U_{\delta,\xi}\phi^{2}\right)\psi_{\delta,\xi}^{0}}_{=:I_{3}}\
$

First of all, let us prove that

$I_{1}=\delta^{2}\mathfrak{C}^{2}H(\xi,\xi)+\varepsilon\delta^{2}\ln\delta\mathfrak{c}^{2}\omega
V(\xi)+o(\delta^{2}).$ (4.5)

We observe that

$\displaystyle I_{1}$
$\displaystyle=\int_{\Omega}\left[U_{\delta,\xi}^{3}-PU_{\delta,\xi}^{3}-\varepsilon
V(x)PU_{\delta,\xi}\right]\psi_{\delta,\xi}^{0}$
$\displaystyle=\int_{\Omega}(U_{\delta,\xi}^{3}-PU_{\delta,\xi}^{3})\psi_{\delta,\xi}^{0}-\varepsilon\int_{\Omega}V(x)PU_{\delta,\xi}\psi_{\delta,\xi}^{0}.$

Now

$\displaystyle\int_{\Omega}(U_{\delta,\xi}^{3}-PU_{\delta,\xi}^{3})$
$\displaystyle\psi_{\delta,\xi}^{0}=3\int_{B(\xi,\rho)}U_{\delta,\xi}^{2}\underbrace{(U_{\delta,\xi}-PU_{\delta,\xi})}_{\delta\mathfrak{C}H(x,\xi)+\mathcal{O}(\delta^{3})}\psi_{\delta,\xi}^{0}+\mathcal{O}(\delta^{3})$
$\displaystyle=3\delta\mathfrak{C}\int_{\Omega}H(x,\xi)U_{\delta,\xi}^{2}\psi_{\delta,\xi}^{0}+\mathcal{O}(\delta^{3})$
$\displaystyle=3\delta^{4}\mathfrak{C}\mathfrak{c}^{3}\int_{B(\xi,\rho)}H(x,\xi)\frac{|x-\xi|^{2}-\delta^{2}}{\left(\delta^{2}+|x-\xi|^{2}\right)^{4}}+o(\delta^{2})$
$\displaystyle=3\delta^{2}\mathfrak{C}\mathfrak{c}^{3}\int_{B(0,\rho/\delta)}\underbrace{H(\xi+\delta
t,\xi)}_{H(\xi,\xi)+\mathcal{O}(\delta)}\frac{|t|^{2}-1}{(1+|t|^{2})^{4}}+o(\delta^{2})$
$\displaystyle=\delta^{2}\mathfrak{C}^{2}H(\xi,\xi)+o(\delta^{2})$

because

$\int_{\mathbb{R}^{4}}\frac{|t|^{2}-1}{(1+|t|^{2})^{4}}=\frac{\omega}{12},$

and

$\displaystyle\varepsilon\int_{\Omega}V(x)PU_{\delta,\xi}\psi_{\delta,\xi}^{0}$
$\displaystyle=\varepsilon\delta^{2}\mathfrak{c}^{2}\int_{B(\xi,\rho)}V(x)\frac{|x-\xi|^{2}-\delta^{2}}{(\delta^{2}+|x-\xi|^{2})^{3}}+\mathcal{O}(\delta^{3})$
$\displaystyle=\varepsilon\delta^{2}\mathfrak{c}^{2}\int_{B(0,\rho/\delta)}V(\delta
t+\xi)\frac{|y|^{2}-1}{(1+|y|^{2})^{3}}+\mathcal{O}(\delta^{3})$
$\displaystyle=-\varepsilon\delta^{2}\ln\delta\mathfrak{c}^{2}\omega
V(\xi)+o(\delta^{2}\varepsilon\ln\delta).$

Let us estimate the other terms. It is important to point out the estimate

$\int_{\partial\Omega}|\partial_{\nu}\phi|^{2}=o(\delta^{2})$

proved in [19]. Then, recalling that
$-\Delta\psi_{\delta,\xi}^{0}=3U_{i}^{2}\psi_{\delta,\xi}^{0}$, for all
$i=1\cdots,k$, we have

$\displaystyle\int_{\Omega}(-\Delta\phi)\psi_{\delta,\xi}^{0}$
$\displaystyle=\underbrace{\int_{\Omega}\phi(-\Delta\psi_{\delta,\xi}^{0})}_{=3\int_{\Omega}\phi
U_{i}^{2}\psi_{\delta,\xi}^{0}}+\int_{\partial\Omega}\underbrace{\phi}_{=0}\nabla\psi_{\delta,\xi}^{0}\cdot\nu-\int_{\partial\Omega}\underbrace{\psi_{\delta,\xi}^{0}}_{=\mathcal{O}(\delta)}\nabla\phi\cdot\nu$
$\displaystyle=3\int_{\Omega}\phi
U_{i}^{2}\psi_{\delta,\xi}^{0}+o(|\delta|^{2}).$

Therefore,

$\displaystyle|I_{2}|$ $\displaystyle\lesssim
3\int_{\Omega}\phi|PU_{\delta,\xi}^{2}-U_{\delta,\xi}^{2}||\psi_{\delta,\xi}^{0}|+\varepsilon\int_{\Omega}|\phi||\psi_{\delta,\xi}^{0}|+o(|\delta|^{2})$
(4.6)
$\displaystyle\lesssim\|\phi\|_{4}\|(PU_{\delta,\xi}^{2}-U_{i}^{2})\psi_{\delta,\xi}^{0}\|_{4/3}+\varepsilon\|\phi\|_{4}\|\psi_{\delta,\xi}^{0}\|_{4/3}=o\left(\delta^{2}\right).$

and

$\displaystyle|I_{3}|$
$\displaystyle\leqslant\int_{\Omega}\left[|\phi|^{3}+3|\phi|^{2}|PU_{\delta,\xi}|\right]|\psi_{\delta,\xi}^{0}|$
(4.7)
$\displaystyle\lesssim\|\phi\|^{3}\|\psi_{\delta,\xi}^{0}\|_{4}+\|\phi\|^{2}\|PU_{\delta,\xi}\psi_{\delta,\xi}^{0}\|_{2}=o(\delta^{2}).$

Finally, (4.4) follows by (4.5), (4.6) and (4.7). ∎

Finally, we write the first order term of (4.3).

######  Lemma 4.3.

For any $\rho>0,$ there exists $\eta_{i}>0$ such that

$\int\limits_{B(\xi_{i},\eta_{i})}(-\Delta u-u^{3}-\epsilon
Vu)\partial_{j}u=-\frac{1}{2}\delta^{2}\left[\mathfrak{C}^{2}{\partial_{j}}\tau(\xi)+(\varepsilon\ln\delta)\omega\mathfrak{c}^{2}{\partial_{j}}V(\xi)\right]+o\left(\delta^{2}\right),$

as $\varepsilon,\delta\to 0$, uniformly respect to $\xi\in\mathcal{O}_{\rho}$.

###### Proof.

First of all, we point out that

$\int\limits_{B(\xi_{i},\eta_{i})}(-\Delta
u-u^{3})\partial_{j}u=\int\limits_{\partial
B(\xi_{i},\eta_{i})}\left(-\partial_{\nu}u\partial_{j}u+\frac{1}{2}|\nabla
u|^{2}\nu_{j}-\frac{1}{4}u^{4}\nu_{j}\right)$ (4.8)

By Lemma 4.4 and (3.1), we choose $\eta_{i}$ such that

$\int_{\partial
B(\xi_{i},\eta_{i})}\left(|\nabla\phi|^{2}+|\phi|^{4}\right)\lesssim\delta^{4}.$
(4.9)

Now, by (2.2)

$U_{\delta,\xi}(x):={\delta\mathfrak{c}\over|x-\xi|^{2}}-\mathfrak{C}\delta
H(x,\xi)+\mathcal{O}\left(|\delta|^{2}\right)$ (4.10)

$C^{1}-$uniformly on $\partial B(\xi,\eta).$ It is crucial to point out that
the function $H(x,\xi)$ is harmonic. Therefore, by (4.8), (4.9) and (4.10)

$\displaystyle\int\limits_{B(\xi,\eta)}(-\Delta u-u^{3})\partial_{j}u$
$\displaystyle=$ $\displaystyle\delta^{2}\int\limits_{\partial
B(\xi,\eta)}\left[-\partial_{\nu}\left({\mathfrak{c}\over|x-\xi|^{2}}+\mathfrak{C}H(x,\xi)\right)\partial_{j}\left({\mathfrak{c}\over|x-\xi|^{2}}+\mathfrak{C}H(x,\xi)\right)\right.$
$\displaystyle+\delta^{2}\frac{1}{2}\left.\left|\nabla\left({\mathfrak{c}\over|x-\xi_{i}|^{2}}+\mathfrak{C}H(x,\xi)\right)\right|^{2}\nu_{j}\right]+o\left(|\delta|^{2}\right)$
$\displaystyle=$ $\displaystyle-\delta^{2}\int_{\partial
B(\xi,\eta)}\nabla\left(\mathfrak{c}\over|x-\xi|^{2}\right)\cdot\nu\partial_{x_{j}}H(x,\xi)+o(\delta^{2})$
$\displaystyle=$ $\displaystyle
2\mathfrak{c}\mathfrak{C}\delta^{2}\frac{1}{|\eta|^{3}}\int\limits_{\partial
B(\xi,\eta)}\partial_{j}H(x,\xi)+o\left(\delta^{3}\right)={\mathfrak{C}^{2}}\delta^{2}\partial_{j}H(x,\xi)\big{|}_{x=\xi}+o\left(\delta^{2}\right)$

because $H(x,\xi)$ is harmonic and also $2\mathfrak{c}\omega=\mathfrak{C}$ .
Finally, as $\tau(x)=H(x,x)$ denotes the Robin’s function, by

$\partial_{\xi_{j}}\tau(\xi)=(\partial_{x_{j}}H(x,y)+\partial_{y_{j}}H(x,y))|_{(x,y)=(\xi,\xi)}=2\partial_{x_{j}}H(x,y)|_{(x,y)=(\xi,\xi)}$

follows that

$\int_{B(\xi,\eta)}(-\Delta
u-u^{3})\partial_{j}u=\frac{1}{2}\mathfrak{C}^{2}\delta^{2}\partial_{j}\tau(\xi)+o(\delta^{2}).$
(4.11)

Arguing in a similar way, we also have

$\displaystyle-\varepsilon\int_{B(\xi,\eta)}V(x)u\partial_{j}u$
$\displaystyle=\varepsilon\frac{1}{2}\int_{B(\xi,\eta)}\left(\partial_{j}V(x)\right)u^{2}-\varepsilon\frac{1}{2}\int_{\partial
B(\xi,\eta)}V(x)u^{2}\nu_{j}$ (4.12)
$\displaystyle=-\frac{1}{2}\varepsilon\delta^{2}\omega\mathfrak{c}^{2}\partial_{j}V(\xi)\ln\delta+o(\delta^{2})$

By (4.11) and (4.12) the claim follows. ∎

######  Lemma 4.4.

If there exists $C_{1}>0$ such that

$\int_{B(\xi,\eta_{1})\backslash B(\xi,\eta_{2})}|f(x)|dx\leqslant C_{1},$
(4.13)

then there exist $C_{2}>0$ and $\bar{\eta}\in(\eta_{1},\eta_{2})$ such that

$\int_{\partial B(\xi,\bar{\eta})}|f(x)|dx\leqslant C_{2}.$

###### Proof.

Assume that for any $\eta\in(\eta_{1},\eta_{2})$ and $C>0$

$\int_{\partial B(\xi,\eta)}|f(x)|dx>C.$

Then by the coarea formula

$\int_{B(\xi,\eta_{1})\backslash
B(\xi,\eta_{2})}|f(x)|dx=\int_{\eta_{1}}^{\eta_{2}}\left(\int_{\partial
B(\xi,\eta)}|f(x)|dx\right)d\eta>C(\eta_{2}-\eta_{1}),$

and a contradiction arises. ∎

###  4.1. Proof of Theorem 1.1: completed

If we set $\delta=e^{-\frac{t}{\varepsilon}}$, with $t>0$, by Lemma 4.1, Lemma
4.2 and Lemma 4.3, the problem reduces to find $t>0$ and $\xi\in\Omega$ such
that

$\begin{cases}c\tau(\xi)-tV(\xi)+o(1)=0\ ,\\\ c\nabla\tau(\xi)-t\nabla
V(\xi)+o(1)=0.\end{cases}$ (4.14)

with $c:=\frac{\mathfrak{C}^{2}}{\omega\mathfrak{c}^{2}}.$ Let
$\xi_{0}\in\Omega$ be a non-degenerate critical point of the function
$f(\xi):={\tau(\xi)\over V(\xi)}$ with $V(\xi_{0})>0$. Then the point
$(t_{0},\xi_{0})$, $t_{0}:=c\frac{\tau(\xi_{0})}{V(\xi_{0})}$, is an isolated
zero of the function
$F(t,\xi):(0,+\infty)\times\Omega\to\mathbb{R}\times\mathbb{R}^{4}$ defined by

$F(t,\xi):=\Big{(}c\tau(\xi)-tV(\xi),c\nabla\tau(\xi)-t\nabla V(\xi)\Big{)}.$

We claim that the local degree

$\mathtt{degloc}\Big{(}F,(t_{0},\xi_{0})\Big{)}\not=0.$ (4.15)

Indeed,

$F^{\prime}(t_{0},x_{0})=\left(\begin{matrix}&-t_{0}V(\xi_{0})&-t_{0}\nabla
V(\xi_{0})\\\ &c\nabla\tau(\xi_{0})-t_{0}\nabla
V(\xi_{0})&c\mathcal{D}^{2}\tau(\xi_{0})-t_{0}\mathcal{D}^{2}V(\xi_{0})\\\
\end{matrix}\right)$

where

$c\nabla\tau(\xi_{0})-t_{0}\nabla V(\xi_{0})=-c\tau(\xi_{0})\nabla
f(\xi_{0})=0,$

so

$\displaystyle\mathtt{det}\ F^{\prime}(t_{0},x_{0})$
$\displaystyle=c^{2}\frac{\tau(\xi_{0})}{V(\xi_{0})}\left(-V(\xi_{0})\mathcal{D}^{2}\tau(\xi_{0})+\tau(\xi_{0})\mathcal{D}^{2}V(\xi_{0})\right)$
$\displaystyle=-c^{2}{V^{3}(\xi_{0})\tau(\xi_{0})}\mathtt{det}\
\mathcal{D}^{2}f(\xi_{0})\not=0.$

Finally, (4.15) implies that if $\epsilon$ is small enough the system (4.14)
has a solution $t_{\epsilon},\xi_{\epsilon}$ such that $t_{\epsilon}\to t_{0}$
and $\xi_{\varepsilon}\to\xi_{0}$ as $\epsilon\to 0.$ That concludes the
proof.

## References

  * [1] A. L. Amadori, F. Gladiali, M. Grossi, A. Pistoia, and G. Vaira.  A complete scenario on nodal radial solutions to the Brezis Nirenberg problem in low dimensions.  Nonlinearity, 34(11):8055–8093, 2021.
  * [2] T. Aubin.  Problèmes isopérimétriques et espaces de Sobolev.  J. Differential Geometry, 11(4):573–598, 1976.
  * [3] G. Bianchi and H. Egnell.  A note on the Sobolev inequality.  J. Funct. Anal., 100(1):18–24, 1991.
  * [4] H. Brézis and L. Nirenberg.  Positive solutions of nonlinear elliptic equations involving critical Sobolev exponents.  Comm. Pure Appl. Math., 36(4):437–477, 1983.
  * [5] H. Brezis and L. A. Peletier.  Asymptotics for elliptic equations involving critical growth.  In Partial differential equations and the calculus of variations, Vol. I, volume 1 of Progr. Nonlinear Differential Equations Appl., pages 149–192. Birkhäuser Boston, Boston, MA, 1989.
  * [6] M. del Pino, J. Dolbeault, and M. Musso.  The Brezis-Nirenberg problem near criticality in dimension 3.  J. Math. Pures Appl. (9), 83(12):1405–1456, 2004.
  * [7] O. Druet.  Elliptic equations with critical Sobolev exponents in dimension 3.  Ann. Inst. H. Poincaré C Anal. Non Linéaire, 19(2):125–142, 2002.
  * [8] P. Esposito, N. Ghoussoub, A. Pistoia, and G. Vaira.  Sign-changing solutions for critical equations with Hardy potential.  Anal. PDE, 14(2):533–566, 2021.
  * [9] Z.-C. Han.  Asymptotic approach to singular solutions for nonlinear elliptic equations involving critical Sobolev exponent.  Ann. Inst. H. Poincaré C Anal. Non Linéaire, 8(2):159–174, 1991.
  * [10] T. König and P. Laurain.  Fine multibubble analysis in the higher-dimensional brezis-nirenberg problem, 2022.
  * [11] T. König and P. Laurain.  Multibubble blow-up analysis for the brezis-nirenberg problem in three dimensions, 2023.
  * [12] A. M. Micheletti and A. Pistoia.  On the effect of the domain geometry on the existence of sign changing solutions to elliptic problems with critical and supercritical growth.  Nonlinearity, 17(3):851–866, 2004.
  * [13] R. Molle and A. Pistoia.  Concentration phenomena in elliptic problems with critical and supercritical growth.  Adv. Differential Equations, 8(5):547–570, 2003.
  * [14] F. Morabito, A. Pistoia, and G. Vaira.  Towering phenomena for the Yamabe equation on symmetric manifolds.  Potential Anal., 47(1):53–102, 2017.
  * [15] M. Musso and A. Pistoia.  Multispike solutions for a nonlinear elliptic problem involving the critical Sobolev exponent.  Indiana Univ. Math. J., 51(3):541–579, 2002.
  * [16] M. Musso and D. Salazar.  Multispike solutions for the Brezis-Nirenberg problem in dimension three.  J. Differential Equations, 264(11):6663–6709, 2018.
  * [17] A. Pistoia and G. Vaira.  Clustering phenomena for linear perturbation of the Yamabe equation.  In Partial differential equations arising from physics and geometry, volume 450 of London Math. Soc. Lecture Note Ser., pages 311–331. Cambridge Univ. Press, Cambridge, 2019.
  * [18] B. Premoselli.  Towers of bubbles for Yamabe-type equations and for the Brézis-Nirenberg problem in dimensions $n\geq 7$.  J. Geom. Anal., 32(3):Paper No. 73, 65, 2022.
  * [19] O. Rey.  The role of the Green’s function in a nonlinear elliptic equation involving the critical Sobolev exponent.  J. Funct. Anal., 89(1):1–52, 1990.
  * [20] G. Talenti.  Best constant in Sobolev inequality.  Ann. Mat. Pura Appl. (4), 110:353–372, 1976.