1 Introduction and main results

Consider the following semilinear elliptic equations involving Grushin operators:

$$\begin{aligned} G_\alpha u+(\alpha +1)^2|y|^{2\alpha }K(x)u=(\alpha +1)^2u^{\frac{Q_\alpha +2}{Q_\alpha -2}},\quad u>0,\ x=(y,z)\in \mathbb {R}^{m_1}\times \mathbb {R}^{m_2}, \end{aligned}$$
(1.1)

where \(\alpha \ge 0\), \(\{m_1,m_2\}\subset \mathbb {N}^+\), K(x) is a function defined in \(\mathbb {R}^{m_1+m_2}\),

$$\begin{aligned} G_\alpha :=-\Delta _y-(\alpha +1)^2|y|^{2\alpha }\Delta _z \end{aligned}$$

is the so-called Grushin operator, \(Q_\alpha :=m_1+(\alpha +1)m_2\) is the appropriate homogeneous dimension, and the power \(\frac{Q_\alpha +2}{Q_\alpha -2}\) is the corresponding critical exponent. For general case \(\alpha >0\), the entire positive solutions of (1.1) with \(K(x)=0\) were discussed by Monti and Morbidelli [21]. For \(\alpha =1\), the problem (1.1) becomes into

$$\begin{aligned} (-\Delta _y-4|y|^2\Delta _z)u+4|y|^2K(x)u=4u^{\frac{m_1+2m_2+2}{m_1+2m_2-2}},\quad u>0,\ x=(y,z)\in \mathbb {R}^{m_1}\times \mathbb {R}^{m_2}. \end{aligned}$$
(1.2)

This case appeared very early in connection with the Cauchy–Riemann Yamabe problem solved by Jerison and Lee [14]. The Sharp Sobolev estimates for the Grushin operator \(G_1=-\Delta _y-4|y|^2\Delta _z\) were calculated by Beckner [1] in low dimension by using hyperbolic symmetry and conformal geometry.

However, as far as we know, there is very little literature on the existence of infinitely many solutions for (1.1) and even (1.2), which is one of our motivations to study this type of problem. In this paper, we shall study the critical Grushin-type problem (1.2) by using Lyapunov–Schmidt reduction argument and local Pohozaev identities. Under appropriate assumptions on K(x), we will prove that (1.2) possesses infinitely many multi-bubbling solutions with arbitrarily large energy and cylindrical symmetry. Our first technique is to transform (1.2) into a Hardy–Sobolev-type problem by a tricky change of variable.

If \(K(x)=K(|y|,z)\), and \(u=\psi (|y|,z)\) solves (1.2), then for \(|y|=\gamma \),

$$\begin{aligned} -\psi _{\gamma \gamma }(\gamma ,z)-\frac{m_1-1}{\gamma }\psi _\gamma (\gamma ,z)-4\gamma ^2\Delta _z\psi (\gamma ,z)+4\gamma ^2K(\gamma ,z)\psi (\gamma ,z)=4\psi ^{\frac{m_1+2m_2+2}{m_1+2m_2-2}}(\gamma ,z). \end{aligned}$$

We define \(v(\gamma ,z)=\psi (\sqrt{\gamma },z),\) then

$$\begin{aligned} \psi _\gamma (\sqrt{\gamma },z)=2\sqrt{\gamma } v_\gamma (\gamma ,z),\quad \psi _{\gamma \gamma }(\sqrt{\gamma },z)=4\gamma v_{\gamma \gamma }(\gamma ,z)+2v_\gamma (\gamma ,z). \end{aligned}$$

Hence v satisfies

$$\begin{aligned} -v_{\gamma \gamma }(\gamma ,z)-\frac{m_1}{2\gamma }v_\gamma (\gamma ,z)-\Delta _z v(\gamma ,z)+K(\sqrt{\gamma },z)v(\gamma ,z)=\frac{1}{\gamma }v^{\frac{m_1+2m_2+2}{m_1+2m_2-2}}(\gamma ,z). \end{aligned}$$

If we denote \(V(x)=K(\sqrt{|y|},z)\), \(k=\frac{m_1+2}{2}\) when \(m_1\) is even, \(N=k+m_2\) and \(2^\star :=\frac{2(N-1)}{N-2}\), then \(u=v(|y|,z)\) solves

$$\begin{aligned} -\Delta u(x)+V(x)u(x)=\frac{u^{2^\star -1}(x)}{|y|},\quad u>0,\ x=(y,z)\in \mathbb {R}^k\times \mathbb {R}^{N-k}. \end{aligned}$$
(1.3)

That is to say, the problem (1.3) is a special form of (1.2).

In the main part of this paper, we will construct infinitely many cylindrically symmetric multi-bubbling solutions for problem (1.3) by applying Lyapunov–Schmidt reduction argument and local Pohozaev identities. As derived and analyzed above, these solutions of (1.3) will produce infinitely many multi-bubbling solutions with cylindrical symmetry for the Grushin-type problem (1.2) when V(x) is cylindrically symmetric. This is also a motivation for us to study the problem (1.3).

Our main idea is motivated by the paper of Peng et al. [22], where the case of \(\alpha =0\) for problem (1.1) is concerned, which leads to study the following problem

$$\begin{aligned} -\Delta u+K(x)u=u^{\frac{N+2}{N-2}},\quad u>0\quad \text {in }\mathbb {R}^N. \end{aligned}$$
(1.4)

The first existence result for (1.4) as far as our knowledge is due to Benci and Cerami [3], where they assumed that K was nonnegative and \(\Vert K\Vert _{L^{N/2}(\mathbb {R}^N)}\) was suitably small. Later in Chen et al. [8], it was proved that (1.4) had infinitely many nonradial solutions under the assumptions that nonnegative K was radially symmetric and \(r^2K(r)\) had a local maximum point or a local minimum point. Recently, Vétois and Wang [25] extended the result in [8] to the optimal dimension four.

In recent work by Peng et al. (see [22]), the authors relaxed the conditions on K and they assumed that \(K(x)=K(|x'|,x'')\) for \((x',x'')\in \mathbb {R}^2\times \mathbb {R}^{N-2}\) and \(r^2K(r,x'')\) has a stable critical point \((r_0,x_0'')\) with \(r_0>0,\ K(r_0,x_0'')>0\) and \(\deg \big (\nabla (r^2K(r,x'')),(r_0,x_0'')\big )\ne 0\), under which they constructed infinitely many positive bubbling solutions for (1.4) that all concentrate at \((r_0,x_0'')\). It is worthwhile to point out that in that paper for the first time, the authors used local Pohozaev identities in determining the location of the bubbles. Motivated by their idea, in the present paper, we will also construct infinitely many bubbling solutions which all concentrate at one point which can be a saddle point of some function related to the potential in our problem. We shall present the details at the end of this section.

For other results about the existence of multiple solutions to noncompact elliptic problems, please refer to [7, 12, 16, 18, 19, 24, 26,27,28,29,30,31,32].

Now we state our assumptions on V(x), which are as below.

(V): \(V(x)=V(|z'|,z'')\in C^1(\mathbb {R}^N)\) is nonnegative bounded for \(x=(y,z',z'') \in \mathbb {R}^k\times \mathbb {R}^2\times \mathbb {R}^{N-k-2}\). Setting \(r=|z'|\), the function \(r^2V(r,z'')\) has a stable critical point \((r_0,z_0'')\) satisfying \(r_0>0,\ V(r_0,z_0'')>0\), and \(\deg \big (\nabla (r^2V(r,z'')),B_\varrho (r_0,z_0'')\big )\ne 0\) for some small constant \(\varrho >0\).

The potential function V(x) is attached some cylindrical symmetry. In fact, the solutions we are going to construct are also of cylindrical symmetry. For similar problems, we refer readers to [2, 6]. In [2], Badiale and Tarantello analyzed the existence and non-existence of cylindrical solutions for a nonlinear elliptic equation in \(\mathbb {R}^3\). In [6], Cao, Peng and Yan constructed some multipeak solutions for a kind of Webster scalar curvature problems on the CR sphere with cylindrically symmetric curvature.

Our main results in this paper are as follows.

Theorem 1.1

Assume that \(N\ge 5,\ \frac{N+1}{2}\le k<N-1\) and (V) holds. Then the problem (1.3) has infinitely many cylindrically symmetric multi-bubbling solutions, whose energy can be made arbitrarily large.

Corollary 1.2

Under the assumptions of Theorem 1.1, if \(m_1=2k-2,\ m_2=N-k\) and \(K(x)=V(|z'|,z'')\), then the critical Grushin-type problem (1.2) has infinitely many multi-bubbling solutions with arbitrarily large energy and cylindrical symmetry.

Remark 1.3

\((\mathrm{i})\) The condition \(\frac{N+1}{2}\le k<N-1\) is equivalent to \(1<h:=N-k\le k-1\), which is needed for the proof of Lemma B.2. \((\mathrm{ii})\) For the case of \(N=4\), the estimate is totally different and the logarithm will appear (see [25]). Thus we only consider the case of \(N\ge 5\). \((\mathrm{iii})\) In order to estimate the local Pohozaev identity (3.1), we have to constrain the potential V independent of the first layer variables y (see Lemma 3.4).

Remark 1.4

In this paper, we only consider the case of \(\alpha =1\) for the Grushin operator \(G_\alpha \) since our methods to obtain the main results require the nondegeneracy of solution to the limit equation and that is true for \(\alpha =1\) (see [4]), while for \(\alpha >0\) and \(\alpha \ne 1\), the nondegeneracy is unknown and thus the corresponding infinitely many solutions problem for general Grushin operator \(G_\alpha \) cannot be obtained here and which is valuable to consider in the future.

There is no doubt that the reduction argument has become more and more useful in studying the existence and properties of solutions for noncompact elliptic problems. The studies [5, 6, 8,9,10,11,12,13, 15,16,17,18,19, 22,23,24, 26,27,31] and so on provide a powerful support. Particularly, Wei and Yan [30, 31] created a new idea, which is to use the number of the bumps or bubbles as the parameter in constructing solutions. That is also what we do in this paper. Next we shall sketch the proof of Theorem 1.1.

It is well known from [4, 20] that

$$\begin{aligned} U_{\xi ,\lambda }(x)=\big [(N-2)(k-1)\big ]^{\frac{N-2}{2}}\left( \frac{\lambda }{(1+\lambda |y|)^2+\lambda ^2|z-\xi |^2}\right) ^{\frac{N-2}{2}},\quad \lambda >0, \xi \in \mathbb {R}^{N-k} \end{aligned}$$

are the only solutions to problem:

$$\begin{aligned} -\Delta u(x)=\frac{u^{2^\star -1}(x)}{|y|},\quad u>0,\ x=(y,z)\in \mathbb {R}^k\times \mathbb {R}^{N-k}, \end{aligned}$$
(1.5)

and \(U_{\xi ,\lambda }(x)\) is nondegenerate in \(D^{1,2}(\mathbb {R}^N)\), where

$$\begin{aligned} D^{1,2}(\mathbb {R}^N):=\left\{ \int _{\mathbb {R}^N}\frac{|u(x)|^{2^\star }}{|y|}\text {d}x<+\infty :|\nabla u|\in L^2(\mathbb {R}^N)\right\} \end{aligned}$$

endowed with the norm \(\Vert u\Vert ^2=\int _{\mathbb {R}^N}|\nabla u|^2\text {d}x\), which is induced by the inner product \((u,v)=\int _{\mathbb {R}^N}\nabla u\nabla v\text {d}x\).

We shall use \(U_{\xi ,\lambda }(x)\) to build up the approximate solution for problem (1.3). In what follows, we always assume that \(m>0\) is a large integer, \(\lambda \in \big [T_0m^{\frac{N-2}{N-4}},T_1m^{\frac{N-2}{N-4}}\big ]\) for some constants \(T_1>T_0>0\), and \(|(\bar{r},\bar{z}'')-( r_0,z_0'')|\le \varrho \) with \(\bar{z}'':=(\bar{z}_3, \bar{z}_4,\ldots ,\bar{z}_{N-k})\in \mathbb {R}^{N-k-2}\). Define

$$\begin{aligned} H_s&=\bigg \{u\in D^{1,2}(\mathbb {R}^N): u(y,z)=u(|y|,z), u(y,z_1,-z_2,z'')=u(y,z_1,z_2,z''),\\&\qquad \qquad u(y,r\cos \theta ,r\sin \theta ,z'')=u\Big (y,r\cos \Big (\theta +\frac{2\pi i}{m}\Big ),r\sin \Big (\theta +\frac{2\pi i}{m}\Big ),z''\Big )\bigg \} \end{aligned}$$

and set

$$\begin{aligned} \xi _i=\Big (\bar{r}\cos \frac{2(i-1)\pi }{m},\bar{r}\sin \frac{2(i-1)\pi }{m},\bar{z}''\Big ),\quad \ i=1,2,\ldots ,m. \end{aligned}$$

Let

$$\begin{aligned} \Omega _i=\bigg \{x=(y,z',z'')\in \mathbb {R}^k\times \mathbb {R}^2\times \mathbb {R}^{N-k-2}:\Big <\frac{z'}{|z'|},\frac{\xi _i'}{|\xi _i'|}\Big >\ge \cos \frac{\pi }{m}\bigg \}, \end{aligned}$$

where \(\xi _i':=\big (\bar{r}\cos \frac{2(i-1)\pi }{m},\bar{r}\sin \frac{2(i-1)\pi }{m}\big )\in \mathbb {R}^2,\ i=1,2,\ldots ,m\).

By the assumption (V), we can choose \(\delta >0\) as a small constant such that \(r^2V(r,z'')\ge C>0\) for \(|(r,z'')-(r_0,z_0'')|\le 10\delta \). To accelerate the decay of \(U_{\xi ,\lambda }\), we define a smooth cutoff function \(\eta (x)=\eta (|y|,|z'|,z'')\) satisfying \(\eta (x)=1\) if \(|(|y|,r,z'')-(0,r_0,z_0'')|\le \delta \) and \(\eta (x)=0\) if \(|(|y|,r,z'')-(0,r_0,z_0'')|\ge 2\delta \) with \(0\le \eta (x)\le 1\) in \(\mathbb {R}^N\). Denote

$$\begin{aligned} \bar{U}_{\xi _i,\lambda }=\eta U_{\xi _i,\lambda },\ \ \bar{W}_{\bar{r},\bar{z}'',\lambda }=\sum _{i=1}^m \bar{U}_{\xi _i,\lambda },\ \ W_{\bar{r},\bar{z}'',\lambda }=\sum _{i=1}^m U_{\xi _i,\lambda }. \end{aligned}$$

Theorem 1.1 is a direct consequence of the following result.

Theorem 1.5

Under the assumptions of Theorem 1.1, there exists an integer \(m_0>0\) such that for any integer \(m\ge m_0\), the problem (1.3) has a solution \(u_m\) of the form

$$\begin{aligned} u_m=\bar{W}_{\bar{r}_m,\bar{z}_m'',\lambda _m}+\varphi _m=\sum _{i=1}^m \eta U_{\xi _i,\lambda _m}+\varphi _m, \end{aligned}$$

where \(\varphi _m\in H_s\). Moreover, as \(m\rightarrow +\infty \), \(\lambda _m\in \big [T_0m^{\frac{N-2}{N-4}},T_1m^{\frac{N-2}{N-4}}\big ]\), \((\bar{r}_m,\bar{z}_m'')\rightarrow (r_0,z_0'')\) and \(\lambda _m^{-\frac{N-2}{2}}\Vert \varphi _m\Vert _{L^\infty (\mathbb {R}^N)}\rightarrow 0\).

We would like to remark that the concentration point \((r_0,z_0'')\) can be a saddle point of \(r^2V(r,z'')\). Thus, the procedure to determine the location of the bubbles in [8, 9, 12, 17, 27,28,29,30,31] cannot take effect any more. Because the solutions constructed there all concentrate at some local maximum points or local minimum points of some functions. Thereupon, a natural method is to estimate the derivatives of the reduced functional. However, the computations are very complicated at times (see [23] for example) and it doesn’t work at all in some cases (refer to [22]). So is our case. In fact, we only use one derivative of the reduced functional. The other two are useless for that the error terms destroy the dominant terms. Motivated and inspired by [22], we will employ some local Pohozaev identities to look for the algebraic equations determining the location of the bubbles. To be precise, we turn to prove that if \((\bar{r}, \bar{z}'')\) satisfies the local Pohozaev identities (3.1) and (3.2) (see Proposition 3.1) in a suitable neighborhood \(D_\rho \) of \((0,r_0, z_0'')\ (0\in \mathbb {R}^k)\), then \(\frac{\partial \mathcal F}{\partial \bar{r}}=0\) and \(\frac{\partial \mathcal F}{\partial \bar{z}_j}=0\ (j=3,4,\ldots ,N-k)\), where

$$\begin{aligned} \mathcal F(\bar{r},\bar{z}'',\lambda ):=I(\bar{W}_{\bar{r},\bar{z}'',\lambda }+\varphi _{\bar{r},\bar{z}'',\lambda }) \end{aligned}$$

is the reduced functional and

$$\begin{aligned} I(u):=\frac{1}{2}\int \Big (|\nabla u|^2+V(r,z'')u^2\Big )-\frac{1}{2^\star }\int \frac{|u|^{2^\star }}{|y|},\ u\in D^{1,2}(\mathbb {R}^N) \end{aligned}$$

is the energy functional corresponding to problem (1.3). By using such identities (3.1) and (3.2), we only need to estimate the error term \(\varphi _{\bar{r},\bar{z}'',\lambda }\) away from the concentration points, which simplifies the calculations greatly. Moreover, to deal with the large number of bubbles in the solution, the reduction procedure will be performed in a weighted space instead of the standard Sobolev space. Since the equations we study in this paper include a singular term, the computations become much more delicate and complicated.

The paper in the sequel is organized as follows. In Sect. 2, we will carry out the finite-dimensional reduction procedure and obtain a good estimate for the error term. Precisely, we will use \(\bar{W}_{\bar{r},\bar{z}'',\lambda }\) as the approximation solution and consider the linearization of problem (1.3) around \(\bar{W}_{\bar{r},\bar{z}'',\lambda }\). The unique solvability of the linearized problem (2.1) will be assured by contradiction discussion and Fredholm’s alternative theorem. After that, we continue to study a perturbation problem for (1.3), namely (2.21). The Contraction Mapping Principle will be employed to prove that (2.21) is uniquely solvable, in which the estimates of \(\mathcal {R}(\varphi )\) and \(\mathcal P_m\) (see (2.23)) will play an important role. In Sect. 3, we will study the reduced finite-dimensional problem to obtain a true solution and complete the proof of Theorem 1.5. In other words, we will apply the local Pohozaev identities (3.1) and (3.2) together with one derivative (3.3) of the reduced functional to find the algebraic equations (3.21), (3.22), (3.23) and finally obtain suitable parameters \(\bar{r}, \bar{z}'',\lambda \), which correspond to the true solution. All the technical estimates will be left in “Appendices A, B and C”. For simplicity and without confusion, sometimes \(\int f(x)\) will denote \(\int _{\mathbb {R}^N}f(x)\text {d}x\). Throughout this paper, C signifies various positive constants independent of m and \(\lambda \).

2 Finite-dimensional reduction

In this section, we perform a finite-dimensional reduction. Let

$$\begin{aligned} \Vert u\Vert _*=\sup _{x\in \mathbb {R}^N}\bigg (\sum _{i=1}^m\frac{1}{(1+\lambda |y|+\lambda |z-\xi _i|)^{\frac{N-2}{2}+\tau }}\bigg )^{-1}\lambda ^{-\frac{N-2}{2}}|u(x)| \end{aligned}$$

and

$$\begin{aligned} \Vert f\Vert _{**}=\sup _{x\in \mathbb {R}^N}\bigg (\sum _{i=1}^m\frac{1}{\lambda |y|(1+\lambda |y|+\lambda |z-\xi _i|)^{\frac{N}{2}+\tau }}\bigg )^{-1}\lambda ^{-\frac{N+2}{2}}|f(x)|, \end{aligned}$$

where \(\tau :=\frac{N-4}{N-2}\). Denote

$$\begin{aligned} \bar{U}_{i,1}=\frac{\partial \bar{U}_{\xi _i,\lambda }}{\partial \lambda },\ \ \bar{U}_{i,2}=\frac{\partial \bar{U}_{\xi _i,\lambda }}{\partial \bar{r}},\ \ \bar{U}_{i,j}=\frac{\partial \bar{U}_{\xi _i,\lambda }}{\partial \bar{z}_j},\quad j=3,4,\ldots ,N-k. \end{aligned}$$

For some real numbers \(c_l\), consider

$$\begin{aligned} {\left\{ \begin{array}{ll} \displaystyle -\Delta \varphi +V(r,z'')\varphi -(2^{\star }-1)\frac{\bar{W}_{\bar{r},\bar{z}'',\lambda }^{2^{\star }-2}}{|y|}\varphi =f+\sum _{l=1}^{N-k}c_l\sum _{i=1}^m\frac{\bar{U}_{\xi _i,\lambda }^{2^{\star }-2}}{|y|}\bar{U}_{i,l} \quad \text {in }\mathbb {R}^N,\\ \varphi \in H_s,\quad \displaystyle \int \frac{\bar{U}_{\xi _i,\lambda }^{2^{\star }-2}}{|y|}\bar{U}_{i,l}\varphi =0,\quad l=1,2,\ldots ,N-k,\quad i=1,2,\ldots ,m. \end{array}\right. } \end{aligned}$$
(2.1)

Lemma 2.1

Assume that \(\varphi _m\) solves problem (2.1) for \(f=f_m\). If \(\Vert f_m\Vert _{**}\rightarrow 0\) as \(m\rightarrow +\infty \), then \(\Vert \varphi _m\Vert _*\rightarrow 0\) as \(m\rightarrow +\infty \).

Proof

We argue by contradiction. Suppose that there exist \(m\rightarrow +\infty , \bar{r}_m\rightarrow r_0, \bar{z}_m''\rightarrow z_0''\), \(\lambda _m\in \big [T_0m^{\frac{N-2}{N-4}},T_1m^{\frac{N-2}{N-4}}\big ]\), and \(\varphi _m\) solves (2.1) for \(f=f_m, \lambda =\lambda _m, \bar{r}=\bar{r}_m, \bar{z}''=\bar{z}_m''\) with \(\Vert f_m\Vert _{**}\rightarrow 0\) and \(\Vert \varphi _m\Vert _*\ge C>0\). We may assume \(\Vert \varphi _m\Vert _*=1\) without loss of generality. For simplicity, we drop the subscript m.

By applying Green representation to \(|\varphi |\), we have

$$\begin{aligned} |\varphi (x)|&\le C\int \frac{1}{|\tilde{x}-x|^{N-2}}\frac{\bar{W}_{\bar{r},\bar{z}'',\lambda }^{2^{\star }-2}(\tilde{x})}{|\tilde{y}|}|\varphi (\tilde{x})|\text {d}\tilde{x}\nonumber \\&+C\int \frac{1}{|\tilde{x}-x|^{N-2}}\Big (|f(\tilde{x})|+\Big |\sum _{l=1}^{N-k}c_l\sum _{i=1}^m\frac{\bar{U}_{\xi _i,\lambda }^{2^{\star }-2}(\tilde{x})}{|\tilde{y}|}\bar{U}_{i,l}(\tilde{x})\Big |\Big )\text {d}\tilde{x}. \end{aligned}$$
(2.2)

Due to Lemmas B.1, B.2 and B.3, we can prove

$$\begin{aligned}&\int \frac{1}{|\tilde{x}-x|^{N-2}}\frac{\bar{W}_{\bar{r},\bar{z}'',\lambda }^{2^{\star }-2}(\tilde{x})}{|\tilde{y}|}|\varphi (\tilde{x})|\text {d}\tilde{x}\le C\Vert \varphi \Vert _*\sum _{i=1}^m\frac{\lambda ^{\frac{N-2}{2}}}{(1+\lambda |y|+\lambda |z-\xi _i|)^{\frac{N-2}{2}+\tau +\iota }}, \end{aligned}$$
(2.3)
$$\begin{aligned}&\int \frac{1}{|\tilde{x}-x|^{N-2}}|f(\tilde{x})|\text {d}\tilde{x}\le C\Vert f\Vert _{**}\sum _{i=1}^m\frac{\lambda ^{\frac{N-2}{2}}}{(1+\lambda |y|+\lambda |z-\xi _i|)^{\frac{N-2}{2}+\tau }}, \end{aligned}$$
(2.4)
$$\begin{aligned}&\int \frac{1}{|\tilde{x}-x|^{N-2}}\Big |\sum _{i=1}^m\frac{\bar{U}_{\xi _i,\lambda }^{2^{\star }-2}(\tilde{x})}{|\tilde{y}|}\bar{U}_{i,l}(\tilde{x})\Big |\text {d}\tilde{x}\le C\sum _{i=1}^m\frac{\lambda ^{\frac{N-2}{2}+n_l}}{(1+\lambda |y|+\lambda |z-\xi _i|)^{\frac{N-2}{2}+\tau }}, \end{aligned}$$
(2.5)

where \(n_l:=-1\) if \(l=1\) and \(n_l:=1\) if \(l=2,3,\ldots ,N-k\).

Next we estimate \(c_l\ (l=1,2,\ldots ,N-k)\). Multiplying (2.1) by \(\bar{U}_{1,l}\ (l=1,2,\ldots ,N-k)\) and integrating, we find that

$$\begin{aligned} \sum _{t=1}^{N-k}c_t\sum _{i=1}^m\int \frac{\bar{U}_{\xi _i,\lambda }^{2^{\star }-2}}{|y|}\bar{U}_{i,t}\bar{U}_{1,l}=\int \Big (-\Delta \varphi +V(r,z'')\varphi -(2^{\star }-1)\frac{\bar{W}_{\bar{r},\bar{z}'',\lambda }^{2^{\star }-2}}{|y|}\varphi -f\Big )\bar{U}_{1,l}. \end{aligned}$$
(2.6)

Taking into account the facts that

$$\begin{aligned}&\frac{1}{\lambda }\le \frac{C}{1+\lambda |y|+\lambda |z-\xi _i|},\quad \forall x\in \text {supp}\eta ,\quad i=1,2,\ldots ,m, \end{aligned}$$
(2.7)
$$\begin{aligned}&N-1-\varepsilon \ge \frac{N+2}{2}+\tau \quad \text {if } \varepsilon >0 \text { small}, \end{aligned}$$
(2.8)

and

$$\begin{aligned} \sum _{i=2}^m\frac{1}{(\lambda |\xi _i-\xi _1|)^\tau }\le C, \end{aligned}$$
(2.9)

we can deduce from Lemma B.1 that

$$\begin{aligned}&\Big |\int f(x)\bar{U}_{1,l}(x)\Big |\nonumber \\&\quad \le C\lambda ^{n_l}\Vert f\Vert _{**}\int \frac{\eta \lambda ^{\frac{N-2}{2}}}{(1+\lambda |y|+\lambda |z-\xi _1|)^{N-2}}\sum _{i=1}^m\frac{\lambda ^{\frac{N+2}{2}}}{\lambda |y|(1+\lambda |y|+\lambda |z-\xi _i|)^{\frac{N}{2}+\tau }}\nonumber \\&\quad \le C\lambda ^{n_l}\Vert f\Vert _{**}\Big (C+\sum _{i=2}^m\frac{C}{(\lambda |\xi _i-\xi _1|)^\tau }\Big )\nonumber \\&\quad \le C\lambda ^{n_l}\Vert f\Vert _{**} \end{aligned}$$
(2.10)

and

$$\begin{aligned}&\Big |\int V(r,z'')\varphi (x)\bar{U}_{1,l}(x)\Big |\nonumber \\&\quad \le C\lambda ^{n_l}\Vert \varphi \Vert _{*}\int \frac{\eta \lambda ^{\frac{N-2}{2}}}{(1+\lambda |y|+\lambda |z-\xi _1|)^{N-2}}\sum _{i=1}^m\frac{\lambda ^{\frac{N-2}{2}}}{(1+\lambda |y|+\lambda |z-\xi _i|)^{\frac{N-2}{2}+\tau }}\nonumber \\&\quad \le C\frac{\lambda ^{n_l}}{\lambda ^{1+\varepsilon }}\Vert \varphi \Vert _{*}\int \frac{\eta \lambda ^{\frac{N+2}{2}}}{(1+\lambda |y|+\lambda |z-\xi _1|)^{N-1-\varepsilon }}\sum _{i=1}^m\frac{\lambda ^{\frac{N-2}{2}}}{(1+\lambda |y|+\lambda |z-\xi _i|)^{\frac{N-2}{2}+\tau }}\nonumber \\&\quad \le C\frac{\lambda ^{n_l}}{\lambda ^{1+\varepsilon }}\Vert \varphi \Vert _{*}\int \frac{\eta \lambda ^{\frac{N+2}{2}}}{(1+\lambda |y|+\lambda |z-\xi _1|)^{\frac{N+2}{2}+\tau }}\sum _{i=1}^m\frac{\lambda ^{\frac{N-2}{2}}}{(1+\lambda |y|+\lambda |z-\xi _i|)^{\frac{N-2}{2}+\tau }}\nonumber \\&\quad \le C\lambda ^{n_l-1-\varepsilon }\Vert \varphi \Vert _{*}. \end{aligned}$$
(2.11)

Noting that

$$\begin{aligned} \sum _{i=2}^m\frac{1}{(\lambda |\xi _i-\xi _1|)^\sigma }\le C\left( \frac{m}{\lambda }\right) ^\sigma \quad \text {if }\sigma >1 \end{aligned}$$
(2.12)

and by direct computation, we can prove

$$\begin{aligned} \Big |\int \Big (-\Delta \varphi -(2^{\star }-1)\frac{\bar{W}_{\bar{r},\bar{z}'',\lambda }^{2^{\star }-2}}{|y|}\varphi \Big )\bar{U}_{1,l}\Big |\le C\lambda ^{n_l-1-\varepsilon }\Vert \varphi \Vert _{*}. \end{aligned}$$
(2.13)

Combining (2.10), (2.11) and (2.13), we have

$$\begin{aligned} \int \Big (-\Delta \varphi +V(r,z'')\varphi -(2^{\star }-1)\frac{\bar{W}_{\bar{r},\bar{z}'',\lambda }^{2^{\star }-2}}{|y|}\varphi -f\Big )\bar{U}_{1,l}=O\Big (\lambda ^{n_l}\Big (\frac{\Vert \varphi \Vert _*}{\lambda ^{1+\varepsilon }}+\Vert f\Vert _{**}\Big )\Big ). \end{aligned}$$
(2.14)

Additionally, direct calculation shows that for some \(\bar{c}>0\),

$$\begin{aligned} \sum _{i=1}^m\int \frac{\bar{U}_{\xi _i,\lambda }^{2^{\star }-2}}{|y|}\bar{U}_{i,t}\bar{U}_{1,l}=\bar{c}\delta _{t,l}\lambda ^{2n_l}+o(1)\lambda ^{n_l},\quad t=1,2,\ldots ,N-k. \end{aligned}$$
(2.15)

After (2.14) and (2.15) substituted in (2.6), we get

$$\begin{aligned} c_l=\frac{1}{\lambda ^{n_l}}\Big [o(\Vert \varphi \Vert _*)+O(\Vert f\Vert _{**})\Big ]. \end{aligned}$$
(2.16)

Now inserting (2.3), (2.4), (2.5), (2.16) into (2.2), we obtain

$$\begin{aligned} \Vert \varphi \Vert _*\le C\left( \frac{\sum _{i=1}^m\frac{1}{(1+\lambda |y|+\lambda |z-\xi _i|)^{\frac{N-2}{2}+\tau +\iota }}}{\sum _{i=1}^m\frac{1}{(1+\lambda |y|+\lambda |z-\xi _i|)^{\frac{N-2}{2}+\tau }}}+\Vert f\Vert _{**}+o(\Vert \varphi \Vert _*)\right) . \end{aligned}$$
(2.17)

Since \(\Vert \varphi \Vert _*=1\), it follows from (2.17) that there exists \(R>0\) such that

$$\begin{aligned} \Vert \lambda ^{-\frac{N-2}{2}}\varphi \Vert _{L^\infty \big (B_{\frac{R}{\lambda }}(0,\xi _1)\big )}\ge C>0,\quad \text {where } (0,\xi _1)\in \mathbb {R}^N. \end{aligned}$$
(2.18)

The function \(\tilde{\varphi }(x):=\lambda ^{-\frac{N-2}{2}}\varphi \big (\lambda ^{-1}x+(0,\xi _1)\big )\) converges uniformly in any compact set to a solution u(x) of

$$\begin{aligned} -\Delta u(x)=(2^{\star }-1)\frac{U_{0,1}^{2^{\star }-2}(x)}{|y|}u(x)\quad \text {in }\mathbb {R}^N. \end{aligned}$$
(2.19)

Owing to that \(\varphi \in H_s\) solves (2.1), then u is even in \(z_2\) and u is perpendicular to the kernel of (2.19). Thus \(u=0\), which contradicts to (2.18). \(\square \)

By Lemma 2.1 and applying similar arguments to the proof of Proposition 4.1 in [10], we can prove the following result.

Lemma 2.2

There exist \(m_0>0\) and a constant \(C>0\) independent of m such that for any \(m\ge m_0\) and all \(f\in L^\infty (\mathbb {R}^N)\), problem (2.1) has a unique solution \(\varphi :=L_m(f)\). Moreover,

$$\begin{aligned} \Vert L_m(f)\Vert _*\le C\Vert f\Vert _{**},\quad |c_l|\le \frac{C}{\lambda ^{n_l}}\Vert f\Vert _{**}. \end{aligned}$$
(2.20)

Now we consider

$$\begin{aligned} {\left\{ \begin{array}{ll} -\Delta (\bar{W}_{\bar{r},\bar{z}'',\lambda }+\varphi )+V(r,z'')(\bar{W}_{\bar{r},\bar{z}'',\lambda }+\varphi )\\ \quad =\displaystyle \frac{(\bar{W}_{\bar{r},\bar{z}'',\lambda }+\varphi )_+^{2^{\star }-1}}{|y|}+\sum _{l=1}^{N-k}c_l\sum _{i=1}^m\frac{\bar{U}_{\xi _i,\lambda }^{2^{\star }-2}}{|y|}\bar{U}_{i,l} \quad \text {in }\mathbb {R}^N,\\ \varphi \in H_s,\ \displaystyle \int \frac{\bar{U}_{\xi _i,\lambda }^{2^{\star }-2}}{|y|}\bar{U}_{i,l}\varphi =0,\quad l=1,2,\ldots ,N-k,\quad i=1,2,\ldots ,m. \end{array}\right. } \end{aligned}$$
(2.21)

The main result in this section is

Proposition 2.3

There exists \(m_0>0\) such that for each \(m\ge m_0\), \(\lambda \in \big [T_0m^{\frac{N-2}{N-4}},T_1m^{\frac{N-2}{N-4}}\big ]\), \(\bar{r}\in [r_0-\varrho ,r_0+\varrho ]\) and \(\bar{z}''\in B_\varrho (z_0'')\), problem (2.21) has a unique solution \(\varphi =\varphi _{\bar{r},\bar{z}'',\lambda }\) satisfying

$$\begin{aligned} \Vert \varphi \Vert _*\le C\Big (\frac{1}{\lambda }\Big )^{1+\varepsilon },\quad |c_l|\le C\Big (\frac{1}{\lambda }\Big )^{1+\varepsilon +n_l}, \end{aligned}$$
(2.22)

where \(\varepsilon >0\) is a small constant.

Rewrite (2.21) as

$$\begin{aligned} {\left\{ \begin{array}{ll} -\Delta \varphi +V(r,z'')\varphi -(2^{\star }-1)\frac{\bar{W}^{2^{\star }-2}_{\bar{r},\bar{z}'',\lambda }}{|y|}\varphi \\ \quad \displaystyle =\mathcal {R}(\varphi )+\mathcal P_m+\sum _{l=1}^{N-k}c_l\sum _{i=1}^m\frac{\bar{U}_{\xi _i,\lambda }^{2^{\star }-2}}{|y|}\bar{U}_{i,l} \quad \text { in }\mathbb {R}^N,\\ \varphi \in H_s,\ \displaystyle \int \frac{\bar{U}_{\xi _i,\lambda }^{2^{\star }-2}}{|y|}\bar{U}_{i,l}\varphi =0,\ l=1,2,\ldots ,N-k,\ i=1,2,\ldots ,m, \end{array}\right. } \end{aligned}$$
(2.23)

where

$$\begin{aligned} \mathcal {R}(\varphi ):=\frac{1}{|y|}\Big [\big (\bar{W}_{\bar{r}, \bar{z}'',\lambda }+\varphi \big )_+^{2^{\star }-1}-\bar{W}_{\bar{r}, \bar{z}'',\lambda }^{2^{\star }-1}-(2^{\star }-1)\bar{W}_{\bar{r}, \bar{z}'',\lambda }^{2^{\star }-2}\varphi \Big ] \end{aligned}$$

and

$$\begin{aligned} \mathcal P_m&:=\frac{1}{|y|}\Big (\bar{W}_{\bar{r}, \bar{z}'',\lambda }^{2^{\star }-1}-\sum _{i=1}^m\eta U_{\xi _i,\lambda }^{2^{\star }-1}\Big )-V(r,z'')\bar{W}_{\bar{r}, \bar{z}'',\lambda }+\Delta \eta W_{\bar{r}, \bar{z}'',\lambda }+2\nabla \eta \nabla W_{\bar{r}, \bar{z}'',\lambda }\\&:=J_0+J_1+J_2+J_3. \end{aligned}$$

In order to apply the Contraction Mapping Principle to prove that (2.23) is uniquely solvable, we need to estimate \(\mathcal {R}(\varphi )\) and \(\mathcal P_m\), respectively.

Lemma 2.4

If \(N\ge 5\), then \(\Vert \mathcal {R}(\varphi )\Vert _{**}\le C\Vert \varphi \Vert _*^{2^{\star }-1}\).

Proof

For \(\theta \in (1,2]\),

$$\begin{aligned} (1+t)_+^\theta -1-\theta t=O(|t|^\theta ),\quad \forall t\in \mathbb {R}. \end{aligned}$$

Thus \(|\mathcal {R}(\varphi )|\le C\frac{|\varphi |^{2^{\star }-1}}{|y|}\), which together with Hölder inequality and (2.9) ensures that

$$\begin{aligned} |\mathcal {R}(\varphi )|&\le C\frac{\Vert \varphi \Vert _*^{2^{\star }-1}}{|y|}\bigg (\sum _{i=1}^m\frac{\lambda ^{\frac{N-2}{2}}}{(1+\lambda |y|+\lambda |z-\xi _i|)^{\frac{N-2}{2}+\tau }}\bigg )^{2^{\star }-1}\\&\le C\frac{\Vert \varphi \Vert _*^{2^{\star }-1}}{|y|}\sum _{i=1}^m\frac{\lambda ^{\frac{N}{2}}}{(1+\lambda |y|+\lambda |z-\xi _i|)^{\frac{N}{2}+\tau }}\bigg (\sum _{i=1}^m\frac{1}{(1+\lambda |y|+\lambda |z-\xi _i|)^\tau }\bigg )^{\frac{2}{N-2}}\\&\le C\frac{\Vert \varphi \Vert _*^{2^{\star }-1}}{|y|}\sum _{i=1}^m\frac{\lambda ^{\frac{N}{2}}}{(1+\lambda |y|+\lambda |z-\xi _i|)^{\frac{N}{2}+\tau }}\bigg (1+\sum _{i=2}^m\frac{C}{(\lambda |\xi _i-\xi _1|)^\tau }\bigg )^{\frac{2}{N-2}}\\&\le C\Vert \varphi \Vert _*^{2^{\star }-1}\sum _{i=1}^m\frac{\lambda ^{\frac{N+2}{2}}}{\lambda |y|(1+\lambda |y|+\lambda |z-\xi _i|)^{\frac{N}{2}+\tau }}, \end{aligned}$$

where we use the symmetry and

$$\begin{aligned} |z-\xi _i|\ge \frac{1}{2}|\xi _i-\xi _1|,\quad \forall x=(y,z)\in \Omega _1. \end{aligned}$$
(2.24)

Therefore, \(\Vert \mathcal {R}(\varphi )\Vert _{**}\le C\Vert \varphi \Vert _*^{2^{\star }-1}\). \(\square \)

Lemma 2.5

If \(N\ge 5\), then there exists a small constant \(\varepsilon >0\) such that \(\Vert \mathcal P_m\Vert _{**}\le C\left( \frac{1}{\lambda }\right) ^{1+\varepsilon }\).

Proof

Firstly, we estimate \(J_0\). By symmetry, we may assume \(x\in \Omega _1\). Then it holds

$$\begin{aligned} |z-\xi _i|\ge |z-\xi _1|,\quad \forall x=(y,z)\in \Omega _1. \end{aligned}$$
(2.25)

If \(|(|y|,r,z'')-(0,r_0,z_0'')|\ge 2\delta \), then \(J_0=0\).

If \(\delta<|(|y|,r,z'')-(0,r_0,z_0'')|<2\delta \), then

$$\begin{aligned} C_1\lambda \le 1+\lambda |y|+\lambda |z-\xi _i|\le C_2\lambda ,\quad i=1,2,\ldots ,m. \end{aligned}$$

Thus we have

$$\begin{aligned} |J_0|\le C\Big (\frac{m}{\lambda ^{\frac{N-2}{2}}}\Big )^{2^\star -1}+Cm\Big (\frac{1}{\lambda ^{\frac{N-2}{2}}}\Big )^{2^\star -1}\le C \Big (\frac{1}{\lambda }\Big )^{1+\varepsilon } \end{aligned}$$

and

$$\begin{aligned} \bigg (\sum _{i=1}^m\frac{\lambda ^{\frac{N+2}{2}}}{\lambda |y|(1+\lambda |y|+\lambda |z-\xi _i|)^{\frac{N}{2}+\tau }}\bigg )^{-1}\le C \Big (\frac{m}{\lambda ^\tau }\Big )^{-1}\le C, \end{aligned}$$

which imply that

$$\begin{aligned} \bigg (\sum _{i=1}^m\frac{\lambda ^{\frac{N+2}{2}}}{\lambda |y|(1+\lambda |y|+\lambda |z-\xi _i|)^{\frac{N}{2}+\tau }}\bigg )^{-1}|J_0|\le C \Big (\frac{1}{\lambda }\Big )^{1+\varepsilon }. \end{aligned}$$

If \(|(|y|,r,z'')-(0,r_0,z_0'')|\le \delta \), then

$$\begin{aligned} |J_0|&\le \frac{C}{|y|}\bigg [U_{\xi _1,\lambda }^{2^{\star }-2}\sum _{i=2}^m U_{\xi _i,\lambda }+\Big (\sum _{i=2}^mU_{\xi _i,\lambda }\Big )^{2^{\star }-1}\bigg ]\\&\le \frac{C}{|y|}\bigg [\frac{\lambda }{(1+\lambda |y|+\lambda |z-\xi _1|)^2}\sum _{i=2}^m\frac{\lambda ^{\frac{N-2}{2}}}{(1+\lambda |y|+\lambda |z-\xi _i|)^{N-2}}\\&\quad +\,\Big (\sum _{i=2}^m\frac{\lambda ^{\frac{N-2}{2}}}{(1+\lambda |y|+\lambda |z-\xi _i|)^{N-2}}\Big )^{2^{\star }-1}\bigg ]. \end{aligned}$$

By (2.25), (2.12) and choosing \(\frac{N-2}{2}<\sigma \le \frac{N-2\tau }{2}\), we have

$$\begin{aligned}&\frac{1}{|y|}\frac{\lambda }{(1+\lambda |y|+\lambda |z-\xi _1|)^2}\sum _{i=2}^m\frac{\lambda ^{\frac{N-2}{2}}}{(1+\lambda |y|+\lambda |z-\xi _i|)^{N-2}}\\&\quad \le C\sum _{i=2}^m\frac{\lambda ^{\frac{N}{2}}}{|y|}\frac{1}{(1+\lambda |y|+\lambda |z-\xi _1|)^{\frac{N}{2}}}\frac{1}{(1+\lambda |y|+\lambda |z-\xi _i|)^{\frac{N}{2}}}\\&\quad \le C\sum _{i=2}^m\frac{1}{(\lambda |\xi _i-\xi _1|)^\sigma }\frac{\lambda ^{\frac{N}{2}}}{|y|(1+\lambda |y|+\lambda |z-\xi _1|)^{N-\sigma }}\\&\quad \le \frac{C}{\lambda ^{1+\varepsilon }}\sum _{i=1}^m\frac{\lambda ^{\frac{N+2}{2}}}{\lambda |y|(1+\lambda |y|+\lambda |z-\xi _i|)^{\frac{N}{2}+\tau }}. \end{aligned}$$

By Hölder inequality, (2.24) and (2.12), we get

$$\begin{aligned}&\frac{1}{|y|}\bigg (\sum _{i=2}^m\frac{\lambda ^{\frac{N-2}{2}}}{(1+\lambda |y|+\lambda |z-\xi _i|)^{N-2}}\bigg )^{2^{\star }-1}\\&\quad \le \frac{1}{|y|}\sum _{i=2}^m\frac{\lambda ^{\frac{N}{2}}}{(1+\lambda |y|+\lambda |z-\xi _i|)^{\frac{N}{2}+\tau }}\bigg (\sum _{i=2}^m\frac{1}{(1+\lambda |y|+\lambda |z-\xi _i|)^{\big (\frac{N-2}{2}-\frac{N-2}{N}\tau \big )\frac{N}{2}}}\bigg )^{\frac{2}{N-2}}\\&\quad \le \frac{C}{|y|}\sum _{i=2}^m\frac{\lambda ^{\frac{N}{2}}}{(1+\lambda |y|+\lambda |z-\xi _i|)^{\frac{N}{2}+\tau }}\bigg (\sum _{i=2}^m\frac{1}{(\lambda |\xi _i-\xi _1|)^{\big (\frac{N-2}{2}-\frac{N-2}{N}\tau \big )\frac{N}{2}}}\bigg )^{\frac{2}{N-2}}\\&\quad \le \frac{C}{\lambda ^{1+\varepsilon }}\sum _{i=1}^m\frac{\lambda ^{\frac{N+2}{2}}}{\lambda |y|(1+\lambda |y|+\lambda |z-\xi _i|)^{\frac{N}{2}+\tau }}. \end{aligned}$$

Hence

$$\begin{aligned} \Vert J_0\Vert _{**}\le C\Big (\frac{1}{\lambda }\Big )^{1+\varepsilon }. \end{aligned}$$
(2.26)

Secondly, we estimate \(J_1\). Due to (2.7) and

$$\begin{aligned} N-2-\varepsilon \ge \frac{N}{2}+\tau \quad \text {if }\varepsilon >0 \text { small}, \end{aligned}$$
(2.27)

we can obtain that

$$\begin{aligned} |J_1|&\le C\sum _{i=1}^m\frac{\eta \lambda ^{\frac{N-2}{2}}}{(1+\lambda |y|+\lambda |z-\xi _i|)^{N-2}} \le \frac{C}{\lambda ^{1+\varepsilon }}\sum _{i=1}^m\frac{\eta \lambda ^{\frac{N+2}{2}}}{\lambda |y|(1+\lambda |y|+\lambda |z-\xi _i|)^{N-2-\varepsilon }}\\&\le \frac{C}{\lambda ^{1+\varepsilon }}\sum _{i=1}^m\frac{\lambda ^{\frac{N+2}{2}}}{\lambda |y|(1+\lambda |y|+\lambda |z-\xi _i|)^{\frac{N}{2}+\tau }}. \end{aligned}$$

Thus

$$\begin{aligned} \Vert J_1\Vert _{**}\le C\Big (\frac{1}{\lambda }\Big )^{1+\varepsilon }. \end{aligned}$$
(2.28)

Thirdly, we estimate \(J_2\). By similar calculation to \(J_1\), we get

$$\begin{aligned} \Vert J_2\Vert _{**}\le C\Big (\frac{1}{\lambda }\Big )^{1+\varepsilon }. \end{aligned}$$
(2.29)

Fourthly, we estimate \(J_3\). Noting the fact that for \(x\in \text {supp}|\nabla \eta |\),

$$\begin{aligned} 1+\lambda |y|+\lambda |z-\xi _i|\ge C\lambda ,\quad i=1,2,\ldots ,m, \end{aligned}$$
(2.30)

and (2.27), we can deduce that

$$\begin{aligned} |J_3|&\le C\sum _{i=1}^m\frac{|\nabla \eta |\lambda ^{\frac{N}{2}}}{(1+\lambda |y|+\lambda |z-\xi _i|)^{N-1}} \le \frac{C}{\lambda ^{1+\varepsilon }}\sum _{i=1}^m\frac{|\nabla \eta |\lambda ^{\frac{N+2}{2}}}{\lambda |y|(1+\lambda |y|+\lambda |z-\xi _i|)^{N-2-\varepsilon }}\\&\le \frac{C}{\lambda ^{1+\varepsilon }}\sum _{i=1}^m\frac{\lambda ^{\frac{N+2}{2}}}{\lambda |y|(1+\lambda |y|+\lambda |z-\xi _i|)^{\frac{N}{2}+\tau }}. \end{aligned}$$

Therefore,

$$\begin{aligned} \Vert J_3\Vert _{**}\le C\Big (\frac{1}{\lambda }\Big )^{1+\varepsilon }. \end{aligned}$$
(2.31)

Combining (2.26), (2.28), (2.29) and (2.31), we finish the proof. \(\square \)

Now we prove Proposition 2.3 by the Contraction Mapping Principle.

Proof of Proposition 2.3

Set

$$\begin{aligned} \mathcal {N}=\Big \{w\in C(\mathbb {R}^N)\cap H_s:\Vert w\Vert _{*}\le \frac{1}{\lambda },\ \int \frac{\bar{U}_{\xi _i,\lambda }^{2^\star -2}}{|y|}\bar{U}_{i,l}w=0, \\ l=1,2,\ldots ,N-k, \ i=1,2,\ldots ,m\Big \}. \end{aligned}$$

Problem (2.23) is equivalent to

$$\begin{aligned} \varphi =\mathcal {A}(\varphi ):=L_m(\mathcal {R}(\varphi ))+L_m(\mathcal P_m), \end{aligned}$$

where \(L_m\) is defined in Lemma 2.2. We will prove that \(\mathcal {A}\) is a contraction mapping from \(\mathcal {N}\) to \(\mathcal {N}\).

On the one hand, by (2.20), Lemmas 2.4 and 2.5, we have

$$\begin{aligned} \Vert \mathcal {A}(\varphi )\Vert _*&\le C\big (\Vert \mathcal {R}(\varphi )\Vert _{**}+\Vert \mathcal P_m\Vert _{**}\big )\le C\Big [\Vert \varphi \Vert _*^{2^\star -1}+\Big (\frac{1}{\lambda }\Big )^{1+\varepsilon }\Big ]\\&\le C\Big [\Big (\frac{1}{\lambda }\Big )^{2^\star -1}+\Big (\frac{1}{\lambda }\Big )^{1+\varepsilon }\Big ]\le \frac{C}{\lambda ^{1+\varepsilon }}\le \frac{1}{\lambda },\quad \forall \varphi \in \mathcal {N}. \end{aligned}$$

Hence \(\mathcal {A}\) maps \(\mathcal {N}\) to \(\mathcal {N}\).

On the other hand,

$$\begin{aligned} \Vert \mathcal {A}(\varphi _1)-\mathcal {A}(\varphi _2)\Vert _*=\Vert L_m(\mathcal {R}(\varphi _1))-L_m(\mathcal {R}(\varphi _2))\Vert _*\le C\Vert \mathcal {R}(\varphi _1)-\mathcal {R}(\varphi _2)\Vert _{**}. \end{aligned}$$

Since \(|\mathcal {R}'(\varphi )|\le C\frac{|\varphi |^{2^\star -2}}{|y|}\) and

$$\begin{aligned} \bigg (\sum _{i=1}^m\frac{\lambda ^{\frac{N-2}{2}}}{(1+\lambda |y|+\lambda |z-\xi _i|)^{\frac{N-2}{2}+\tau }}\bigg )^{2^\star -1}\le \sum _{i=1}^m\frac{C\lambda ^{\frac{N}{2}}}{(1+\lambda |y|+\lambda |z-\xi _i|)^{\frac{N}{2}+\tau }}, \end{aligned}$$
(2.32)

we get that there is \(0<\vartheta <1\) such that

$$\begin{aligned}&|\mathcal {R}(\varphi _1)-\mathcal {R}(\varphi _2)|\le |\mathcal {R}'(\varphi _1+\vartheta (\varphi _2-\varphi _1))|\cdot |\varphi _1-\varphi _2|\\&\quad \le C\frac{|\varphi _1|^{2^\star -2}+|\varphi _2|^{2^\star -2}}{|y|}|\varphi _1-\varphi _2|\\&\quad \le C\big (\Vert \varphi _1\Vert _*^{2^\star -2}+\Vert \varphi _2\Vert _*^{2^\star -2}\big )\Vert \varphi _1-\varphi _2\Vert _* \cdot \frac{1}{|y|}\bigg (\sum _{i=1}^m\frac{\lambda ^{\frac{N-2}{2}}}{(1+\lambda |y|+\lambda |z-\xi _i|)^{\frac{N-2}{2}+\tau }}\bigg )^{2^\star -1}\\&\quad \le C\Big (\frac{1}{\lambda }\Big )^{2^\star -2}\Vert \varphi _1-\varphi _2\Vert _*\sum _{i=1}^m\frac{\lambda ^{\frac{N+2}{2}}}{\lambda |y|(1+\lambda |y|+\lambda |z-\xi _i|)^{\frac{N}{2}+\tau }}, \end{aligned}$$

which implies that

$$\begin{aligned} \Vert \mathcal {R}(\varphi _1)-\mathcal {R}(\varphi _2)\Vert _{**}\le C\Big (\frac{1}{\lambda }\Big )^{2^\star -2}\Vert \varphi _1-\varphi _2\Vert _*=o(1)\Vert \varphi _1-\varphi _2\Vert _*. \end{aligned}$$

Thus \(\mathcal {A}\) is a contraction mapping on \(\mathcal {N}\).

In virtue of the Contraction Mapping Principle, there exists an unique \(\varphi =\varphi _{\bar{r},\bar{z}'',\lambda }\in \mathcal {N}\) satisfying problem (2.23). Moreover, by (2.20), Lemmas 2.4 and 2.5 again, we achieve (2.22). \(\square \)

3 Proof of the main result

In this section, we prove the existence of infinitely many positive multi-bubbling solutions for problem (1.3). Namely, we choose suitable \((\bar{r}, \bar{z}'',\lambda )\) satisfying \(c_l=0\ ( l=1,2,\ldots ,N-k)\), and thus \(\bar{W}_{\bar{r},\bar{z}'',\lambda }+\varphi _{\bar{r},\bar{z}'',\lambda }\) is a solution of problem (1.3).

For this purpose, the following result is needed. We point out in advance that the left-hand side of (3.1) is the local Pohozaev identity generating from scaling, while the left-hand side of (3.2) is the local Pohozaev identities generating from translations.

Proposition 3.1

Suppose that \((\bar{r},\bar{z}'',\lambda )\) satisfies

$$\begin{aligned}&\int _{D_\rho }\Big (-\Delta u_m+V(r,z'')u_m-\frac{(u_m)_+^{2^\star -1}}{|y|}\Big )\left\langle x,\nabla u_m\right\rangle =0, \end{aligned}$$
(3.1)
$$\begin{aligned}&\int _{D_\rho }\Big (-\Delta u_m+V(r,z'')u_m-\frac{(u_m)_+^{2^\star -1}}{|y|}\Big )\frac{\partial u_m}{\partial z_j}=0,\quad j=3,4,\ldots ,N-k, \end{aligned}$$
(3.2)
$$\begin{aligned}&\int _{\mathbb {R}^N}\Big (-\Delta u_m+V(r,z'')u_m-\frac{(u_m)_+^{2^\star -1}}{|y|}\Big )\frac{\partial \bar{W}_{\bar{r},\bar{z}'',\lambda }}{\partial \lambda }=0, \end{aligned}$$
(3.3)

where \(D_\rho :=\big \{(y,z',z'')\in \mathbb {R}^N:|(|y|,|z'|,z'')-(0,r_0,z_0'')|\le \rho \big \}\) with \(\rho \in (2\delta ,5\delta )\) and \(u_m:=\bar{W}_{\bar{r},\bar{z}'',\lambda }+\varphi _{\bar{r},\bar{z}'',\lambda }\) gotten from Proposition 2.3. Then

$$\begin{aligned} c_l=0,\quad l=1,2,\ldots ,N-k. \end{aligned}$$

Proof

Since \(\text {supp}\eta \subset D_\rho \), we find \(\bar{U}_{\xi _i,\lambda }=0\) in \(\mathbb {R}^N{\setminus } D_\rho \) for \(i=1,2,\ldots ,m\). It follows from (3.1), (3.2), (3.3) that

$$\begin{aligned} \sum _{l=1}^{N-k}c_l\sum _{i=1}^m\int _{\mathbb {R}^N}\frac{\bar{U}_{\xi _i,\lambda }^{2^\star -2}}{|y|}\bar{U}_{i,l}v=0 \end{aligned}$$
(3.4)

for \(v=\left\langle x,\nabla u_m\right\rangle , \ v=\frac{\partial u_m}{\partial z_j} \ (j=3,4,\ldots ,N-k)\) and \(v=\frac{\partial \bar{W}_{\bar{r},\bar{z}'',\lambda }}{\partial \lambda }\), respectively.

By direct computations, we can prove

$$\begin{aligned}&\sum _{i=1}^m\int _{\mathbb {R}^N}\frac{\bar{U}_{\xi _i,\lambda }^{2^\star -2}}{|y|}\bar{U}_{i,2}\left\langle z',\nabla _{z'}\bar{W}_{\bar{r},\bar{z}'',\lambda }\right\rangle =m\lambda ^2[a_1+o(1)], \end{aligned}$$
(3.5)
$$\begin{aligned}&\sum _{i=1}^m\int _{\mathbb {R}^N}\frac{\bar{U}_{\xi _i,\lambda }^{2^\star -2}}{|y|}\bar{U}_{i,j}\frac{\partial \bar{W}_{\bar{r},\bar{z}'',\lambda }}{\partial z_j}=m\lambda ^2[a_2+o(1)],\quad j=3,4,\ldots ,N-k,\end{aligned}$$
(3.6)
$$\begin{aligned}&\sum _{i=1}^m\int _{\mathbb {R}^N}\frac{\bar{U}_{\xi _i,\lambda }^{2^\star -2}}{|y|}\bar{U}_{i,1}\frac{\partial \bar{W}_{\bar{r},\bar{z}'',\lambda }}{\partial \lambda }=\frac{m}{\lambda ^2}[a_3+o(1)] \end{aligned}$$
(3.7)

for some constants \(a_1>0,\ a_2<0,\ a_3>0\).

By means of integration by parts and (2.22), we get

$$\begin{aligned} \sum _{l=1}^{N-k}c_l\sum _{i=1}^m\int _{\mathbb {R}^N}\frac{\bar{U}_{\xi _i,\lambda }^{2^\star -2}}{|y|}\bar{U}_{i,l}v=o(m|c_1|)+o(m\lambda ^2)\sum _{l=2}^{N-k}|c_l| \end{aligned}$$

for \(v=\left\langle x,\nabla \varphi _{\bar{r},\bar{z}'',\lambda }\right\rangle \) and \(v=\frac{\partial \varphi _{\bar{r},\bar{z}'',\lambda }}{\partial z_j} \ (j=3,4,\ldots ,N-k)\), respectively, which together with (3.4) imply that

$$\begin{aligned} \sum _{l=1}^{N-k}c_l\sum _{i=1}^m\int _{\mathbb {R}^N}\frac{\bar{U}_{\xi _i,\lambda }^{2^\star -2}}{|y|}\bar{U}_{i,l}v=o(m|c_1|)+o(m\lambda ^2)\sum _{l=2}^{N-k}|c_l| \end{aligned}$$
(3.8)

for \(v=\left\langle x,\nabla \bar{W}_{\bar{r},\bar{z}'',\lambda }\right\rangle \) and \(v=\frac{\partial \bar{W}_{\bar{r},\bar{z}'',\lambda }}{\partial z_j} \ (j=3,4,\ldots ,N-k)\), respectively.

Noting

$$\begin{aligned} \left\langle x,\nabla \bar{W}_{\bar{r},\bar{z}'',\lambda }\right\rangle =\left\langle y,\nabla _y\bar{W}_{\bar{r},\bar{z}'',\lambda }\right\rangle +\left\langle z',\nabla _{z'}\bar{W}_{\bar{r},\bar{z}'',\lambda }\right\rangle +\left\langle z'',\nabla _{z''}\bar{W}_{\bar{r},\bar{z}'',\lambda }\right\rangle , \end{aligned}$$

we can obtain

$$\begin{aligned}&\sum _{l=1}^{N-k}c_l\sum _{i=1}^m\int _{\mathbb {R}^N}\frac{\bar{U}_{\xi _i,\lambda }^{2^\star -2}}{|y|}\bar{U}_{i,l}\left\langle x,\nabla \bar{W}_{\bar{r},\bar{z}'',\lambda }\right\rangle \nonumber \\&\quad =o(m|c_1|)+c_2\sum _{i=1}^m\int _{\mathbb {R}^N}\frac{\bar{U}_{\xi _i,\lambda }^{2^\star -2}}{|y|}\bar{U}_{i,2}\left\langle z',\nabla _{z'} \bar{W}_{\bar{r},\bar{z}'',\lambda }\right\rangle \nonumber \\&\qquad +\,o(m\lambda ^2)|c_2|+\sum _{l=3}^{N-k}c_l[b_l+o(1)]m\lambda ^2,\quad b_l\in \mathbb {R}. \end{aligned}$$
(3.9)

In addition, for \(j=3,4,\ldots ,N-k\),

$$\begin{aligned}&\sum _{l=1}^{N-k}c_l\sum _{i=1}^m\int _{\mathbb {R}^N}\frac{\bar{U}_{\xi _i,\lambda }^{2^\star -2}}{|y|}\bar{U}_{i,l}\frac{\partial \bar{W}_{\bar{r},\bar{z}'',\lambda }}{\partial z_j}\nonumber \\&\quad =o(m|c_1|)+c_j\sum _{i=1}^m\int _{\mathbb {R}^N}\frac{\bar{U}_{\xi _i,\lambda }^{2^\star -2}}{|y|}\bar{U}_{i,j}\frac{\partial \bar{W}_{\bar{r},\bar{z}'',\lambda }}{\partial z_j}+o(m\lambda ^2)\sum _{l\ne 1,j}^{N-k}|c_l|. \end{aligned}$$
(3.10)

By (3.8), (3.9) and (3.5), we have

$$\begin{aligned} c_2[a_1+o(1)]=\sum _{l=3}^{N-k}c_l[b_l+o(1)]+o\Big (\frac{1}{\lambda ^2}\Big )c_1. \end{aligned}$$
(3.11)

By (3.8), (3.10) and (3.6), we get for \(j=3,4,\ldots ,N-k\),

$$\begin{aligned} c_j[a_2+o(1)]=o(1)\sum _{l\ne 1,j}^{N-k}c_l+o\Big (\frac{1}{\lambda ^2}\Big )c_1. \end{aligned}$$
(3.12)

It follows from (3.11) and (3.12) that

$$\begin{aligned} c_j=o\Big (\frac{1}{\lambda ^2}\Big )c_1,\quad j=2,3,\ldots ,N-k. \end{aligned}$$
(3.13)

In terms of (3.4), (3.7), (3.13) and

$$\begin{aligned} \sum _{i=1}^m\int _{\mathbb {R}^N}\frac{\bar{U}_{\xi _i,\lambda }^{2^\star -2}}{|y|}\bar{U}_{i,j}\frac{\partial \bar{W}_{\bar{r},\bar{z}'',\lambda }}{\partial \lambda }=mo(1),\quad j=2,3,\ldots ,N-k, \end{aligned}$$

we obtain that

$$\begin{aligned} 0&=\sum _{l=1}^{N-k}c_l\sum _{i=1}^m\int _{\mathbb {R}^N}\frac{\bar{U}_{\xi _i,\lambda }^{2^\star -2}}{|y|}\bar{U}_{i,l}\frac{\partial \bar{W}_{\bar{r},\bar{z}'',\lambda }}{\partial \lambda }\\&=c_1\sum _{i=1}^m\int _{\mathbb {R}^N}\frac{\bar{U}_{\xi _i,\lambda }^{2^\star -2}}{|y|}\bar{U}_{i,1}\frac{\partial \bar{W}_{\bar{r},\bar{z}'',\lambda }}{\partial \lambda }+\sum _{l=2}^{N-k}c_l\sum _{i=1}^m\int _{\mathbb {R}^N}\frac{\bar{U}_{\xi _i,\lambda }^{2^\star -2}}{|y|}\bar{U}_{i,l}\frac{\partial \bar{W}_{\bar{r},\bar{z}'',\lambda }}{\partial \lambda }\\&=c_1\frac{m}{\lambda ^2}[a_3+o(1)]+o\Big (\frac{1}{\lambda ^2}\Big )c_1\cdot mo(1) =\Big [\frac{ma_3}{\lambda ^2}+o\Big (\frac{m}{\lambda ^2}\Big )\Big ]c_1, \end{aligned}$$

which implies \(c_1=0\). \(\square \)

Next we will prove there really exists \((\bar{r},\bar{z}'',\lambda )\) satisfying (3.1), (3.2) and (3.3). To achieve it, we have to find the equivalent forms of (3.1), (3.2) and (3.3), respectively.

Now we estimate (3.3).

Lemma 3.2

Equation (3.3) is equivalent to

$$\begin{aligned} m\bigg (-\frac{B_1}{\lambda ^3}V(\bar{r},\bar{z}'')+\frac{B_3m^{N-2}}{\lambda ^{N-1}}+O\Big (\frac{1}{\lambda ^{3+\varepsilon }}\Big )\bigg )=0, \end{aligned}$$

where \(B_1>0,\ B_3>0\) are defined in Lemma A.1.

Proof

Since \(u_m=\bar{W}_{\bar{r},\bar{z}'',\lambda }+\varphi _{\bar{r},\bar{z}'',\lambda }\), we have

$$\begin{aligned}&\int \bigg [-\Delta u_m+V(r,z'')u_m-\frac{(u_m)_+^{2^\star -1}}{|y|}\bigg ]\frac{\partial \bar{W}_{\bar{r},\bar{z}'',\lambda }}{\partial \lambda }\\&\quad =\big<I'(\bar{W}_{\bar{r},\bar{z}'',\lambda }),\frac{\partial \bar{W}_{\bar{r},\bar{z}'',\lambda }}{\partial \lambda }\big>-\int \mathcal {R}(\varphi )\frac{\partial \bar{W}_{\bar{r},\bar{z}'',\lambda }}{\partial \lambda }\\&\qquad +\int \bigg [-\Delta \varphi +V(r,z'')\varphi -(2^\star -1)\frac{\bar{W}_{\bar{r},\bar{z}'',\lambda }^{2^\star -2}}{|y|}\varphi \bigg ]\frac{\partial \bar{W}_{\bar{r},\bar{z}'',\lambda }}{\partial \lambda }\\&\quad :=\big <I'(\bar{W}_{\bar{r},\bar{z}'',\lambda }),\frac{\partial \bar{W}_{\bar{r},\bar{z}'',\lambda }}{\partial \lambda }\big >-I_1+I_2. \end{aligned}$$

Taking into account (2.32), (2.9), (2.22) and for \(\theta \in (1,2]\),

$$\begin{aligned} (1+t)_+^\theta -1-\theta t=O(t^2),\quad \forall t\in \mathbb {R}, \end{aligned}$$

we can obtain

$$\begin{aligned} |I_1|&\le C\int \frac{\bar{W}_{\bar{r},\bar{z}'',\lambda }^{2^\star -3}}{|y|}\varphi ^2\Big |\frac{\partial \bar{W}_{\bar{r},\bar{z}'',\lambda }}{\partial \lambda }\Big |\\&\le C\frac{\Vert \varphi \Vert _*^2}{\lambda }\int \frac{1}{|y|}\bigg (\sum _{i=1}^m\frac{\lambda ^{\frac{N-2}{2}}}{(1+\lambda |y|+\lambda |z-\xi _i|)^{\frac{N-2}{2}+\tau }}\bigg )^{2^\star }\\&\le C\frac{\Vert \varphi \Vert _*^2}{\lambda }\int \frac{1}{|y|}\sum _{i=1}^m\frac{\lambda ^\frac{N}{2}}{(1+\lambda |y|+\lambda |z-\xi _i|)^{\frac{N}{2}+\tau }}\sum _{i=1}^m\frac{\lambda ^{\frac{N-2}{2}}}{(1+\lambda |y|+\lambda |z-\xi _i|)^{\frac{N-2}{2}+\tau }}\\&\le C\frac{m\Vert \varphi \Vert _*^2}{\lambda }\int \frac{1}{|y|}\frac{\lambda ^\frac{N}{2}}{(1+\lambda |y|+\lambda |z-\xi _1|)^{\frac{N}{2}+\tau }}\sum _{i=1}^m\frac{\lambda ^{\frac{N-2}{2}}}{(1+\lambda |y|+\lambda |z-\xi _i|)^{\frac{N-2}{2}+\tau }}\\&\le C\frac{m\Vert \varphi \Vert _*^2}{\lambda }\Big (C+\sum _{i=2}^m\frac{C}{(\lambda |\xi _i-\xi _1|)^\tau }\Big )\le C\frac{m}{\lambda ^{3+\varepsilon }}. \end{aligned}$$

By (2.11), (2.13) and (2.22), we get

$$\begin{aligned} I_2=O\Big (\frac{m\Vert \varphi \Vert _*}{\lambda ^{2+\varepsilon }}\Big )=O\Big (\frac{m}{\lambda ^{3+\varepsilon }}\Big ). \end{aligned}$$

Therefore,

$$\begin{aligned} \big<I'(\bar{W}_{\bar{r},\bar{z}'',\lambda }+\varphi ),\frac{\partial \bar{W}_{\bar{r},\bar{z}'',\lambda }}{\partial \lambda }\big>=\big <I'(\bar{W}_{\bar{r},\bar{z}'',\lambda }),\frac{\partial \bar{W}_{\bar{r},\bar{z}'',\lambda }}{\partial \lambda }\big >+O\Big (\frac{m}{\lambda ^{3+\varepsilon }}\Big ). \end{aligned}$$

By Lemma A.1 in “Appendix A”, the proof is completed. \(\square \)

Now we estimate (3.1) and (3.2).

Lemma 3.3

Equation (3.2) is equivalent to

$$\begin{aligned} \int _{D_\rho }\frac{\partial V(r,z'')}{\partial z_j}u_m^2=O\bigg (\int _{\partial D_\rho }\Big (|\nabla \varphi |^2+\varphi ^2+\frac{|\varphi |^{2^\star }}{|y|}\Big )\bigg ),\ j=3,4,\ldots ,N-k. \end{aligned}$$

Lemma 3.4

Equation (3.1) is equivalent to

$$\begin{aligned} \int _{D_\rho }\frac{1}{2r}\frac{\partial (r^2V(r,z''))}{\partial r}u_m^2=o\Big (\frac{m}{\lambda ^2}\Big )+O\bigg (\int _{\partial D_\rho }\Big (|\nabla \varphi |^2+\varphi ^2+\frac{|\varphi |^{2^\star }}{|y|}\Big )\bigg ). \end{aligned}$$

Proof

Firstly, we have

$$\begin{aligned} \int _{D_\rho }\Big (V(z)+\frac{1}{2}\big <z,\nabla _z V(z)\big >\Big )u_m^2&= -\frac{N-2}{2}\int _{D_\rho }\Big (|\nabla u_m|^2+V(r,z'')u_m^2-\frac{(u_m)_+^{2^\star }}{|y|}\Big )\nonumber \\&\quad +\,O\Big (\int _{\partial D_\rho }\Big (|\nabla \varphi |^2+\varphi ^2+\frac{|\varphi |^{2^\star }}{|y|}\Big )\Big ). \end{aligned}$$
(3.14)

Secondly, by (2.21), we get

$$\begin{aligned}&\int _{D_\rho }\Big (|\nabla u_m|^2+V(r,z'')u_m^2-\frac{(u_m)_+^{2^\star }}{|y|}\Big )\nonumber \\&\quad =\sum _{l=1}^{N-k}c_l\sum _{i=1}^m\int _{\mathbb {R}^N}\frac{\bar{U}_{\xi _i,\lambda }^{2^{\star }-2}}{|y|}\bar{U}_{i,l}u_m+O\Big (\int _{\partial D_\rho }\big (|\nabla \varphi |^2+\varphi ^2\big )\Big )\nonumber \\&\quad =\sum _{l=1}^{N-k}c_l\sum _{i=1}^m\int _{\mathbb {R}^N}\frac{\bar{U}_{\xi _i,\lambda }^{2^{\star }-2}}{|y|}\bar{U}_{i,l}\bar{W}_{\bar{r},\bar{z}'',\lambda }+O\Big (\int _{\partial D_\rho }\big (|\nabla \varphi |^2+\varphi ^2\big )\Big ). \end{aligned}$$
(3.15)

Inserting (3.15) into (3.14), we obtain

$$\begin{aligned} \int _{D_\rho }\Big (V(z)+\frac{1}{2}\big <z,\nabla _z V(z)\big >\Big )u_m^2&= -\frac{N-2}{2}\sum _{l=1}^{N-k}c_l\sum _{i=1}^m\int _{\mathbb {R}^N}\frac{\bar{U}_{\xi _i,\lambda }^{2^{\star }-2}}{|y|}\bar{U}_{i,l}\bar{W}_{\bar{r},\bar{z}'',\lambda }\nonumber \\&\quad +\,O\Big (\int _{\partial D_\rho }\Big (|\nabla \varphi |^2+\varphi ^2+\frac{|\varphi |^{2^\star }}{|y|}\Big )\Big ). \end{aligned}$$
(3.16)

It follows from (3.16), (2.22) and

$$\begin{aligned} \sum _{i=1}^m\int _{\mathbb {R}^N}\frac{\bar{U}_{\xi _i,\lambda }^{2^{\star }-2}}{|y|}\bar{U}_{i,l}\bar{W}_{\bar{r},\bar{z}'',\lambda }&=m\bigg (\int _{\mathbb {R}^N}\frac{\bar{U}_{\xi _1,\lambda }^{2^{\star }-1}}{|y|}\bar{U}_{1,l}+\int _{\mathbb {R}^N}\frac{\bar{U}_{\xi _1,\lambda }^{2^{\star }-2}}{|y|}\bar{U}_{1,l}\sum _{i=2}^m\bar{U}_{\xi _i,\lambda }\bigg )\\&=m\Big (O\Big (\frac{\lambda ^{n_l}}{\lambda ^{N-1}}\Big )+O\Big (\frac{\lambda ^{n_l}}{\lambda ^2}\Big )\Big )=O\Big (\frac{m\lambda ^{n_l}}{\lambda ^2}\Big ) \end{aligned}$$

that for some small \(\varepsilon >0\),

$$\begin{aligned} \int _{D_\rho }\Big (V(z)+\frac{1}{2}\big <z,\nabla _z V(z)\big >\Big )u_m^2 =O\Big (\frac{m}{\lambda ^{2+\varepsilon }}\Big )+O\Big (\int _{\partial D_\rho }\Big (|\nabla \varphi |^2+\varphi ^2+\frac{|\varphi |^{2^\star }}{|y|}\Big )\Big ).\qquad \end{aligned}$$
(3.17)

Noting

$$\begin{aligned} \left\langle z,\nabla _z V(z)\right\rangle =r\frac{\partial V(r,z'')}{\partial r}+\sum _{j=3}^{N-k}z_j\frac{\partial V(r,z'')}{\partial z_j} \end{aligned}$$

and by Lemma 3.3 and (3.17), we have

$$\begin{aligned} \int _{D_\rho }\Big (V(r,z'')+\frac{1}{2}r\frac{\partial V(r,z'')}{\partial r}\Big )u_m^2=o\Big (\frac{m}{\lambda ^2}\Big )+O\Big (\int _{\partial D_\rho }\Big (|\nabla \varphi |^2+\varphi ^2+\frac{|\varphi |^{2^\star }}{|y|}\Big )\Big ). \end{aligned}$$

\(\square \)

To simplify the results in Lemmas 3.3 and 3.4 and obtain the final equivalent forms of (3.1) and (3.2), we need the following several results.

Lemma 3.5

There holds

$$\begin{aligned} \int _{\mathbb {R}^N}\Big (|\nabla \varphi |^2+V(r,z'')\varphi ^2\Big )=O\Big (\frac{m}{\lambda ^{2+\varepsilon }}\Big ). \end{aligned}$$

Proof

By (2.21), we have

$$\begin{aligned}&\int _{\mathbb {R}^N}\Big (|\nabla \varphi |^2+V(r,z'')\varphi ^2\Big )\\&\quad =\int _{\mathbb {R}^N}\frac{1}{|y|}\Big [\big (\bar{W}_{\bar{r},\bar{z}'',\lambda }+\varphi \big )_+^{2^\star -1}-\bar{W}_{\bar{r},\bar{z}'',\lambda }^{2^\star -1}\Big ]\varphi \\&\qquad -\int _{\mathbb {R}^N} V(r,z'')\bar{W}_{\bar{r},\bar{z}'',\lambda }\varphi +\int _{\mathbb {R}^N}\bigg (\frac{\bar{W}_{\bar{r},\bar{z}'',\lambda }^{2^\star -1}}{|y|}+\Delta \bar{W}_{\bar{r},\bar{z}'',\lambda }\bigg )\varphi \\&\quad :=\tilde{J}_1+\tilde{J}_2+\tilde{J}_3. \end{aligned}$$

Similarly to the arguments of \(I_1\) in Lemma 3.2 and applying (2.22), we have

$$\begin{aligned} |\tilde{J}_1|&\le C\int \frac{1}{|y|}\Big (\bar{W}_{\bar{r},\bar{z}'',\lambda }^{2^\star -2}\varphi ^2+|\varphi |^{2^\star }\Big )\\&\le C\big (\Vert \varphi \Vert _*^2+\Vert \varphi \Vert _*^{2^\star }\big )\int \frac{1}{|y|}\bigg (\sum _{i=1}^m\frac{\lambda ^{\frac{N-2}{2}}}{(1+\lambda |y|+\lambda |z-\xi _i|)^{\frac{N-2}{2}+\tau }}\bigg )^{2^\star }\\&\le C\frac{m}{\lambda ^{2+\varepsilon }}. \end{aligned}$$

As similar as the arguments of (2.11) and by (2.22), we get

$$\begin{aligned} |\tilde{J}_2|&\le C\Vert \varphi \Vert _*\int \sum _{i=1}^m\frac{\eta \lambda ^{\frac{N-2}{2}}}{(1+\lambda |y|+\lambda |z-\xi _i|)^{N-2}}\sum _{i=1}^m\frac{\lambda ^{\frac{N-2}{2}}}{(1+\lambda |y|+\lambda |z-\xi _i|)^{\frac{N-2}{2}+\tau }}\\&\le C\frac{m\Vert \varphi \Vert _*}{\lambda ^{1+\varepsilon }}\int \frac{\eta \lambda ^{\frac{N+2}{2}}}{(1+\lambda |y|+\lambda |z-\xi _1|)^{\frac{N+2}{2}+\tau }}\sum _{i=1}^m\frac{\lambda ^{\frac{N-2}{2}}}{(1+\lambda |y|+\lambda |z-\xi _i|)^{\frac{N-2}{2}+\tau }}\\&\le C\frac{m}{\lambda ^{2+\varepsilon }}. \end{aligned}$$

By using the estimates of \(J_0,J_2,J_3\) in Lemma 2.5, (2.9) and (2.22), we obtain

$$\begin{aligned} |\tilde{J}_3|&= \bigg |\int \bigg [\bigg (\frac{\bar{W}_{\bar{r},\bar{z}'',\lambda }^{2^\star -1}}{|y|}-\sum _{i=1}^m\eta \frac{U_{\xi _i,\lambda }^{2^\star -1}}{|y|}\bigg )+\Delta \eta W_{\bar{r},\bar{z}'',\lambda }+2\nabla \eta \nabla W_{\bar{r},\bar{z}'',\lambda }\bigg ]\varphi \bigg |\\&\le C\frac{\Vert \varphi \Vert _*}{\lambda ^{1+\varepsilon }}\int \sum _{i=1}^m\frac{\lambda ^{\frac{N}{2}}}{|y|(1+\lambda |y|+\lambda |z-\xi _i|)^{\frac{N}{2}+\tau }}\sum _{i=1}^m\frac{\lambda ^{\frac{N-2}{2}}}{(1+\lambda |y|+\lambda |z-\xi _i|)^{\frac{N-2}{2}+\tau }}\\&\le C\frac{m}{\lambda ^{2+\varepsilon }}. \end{aligned}$$

\(\square \)

Lemma 3.6

There holds

$$\begin{aligned} \int _{D_{5\delta }{\setminus } D_{2\delta }}\Big (|\nabla \varphi |^2+\varphi ^2+\frac{|\varphi |^{2^\star }}{|y|}\Big )=O\Big (\frac{m}{\lambda ^{2+\varepsilon }}\Big ). \end{aligned}$$

Proof

By Hardy–Sobolev inequality and Lemma 3.5, we have

$$\begin{aligned} \int _{\mathbb {R}^N}\frac{|\varphi |^{2^\star }}{|y|}\le \Big (\int _{\mathbb {R}^N}|\nabla \varphi |^2\Big )^{\frac{2^\star }{2}}\le C\Big (\frac{m}{\lambda ^{2+\varepsilon }}\Big )^{\frac{2^\star }{2}}\le \frac{Cm}{\lambda ^{2+\varepsilon }}. \end{aligned}$$

Hence

$$\begin{aligned}&\int _{D_{5\delta }{\setminus } D_{2\delta }}\Big (|\nabla \varphi |^2+\varphi ^2+\frac{|\varphi |^{2^\star }}{|y|}\Big )\\&\quad \le C\int _{D_{5\delta }{\setminus } D_{2\delta }}\Big (|\nabla \varphi |^2+V(r,z'')\varphi ^2+\frac{|\varphi |^{2^\star }}{|y|}\Big )\\&\quad \le C\int _{\mathbb {R}^N}\Big (|\nabla \varphi |^2+V(r,z'')\varphi ^2+\frac{|\varphi |^{2^\star }}{|y|}\Big ) \le \frac{Cm}{\lambda ^{2+\varepsilon }}. \end{aligned}$$

\(\square \)

In view of Lemma 3.6, the following result is true.

Lemma 3.7

There exists \(\rho \in (3\delta ,4\delta )\) such that

$$\begin{aligned} \int _{\partial D_\rho }\Big (|\nabla \varphi |^2+\varphi ^2+\frac{|\varphi |^{2^\star }}{|y|}\Big )=o\Big (\frac{m}{\lambda ^2}\Big ). \end{aligned}$$

Lemma 3.8

For any \(g(r,z'')\in C(\mathbb {R}^{N-k-1},\mathbb {R})\), there holds

$$\begin{aligned} \int _{D_\rho }g(r,z'')u_m^2=m\bigg (\frac{1}{\lambda ^2}g(\bar{r},\bar{z}'')\int _{\mathbb {R}^N}U_{0,1}^2+o\Big (\frac{1}{\lambda ^2}\Big )\bigg ). \end{aligned}$$

Proof

Although the proof is similar to that of Lemma 3.4 in [22], we give the details for completeness. Since \(u_m=\bar{W}_{\bar{r},\bar{z}'',\lambda }+\varphi \), we have

$$\begin{aligned} \int _{D_\rho }g(r,z'')u_m^2=\int _{D_\rho }g(r,z'')\bar{W}_{\bar{r},\bar{z}'',\lambda }^2+2\int _{D_\rho }g(r,z'')\bar{W}_{\bar{r},\bar{z}'',\lambda }\varphi +\int _{D_\rho }g(r,z'')\varphi ^2.\qquad \end{aligned}$$
(3.18)

It is similar to the estimate of \(\tilde{J}_2\) in Lemma 3.5 that

$$\begin{aligned} \Big |2\int _{D_\rho }g(r,z'')\bar{W}_{\bar{r},\bar{z}'',\lambda }\varphi \Big |\le C\int _{D_\rho }\bar{W}_{\bar{r},\bar{z}'',\lambda }|\varphi |\le \frac{Cm\Vert \varphi \Vert _*}{\lambda ^{1+\varepsilon }}\le \frac{Cm}{\lambda ^{2+\varepsilon }}. \end{aligned}$$
(3.19)

By Lemma 3.5, we have

$$\begin{aligned} \Big |\int _{D_\rho }g(r,z'')\varphi ^2\Big |\le C\int _{D_\rho }V(r,z'')|\varphi |^2\le C\int _{\mathbb {R}^N}V(r,z'')|\varphi |^2\le \frac{Cm}{\lambda ^{2+\varepsilon }}. \end{aligned}$$
(3.20)

Inserting (3.19) and (3.20) into (3.18), we obtain

$$\begin{aligned} \int _{D_\rho }g(r,z'')u_m^2&= \int _{D_\rho }g(r,z'')\bar{W}_{\bar{r},\bar{z}'',\lambda }^2+O\Big (\frac{m}{\lambda ^{2+\varepsilon }}\Big )\\&= \int _{D_\rho }g(r,z'')\bigg (\sum _{i=1}^m\bar{U}_{\xi _i,\lambda }^2+\sum _{i=1}^m\sum _{j\ne i}^m\bar{U}_{\xi _i,\lambda }\bar{U}_{\xi _j,\lambda }\bigg )+o\Big (\frac{m}{\lambda ^2}\Big ), \end{aligned}$$

where

$$\begin{aligned} \int _{D_\rho }g(r,z'')\bar{U}_{\xi _i,\lambda }^2=\frac{1}{\lambda ^2}g(\bar{r},\bar{z}'')\int _{\mathbb {R}^N}U_{0,1}^2+o\Big (\frac{1}{\lambda ^2}\Big ) \end{aligned}$$

and

$$\begin{aligned} \bigg |\sum _{j\ne i}^m\int _{D_\rho }g(r,z'')\bar{U}_{\xi _i,\lambda }\bar{U}_{\xi _j,\lambda }\bigg |\le \frac{C}{\lambda ^2}\sum _{j\ne i}^m\frac{1}{(\lambda |\xi _j-\xi _i|)^{\frac{N-4}{2}}} =o\Big (\frac{1}{\lambda ^2}\Big ). \end{aligned}$$

\(\square \)

Choosing \(g(r,z'')=\frac{\partial V(r,z'')}{\partial z_j}\) and \(g(r,z'')=\frac{1}{2r}\frac{\partial (r^2 V(r,z''))}{\partial r}\) in Lemma 3.8, respectively, we have the following result.

Remark 3.9

For \(j=3,4,\ldots ,N-k\),

$$\begin{aligned} \int _{D_\rho }\frac{\partial V(r,z'')}{\partial z_j}u_m^2&= m\bigg (\frac{1}{\lambda ^2}\frac{\partial V(\bar{r},\bar{z}'')}{\partial \bar{z}_j}\int _{\mathbb {R}^N}U_{0,1}^2+o\Big (\frac{1}{\lambda ^2}\Big )\bigg ),\\ \int _{D_\rho }\frac{1}{2r}\frac{\partial (r^2 V(r,z''))}{\partial r}u_m^2&= m\bigg (\frac{1}{\lambda ^2}\frac{1}{2\bar{r}}\frac{\partial (\bar{r}^2 V(\bar{r},\bar{z}''))}{\partial \bar{r}}\int _{\mathbb {R}^N}U_{0,1}^2+o\Big (\frac{1}{\lambda ^2}\Big )\bigg ). \end{aligned}$$

Now we are ready to prove Theorem 1.5.

Proof of Theorem 1.5

Combining Lemmas 3.2, 3.3, 3.4, 3.7 and Remark 3.9, we find that there exists \(\rho \in (3\delta ,4\delta )\) such that (3.1), (3.2) and (3.3) are, respectively, equivalent to

$$\begin{aligned}&\frac{\partial (\bar{r}^2V(\bar{r},\bar{z}''))}{\partial \bar{r}}=o(1), \end{aligned}$$
(3.21)
$$\begin{aligned}&\frac{\partial (\bar{r}^2V(\bar{r},\bar{z}''))}{\partial \bar{z}_j}=o(1),\quad j=3,4,\ldots ,N-k, \end{aligned}$$
(3.22)
$$\begin{aligned}&B_1V(\bar{r},\bar{z}'')-\frac{B_3m^{N-2}}{\lambda ^{N-4}}=o(1). \end{aligned}$$
(3.23)

Set \(\lambda =tm^{\frac{1}{\tau }}\), then we get from (3.23) that

$$\begin{aligned} B_1V(\bar{r},\bar{z}'')-\frac{B_3}{t^{N-4}}=o(1),\quad t\in [T_0,T_1]. \end{aligned}$$
(3.24)

Let \(F(t,\bar{r},\bar{z}'')=\big (\nabla _{\bar{r},\bar{z}''}(\bar{r}^2V(\bar{r},\bar{z}'')),B_1V(\bar{r},\bar{z}'')-\frac{B_3}{t^{N-4}}\big )\), then

$$\begin{aligned}&\deg \big (F(t,\bar{r},\bar{z}''),[T_0,T_1]\times B_\varrho (r_0,z_0'')\big )\\&\quad =\deg \big (\nabla _{\bar{r},\bar{z}''}(\bar{r}^2V(\bar{r},\bar{z}'')),B_\varrho (r_0,z_0'')\big )\cdot \deg \Big (B_1V(\bar{r},\bar{z}'')-\frac{B_3}{t^{N-4}},[T_0,T_1]\Big )\\&\quad = \deg \big (\nabla _{\bar{r},\bar{z}''}(\bar{r}^2V(\bar{r},\bar{z}'')),B_\varrho (r_0,z_0'')\big )\ne 0. \end{aligned}$$

Hence there exist \(t_m\in [T_0,T_1]\) and \((\bar{r}_m,\bar{z}_m'')\in B_\varrho (r_0,z_0'')\) satisfying (3.21), (3.22) and (3.24). \(\square \)