1. Introduction and main results
In present paper, we are interested in the following fractional system of Hamiltonian type:
\begin{equation}
\begin{cases}
\begin{aligned}
(-\Delta)^{s} u&=H_v(u,v) \;\;&&\text{in}~\Omega,\\
(-\Delta)^{s} v&=H_u(u,v) &&\text{in}~\Omega,\\
u &= v = 0 &&\text{in} ~ \mathbb{R}^N\setminus\Omega,
\end{aligned}
\end{cases}
\end{equation} where 
$s\in (0,1)$, 
$N \gt 2s$, 
$H \in C^1(\mathbb{R}^2, \mathbb{R})$, and 
$\Omega \subset \mathbb{R}^N$ is a bounded domain. We will not focus on the regularity condition of the domain Ω, and say simply that Ω is smooth, even most results work with 
$C^{1, 1}$ boundary. The study of system (1.1) is mainly motivated by the following classical Hamiltonian system:
\begin{equation}
\begin{cases}
\begin{aligned}
-\Delta u&=H_v(u,v) \,&&\text{in}~\Omega,\\
-\Delta v&=H_u(u,v) &&\text{in}~\Omega,\\
u &= v = 0 &&\text{in} ~ \partial\Omega.\\
\end{aligned}
\end{cases}
\end{equation}Formally the energy functional of (1.2) is
\begin{equation}
\mathcal{L}(u, v)=\int_{\Omega}\nabla u \nabla v dx-\int_{\Omega}H(u, v) dx.
\end{equation} A crucial question is to decide on which space we will consider the functional 
$\mathcal{L}$. Let us look at the famous Lane–Emden system where 
$H(u,v) = \frac{|u|^{p+1}}{p+1} + \frac{|v|^{q+1}}{q+1}$, 
$p, q \gt 0$, i.e.
\begin{equation}
\begin{cases}
\begin{aligned}
-\Delta u& =|v|^{q-1}v \,&&\text{in}~\Omega,\\
-\Delta v& =|u|^{p-1}u \,&&\text{in}~\Omega,\\
u &= v = 0~~&&\text{in} ~\partial\Omega.\\
\end{aligned}
\end{cases}
\end{equation}It is well known that the existence of positive solutions to (1.4) on a ball is decided by the position of (p, q) with respect to the critical hyperbola
\begin{equation}
p, q \gt 0, \quad \frac{1}{p+1}+\frac{1}{q+1}=\frac{N-2}{N}.
\end{equation} A naive choice of the functional space for 
$\mathcal{L}$ is 
$H_0^1(\Omega)^2$. However, the Sobolev embedding greatly restricts the growth of H, and we could only work with 
$\max(p, q) \leq \frac{N+2}{N-2}$, hence many other choices of H were eliminated. Another difficulty is to deal with the strong indefiniteness of the quadratic part in 
$\mathcal{L}$, which is neither bounded from below nor from above on any subspace of 
$H_0^1(\Omega)^2$ with finite codimension.
To include more choices of H as for (1.4), people thought about functional spaces without symmetry in u and v, such as 
$W_0^{1,t}(\Omega)\times W_0^{1,\frac{t}{t-1}}(\Omega)$ with t > 1. But a new difficulty occurs since this choice is not a Hilbert space when t ≠ 2, which prevents us from applying the linking theory due to Benci and Rabinowitz [Reference Benci and Rabinowitz3].
As far as we are aware, Hulshof and van der Vorst [Reference Hulshof and van der Vorst24] first used the interpolation space framework to deal with the system with 
$H(u, v) = F(u)+G(v)$, that is
\begin{equation}
\begin{cases}
\begin{aligned}
-\Delta u & =g(v) \;\;&&\text{in}~\Omega,\\
-\Delta v& =f(u) &&\text{in}~\Omega,\\
u&= v=0 &&\text{in} ~ \partial\Omega,
\end{aligned}
\end{cases}
\end{equation} where 
$F(t)=\displaystyle \int_0^{t}f(\tau) d\tau$, 
$G(t)=\displaystyle \int_0^{t}g(\tau) d\tau$. They replaced the first integral in (1.3) by
\begin{equation*}
\int_{\Omega}\widetilde{A}^{t}u\widetilde{A}^{2-t}v dx \quad \text{with}\;\; t\in (0, 2),
\end{equation*} where 
$\widetilde{A}^t: \widetilde{H}^{t}(\Omega)\to L^2(\Omega)$ is a family of interpolation operators and 
$\widetilde{H}^{t}(\Omega)$ is a family of interpolation spaces between 
$L^2(\Omega)$ and 
$H^2(\Omega)\cap H_0^1(\Omega)$. More precisely, for any
\begin{equation}
p, q \gt 0, \quad \frac{1}{p+1}+\frac{1}{q+1} \gt \frac{N-2}{N}, \end{equation} there exist some 
$t\in (0, 2)$ such that the embedding
\begin{equation*}\widetilde{H}^t(\Omega)\times \widetilde{H}^{2-t}(\Omega)\subset L^{p+1}(\Omega)\times L^{q+1}(\Omega)\end{equation*}is compact, which yields the existence of nontrivial solution to (1.6) under following conditions
• (Subcritical)
$f(u)=O(|u|^p), g(v)=O(|v|^q)$ as
$|u|, |v|\to \infty$ with
$p, q \gt 1$ satisfying (1.7);• (Superlinear)
$f(u)=o(|u|)$,
$g(v)=o(|v|)$ as
$|u|, |v|\to 0$;• (AR condition)
$f, g\in C(\mathbb{R})$,
$f(0)=g(0)=0$, and there exist γ > 2, R > 0 such that
\begin{equation*}
0 \lt \gamma F(u)\le uf(u),\quad 0 \lt \gamma G(v)\le vg(v),\quad \forall\; |u|, |v|\ge R.
\end{equation*}
This generalized clearly the study of Lane–Emden system (1.4).
Later, de Figueiredo and Felmer also applied in [Reference de Figueiredo and Felmer19] the interpolation space method to handle (1.2), with more general coupled H where 
$p, q$ satisfy (1.7),
\begin{equation*}
0 \le H(u, v)\le C(|u|^{p+1}+|v|^{q+1}),
\end{equation*}and H is superlinear at 0, that is
\begin{equation*}
C(|u|^{p+1}+|v|^{q+1}) \le H(u, v) \;\;\text{with}\;\; pq \gt 1, |u|+|v| \lt r \mbox{ for some } r \gt 0.
\end{equation*} In addition, Ambrosetti–Rabinowitz ((AR) for short) type condition (as 
$(H2)$ below) was assumed. They proved that (1.2) admits a nontrivial solution if
\begin{equation*}
(N-4)\max\{p, q\} \lt N+4.
\end{equation*} Remark that to work with more general coupled H, we need to restrict the upper bound of 
$p, q$.
Clément and van der Vorst in [Reference Clément and van der Vorst15] proposed another approach to study (1.2), their idea is to take advantage of the dual method developed by Clarke, Ekeland, and Temam [Reference Clarke and Ekeland16, Reference Ekeland and Temam21]. Here the nonlinearity H is assumed strictly convex, subcritical (see (1.7)) and superlinear at 0. This dual method will transform (1.2) into a problem where the energy functional possesses a mountain pass structure. In fact, consider the Legendre–Fenchel transform of H (see [Reference Ekeland and Temam21], [Reference Struwe37, chapter I, section 6]):
\begin{equation*}
H^*(x)=\underset{w\in \mathbb{R}^2}{\sup}\big\{\langle w, x\rangle -H(x)\big\}\quad \mbox{for any } x \in \mathbb{R}^2,
\end{equation*}one can obtain solutions to (1.2) by critical points of
\begin{equation*}
\mathcal{J}(u, v)=\int_{\Omega}{H}^*(u, v) dx-\int_{\Omega}v\mathcal{A}u dx,\quad \forall \, (u, v)\in L^{1+\frac1p}(\Omega)\times L^{1+\frac1q}(\Omega),
\end{equation*} where 
$\mathcal{A}$ is the inverse of
\begin{align}
-\Delta : W^{2, 1+\frac1p}(\Omega)\cap W_0^{1, 1+\frac1p}(\Omega)\to L^{1+\frac1p}(\Omega).
\end{align} If moreover H satisfies (AR) type condition (similar to 
$(H5)$ and 
$(H7)$ below), the existence of nontrivial solution to (1.2) was shown in [Reference Clément and van der Vorst15]. This dual method was also used to handle system with critical growth, see Hulshof et al. [Reference Hulshof, Mitidieri and van der Vorst25].
Some other approaches exist. de Figueiredo, do Ó and Ruf [Reference de Figueiredo, do Ó and Ruf18] used Orlicz space to obtain nontrivial solutions of (1.6), they replaced 
$W_0^{1,t}(\Omega)\times W_0^{1,\frac{t}{t-1}}(\Omega)$ by Sobolev–Orlicz space 
$W_0^1L_A(\Omega)\times W_0^1L_{\widetilde A}(\Omega)$, where A is an N-function and 
${\widetilde A}$ is its Young-conjugate. Owing to the fact that this setting is not a Hilbert space, they used a finite-dimensional approximation method. Their models also contain nonlinearites with nearly critical growth. The Lyapunov–Schmidt reduction approach was also applied to problem (1.6), see for instance [Reference Ramos and Tavares32, Reference Ramos, Tavares and Zou33]. For more literature on this topic, we refer to [Reference Bonheure, Moreira dos Santos and Tavares6] and references therein.
Coming back to the special case (1.4). In [Reference Clément, de Fegueiredo and Mitedieri14], the existence of positive solution to (1.4) was firstly considered based on a topological method. When 
$p, q \gt 0$ and pq < 1, the uniqueness of positive solution was investigated in [Reference Dalmasso17]. In [Reference Bonheure, Moreira dos Santos and Ramos5], (1.4) was reduced to the following single equation
\begin{equation}
\begin{cases}
\begin{aligned}
\Delta \big(|\Delta u|^{\frac1p -1}\Delta u\big) & =|u|^{q-1}u~~&& \text{in}~\Omega,\\
u & = \Delta u=0 &&\text{in}~\partial\Omega.
\end{aligned}
\end{cases}
\end{equation} If 
$p, q \gt 0$, pq ≠ 1 and subcritical as in (1.7), then (1.4) admits a positive classical ground state solution (it has minimal energy among all solutions), see [Reference Bonheure, Moreira dos Santos and Ramos5]. If Ω is a ball, we can use the Schwartz rearrangement to show that the ground state solution is radially symmetric. Furthermore, (1.4) has no positive solutions, if Ω is star-shaped and 
$p, q \gt 0$ do not satisfy (1.7), see [Reference Mitidieri29].
To the best of our knowledge, for the fractional Laplacian case, although some special cases such as fractional Lane–Emden systems were studied (see [Reference Leite and Montenegro26] and references therein), there exists no study of (1.1) for general coupled subcritical nonlinearities so far. As mentioned above, a key step to handle (1.1) with the variational approach is to establish a suitable functional framework. Furthermore, in the Laplacian case, regardless of interpolation method, dual method, or reduced into a single equation, one needs the isomorphism given by (1.8). In the fractional Laplacian case, we need to find suitable functional space which plays the role of 
$W^{2, 1+\frac1p}(\Omega)\cap W_0^{1, 1+\frac1p}(\Omega)$.
1.1. Weak solution and fractional spaces
Let Ω be a smooth bounded domain and 
$s \in (0, 1)$, we denote
\begin{equation*}
X_0^s(\Omega):=\left\{\!u\in L^2(\mathbb{R}^N)\,{:}\,u=0~ \mbox{a.e.~in}~\mathbb{R}^N \,{\setminus}\, \Omega \;\;\text{and}\; \;\int_{\mathbb{R}^{2N}}\!{\frac{|u(x)-u(y)|^2}{|x-y|^{N+2s}}}dxdy \,{ \lt }\, \infty\!\right\}\!,
\end{equation*}endowed with the norm
\begin{equation*}
\|u\|_{X_0^s(\Omega)}:=\left( \int_{\mathbb{R}^{2N}}{\frac{|u(x)-u(y)|^2}{|x-y|^{N+2s}}}dxdy\right)^{1/2}.
\end{equation*} The embedding 
$X_0^s(\Omega)\hookrightarrow L^r(\Omega)$ is continuous for 
$r\in [1,\frac{2N}{N-2s}]$ and compact for 
$r\in [1,\frac{2N}{N-2s})$, see [Reference Di Nezza, Palatucci and Valdinoci20, Theorems 6.5, 7.1]. 
$(-\Delta)^s$ is an operator from 
$X_0^s(\Omega)$ into its dual space, namely
\begin{equation*}
\langle (-\Delta)^s u, v \rangle=\int_{\mathbb{R}^{2N}}{\frac{(u(x)-u(y))(v(x)-v(y))}{|x-y|^{N+2s}}}dxdy,\quad \forall\, u, v\in X_0^s(\Omega).
\end{equation*}Formally the energy functional associated to (1.1) is
\begin{equation}
\mathcal{K}(u, v)=\int_{\mathbb{R}^{2N}}{\frac{(u(x)-u(y))(v(x)-v(y))}{|x-y|^{N+2s}}}dxdy-\int_{\Omega}H(u, v) dx.
\end{equation} As mentioned before, the crucial question is to find a suitable functional space to work with 
$\mathcal{K}$, and 
$X_0^s(\Omega)^2$ would not be the right one if we hope to handle more general function H.
Next we recall the regularity result for the linear equation
\begin{equation}
(-\Delta)^{s} u\;=f \;\;\text{in}~\Omega,\quad
u =0~\;\;\text{in} ~ \mathbb{R}^N\setminus\Omega.
\end{equation}We denote
\begin{equation*}
\delta(x):=\text{dist}(x, \partial\Omega),\quad x\in\Omega,
\end{equation*}and
\begin{equation*}
C_{\delta}^{\alpha}(\overline{\Omega}):= \Big\{u\in C(\overline\Omega): \frac{u}{\delta^s} \mbox{admits a continuous extension belonging to } C^{\alpha}(\overline{\Omega})\Big\}.
\end{equation*}Ros-Oton and Serra in [Reference Ros-Oton and Serra35, Proposition 1.1, Theorem 1.2] proved that
Lemma 1.1. Let Ω be a bounded 
$C^{1,1}$ domain and 
$f\in L^{\infty}(\Omega)$. There exists a unique 
$u\in X_0^s(\Omega)$ solving (1.11). Moreover 
$u\in C^s(\overline{\Omega})\cap C_{\delta}^{\alpha}(\overline{\Omega})$, and
\begin{equation*}
\|u\|_{C^s(\overline{\Omega})}+ \left\|\frac{u}{\delta^s}\right\|_{C^\alpha(\overline{\Omega})}\le C\|f\|_{\infty}
\end{equation*} for some 
$0 \lt \alpha \lt \min\{s, 1-s\}$. The constants α, 
$C$ depend only on Ω and 
$s$.
Here we shall consider a weaker solution (see definition below) and choose the test function space as
\begin{equation*}
\mathcal{T}_s(\Omega):=\{u\in X_0^s(\Omega): (-\Delta)^s u\in C_c^{\infty}(\Omega)\}.
\end{equation*} Hence 
$\mathcal{T}_s(\Omega) \subset C^s(\overline{\Omega})\cap C_{\delta}^{\alpha}(\overline{\Omega})$. Let
\begin{equation*}
L^1(\Omega; \delta^s dx) :=\Big\{u\in L^1_{loc}(\Omega): \int_{\Omega} |u|\delta^sdx \lt \infty\Big\}
\end{equation*} be endowed with the norm 
$\|u\|_{L^1(\Omega; \delta^s dx)}= |u \delta^s|_1$. In this paper, 
$|\cdot|_r$ denotes always the norm of 
$L^r(\Omega)$.
Definition 1.2. Let 
$s\in (0, 1)$, 
$N \gt 2s$ and 
$f\in L^1(\Omega; \delta^s dx)$. We say that 
$u$ is a 
$L^1$-weak solution to (1.11) if 
$u\in L^1(\Omega)$, and
\begin{align}
\int_{\Omega}u(-\Delta)^s\varphi dx=\int_{\Omega}f \varphi dx,\quad \forall \, \varphi\in \mathcal{T}_s(\Omega).
\end{align} Note that similar notions were given in [Reference Leonori, Peral, Primo and Soria27, Reference Ros-Oton and Serra34]. We will prove in proposition 2.3 that for any 
$f\in L^1(\Omega; \delta^s dx)$, there exists a unique L 1-weak solution to (1.11). We also give the comparison principle and maximum principle for fractional Laplacian in L 1-weak sense, see lemmas 2.1 and 2.2. Accordingly, we define
Definition 1.3. For 
$s\in (0, 1)$ and 
$N \gt 2s$, 
$(u, v)$ is said a 
$L^1$-weak solution to system (1.1) if 
$u, v \in L^1(\Omega)$, 
$H_u(u, v), H_v(u, v)\in L^1(\Omega; \delta^sdx)$ and for any 
$(\varphi, \psi)\in \mathcal{T}_s(\Omega)\times \mathcal{T}_s(\Omega)$,
\begin{equation*}
\int_{\Omega} u(-\Delta)^s\varphi dx=\int_{\Omega} H_v(u, v)\varphi dx,\quad
\int_{\Omega} v(-\Delta)^s\psi dx=\int_{\Omega} H_u(u, v)\psi dx.
\end{equation*} If in addition 
$u, v\in L^{\infty}(\Omega)$, we call 
$(u, v)$ a classical solution.
Another important functional space for us is the set of u such that 
$(-\Delta)^s u\in L^r(\Omega)$, namely for 
$r \geq 1$,
\begin{equation}
\mathcal{W}^{2s, r}(\Omega):=\{u\in L^1(\Omega): \exists\; f \in L^r(\Omega) \mbox{ such that (1.12) holds true}\}
\end{equation}endowed with the norm
\begin{equation*}\|u\|_{\mathcal{W}^{2s, r}(\Omega)}=|(-\Delta)^s u|_r.\end{equation*} We shall use 
$\mathcal{W}^{2s, r}(\Omega)$ to play the role of 
$W^{2, r}(\Omega)\cap W_0^{1, r}(\Omega)$ in the Laplacian case. In particular, we denote
\begin{equation*}
\mathcal{H}^{2s}(\Omega)=\mathcal{W}^{2s, 2}(\Omega).
\end{equation*} Some (compact) embedding properties of 
$\mathcal{W}^{2s, r}(\Omega)$ will be shown in proposition 2.6.
1.2. Fractional Hamiltonian system
Motivated by [Reference de Figueiredo and Felmer19, Reference Hulshof and van der Vorst24], we apply firstly the interpolation space method to study the system (1.1). For 
$0\le\alpha\le 2s$, consider the interpolation space
\begin{equation}
E^\alpha:=\Big\{u=\sum_{j\ge 1}a_j\varphi_j\in L^2(\Omega): \sum_{j\ge 1}\lambda_j^{\frac{\alpha}{s}}a_j^2 \lt \infty \Big\},
\end{equation} where λj is the j-th eigenvalue of 
$(-\Delta)^s$ with corresponding eigenfunction φj, and 
$\{\varphi_j\}$ forms an orthonormal basis of 
$L^2(\Omega)$. For 
$\alpha \in (0, 2s)$, let 
$A^{\alpha}: E^{\alpha}\to L^{2}(\Omega)$ be given by
\begin{equation}
A^{\alpha}u:=\sum_{j\ge 1}\lambda_j^{\frac{\alpha}{2s}}a_j\varphi_j \in L^2(\Omega),\quad \forall \;u=\sum_{j\ge 1}a_j\varphi_j\in E^\alpha.
\end{equation}E α is clearly a Hilbert space with the scalar product
\begin{equation*}
(u,v)_{E^\alpha}= \int_{\Omega}A^{\alpha}uA^{\alpha}v dx = \sum_{j\ge 1}\lambda_j^{\frac{\alpha}{s}}\langle u, \varphi_j\rangle_{L^2}\langle v, \varphi_j\rangle_{L^2}.
\end{equation*}To handle the Hamiltonian system (1.1), we define
\begin{equation*}
\mathbf{E}_{\alpha}:=E^{\alpha}\times E^{2s-\alpha}, \quad \alpha \in (0, 2s).
\end{equation*} Instead of considering 
$\mathcal K$ over 
$X_0^s(\Omega)^2$, we consider the energy functional
\begin{equation}
\mathcal{E}(u,v)=\int_{\Omega}A^{\alpha}uA^{2s-\alpha}v dx-\int_{\Omega}H(u, v) dx,\quad \forall\; (u,v)\in \mathbf{E}_{\alpha}.
\end{equation} The choice of 
$\mathbf{E}_{\alpha}$ originates from three observations.
•
$\mathcal{H}^{2s}(\Omega) = E^{2s}$, seeing proposition 3.1;• E α can be embedded compactly into some
$L^r(\Omega)$, seeing proposition 2.6 and remark 3.2;• Every critical point of
$\mathcal{E}$ is a L 1-weak solution to (1.1), seeing proposition 3.4.
Applying a linking theorem in [Reference Felmer22], we show the following existence result, which extends the study for the Laplacian case in [Reference de Figueiredo and Felmer19].
Theorem 1.4. Let 
$N \gt 2s$ and 
$p, q \gt 0$ satisfy
\begin{equation}
1 \gt \frac1{p+1}+\frac1{q+1} \gt \frac{N-2s}{N},
\end{equation}and
\begin{equation}
(N-4s)\max\{p, q\} \lt N+4s.
\end{equation} Assume that 
$H \in C^1(\mathbb{R}^2, \mathbb{R})$ satisfies
(H 1)
$H \ge 0$ in
$\mathbb{R}^2$;(H 2) There exists
$R \gt 0$ such that
(1.19)\begin{equation}
\frac{1}{p+1}H_u(u,v)u+\frac{1}{q+1}H_v(u,v)v\ge H(u,v) \gt 0, \quad \forall\; |u|+|v|\ge R;
\end{equation}(H 3) There exist
$r \gt 0$,
$C \gt 0$ such that
(1.20)\begin{equation}
|H(u,v)|\le C(|u|^{p+1}+|v|^{q+1}), \quad \forall \; |u|+|v|\le r;
\end{equation}(H 4) There exists
$C \gt 0$ such that for any
$(u, v) \in \mathbb R^2$,
(1.21)\begin{equation}
|H_u(u,v)|\le C\Big(|u|^p+|v|^\frac{p(q+1)}{p+1}+1\Big),\quad |H_v(u,v)|\le C\Big(|v|^q+|u|^\frac{q(p+1)}{q+1}+1\Big).
\end{equation}
Then, there exists a nontrivial classical solution to (1.1).
Note that theorem 1.4 requires the assumption (1.18). For getting the existence of solutions to (1.1) in a more broad range of p, q, we will also apply the dual method. Three major difficulties subsist.
• The first problem is still to construct a suitable functional framework. We choose (see §4) the energy functional as follows
\begin{equation*}
\mathcal{J}(f, g)=\mathcal{H}^*(f, g)-\int_{\Omega}g\mathcal{A}f dx,\quad \forall \, (f, g)\in X:=L^{1+\frac1p}(\Omega)\times L^{1+\frac1q}(\Omega)
\end{equation*}where
$\mathcal{A}$ is the inverse of
$(-\Delta)^s: \mathcal{W}^{2s, 1+\frac1p}(\Omega)\to L^{1+\frac1p}(\Omega)$ and
$\mathcal{H}^*$ denotes the Legendre–Fenchel transform of
\begin{equation*}
\mathcal{H}(u, v)=\int_{\Omega}H(u, v) dx, \quad \forall\; (u, v)\in X^*=L^{p+1}(\Omega)\times L^{q+1}(\Omega).
\end{equation*}• Secondly, we need to verify the differentiability of
$\mathcal{J}$, and the correspondence between critical points of
$\mathcal{J}$ and solutions to (1.1). It should be mentioned that these arguments were not proven explicitly in [Reference Clément and van der Vorst15] for the Laplacian case.• Finally, we need to check the compact embedding properties of
$\mathcal{W}^{2s, 1+\frac1p}(\Omega)$ in order to check the Palais–Smale condition.• Finally, we need to check the compact embedding properties of
$\mathcal{W}^{2s, 1+\frac1p}(\Omega)$ in order to check the Palais–Smale condition.
To deal with these difficulties, we use properties of Legendre–Fenchel transform (see lemma 2.9) to ensure the well-definedness of 
$\mathcal{J}$. For the differentiability of 
$\mathcal{J}$, we will prove in lemma 4.5 that
\begin{equation*}
\mathcal{H}^*(f, g)=\int_{\Omega}H^*(f, g) dx,\quad \forall\; (f, g)\in X,
\end{equation*} where 
$H^*$ is the Legendre–Fenchel transform of 
$H: \mathbb{R}^2\to \mathbb{R}$. We use the fact 
$(\nabla H)^{-1}=\nabla H^*$ (see lemma 4.3) which guarantees that 
$\mathcal{H}^*$ is well defined and of class C 1 over X. Moreover, under the superlinear growth assumption of H, 
$\mathcal{J}$ has a mountain pass geometry, and 
$H_z^*(f, g) := \nabla H^*(f, g)$ provides a weak solution to (1.1) given any critical point 
$(f, g)\in X$ for 
$\mathcal{J}$, see proposition 4.7. The compact embeddings of 
$\mathcal{W}^{2s, 1+\frac1p}(\Omega)$ are given in proposition 2.6.
Theorem 1.5. Let 
$N \gt 2s$, 
$p, q \gt 0$ satisfy (1.17). Assume that 
$H \in C^1(\mathbb{R}^2, \mathbb{R})$ satisfies
(H 5) There exist positive numbers
$C_1, C_2$ such that
\begin{equation*}
\begin{split}
&C_1|u|^{p+1}\le H_u(u,v)u\le C_2\left( |u|^{p+1}+ |u|^{\alpha}|v|^{\beta}\right),\\
&C_1|v|^{q+1}\le H_v(u,v)v\le C_2\left( |v|^{q+1}+ |u|^{\alpha}|v|^{\beta}\right),
\end{split}
\end{equation*}with
(1.22)\begin{equation}
\frac{\alpha}{p+1}+\frac{\beta}{q+1}=1,\;\; \alpha, \beta \gt 1;
\end{equation}(H 6)
$\nabla H$ is strictly monotone, i.e.
\begin{align*}
\langle \nabla H(u_1, v_1)- & \nabla H(u_2, v_2),\, (u_1,v_1)- (u_2,v_2)\rangle \gt 0,\quad \mbox{for any disjoint } \\ & (u_1, v_1),\, (u_2, v_2)\in \mathbb{R}^2;
\end{align*}(H 7) There are
$\theta \in (0, 1)$ and positive numbers
$C_3$,
$C_4$ such that
\begin{equation*}
\theta H_u(u,v)u+ (1- \theta)H_v(u,v)v- H(u,v)\ge C_3\Big(|u|^{p+1}+|v|^{q+1}\Big)-C_4.
\end{equation*}
Then there exists a nontrivial classical solution to (1.1).
Remark 1.6. 
$(H5)$ implies indeed 
$(H1)$ and 
$(H4)$. 
$(H7)$ is useful to prove that every Palais–Smale sequence of 
$\mathcal{J}$ is bounded, see lemma 4.8. Without loss of generality, we can assume 
$H(0, 0)=0$. Otherwise, we replace 
$H(u, v)$ by 
$H(u, v)-H(0, 0)$.
Remark 1.7. The nonlinearities H in theorems 1.4 and 1.5 both have subcritical and superlinear growth. But some differences exist between the two families of assumptions. For example, let H have the form
\begin{equation}
H_\varepsilon(u, v)=|u|^{p+1}+|v|^{q+1}+\varepsilon |u|^{\alpha}|v|^{\beta},
\end{equation} where 
$p, q$ satisfy (1.17), and 
$\alpha, \beta$ satisfy (1.22). When (1.18) holds true, 
$H_{\varepsilon}$ satisfies 
$(H1)$–
$(H4)$ for all 
$\varepsilon \gt 0$. On the other hand, we need not (1.18) in theorem 1.5, but the strictly convexity assumption 
$(H6)$ fails for 
$H_\varepsilon$ when 
$\varepsilon$ is sufficiently large. In other words, theorem 1.4 holds for more broad coupling nonlinearities but requires narrow choices of 
$p, q$; theorem 1.5 can work for all subcritical and superlinear 
$p, q$, meanwhile the strict convexity is more restrictive for the coupling term.
1.3. Fractional Lane–Emden system
As a special example, we revisit the fractional Lane–Emden system
\begin{equation}
\begin{cases}
\begin{aligned}
(-\Delta)^{s} u&=|v|^{q-1}v \,&&\text{in}~\Omega,\\
(-\Delta)^{s} v&=|u|^{p-1}u \,&&\text{in}~\Omega,\\
u &= v = 0 &&\text{in} ~ \mathbb{R}^N\setminus\Omega,\\
\end{aligned}
\end{cases}
\end{equation} where 
$s\in (0, 1)$, 
$p, q\in (0, \infty)$, 
$N \gt 2s$. As for the classical Laplacian case, we consider subcritical exponents 
$p, q$, that is
\begin{equation}
p, q \gt 0,\quad\frac{1}{p+1}+\frac{1}{q+1} \gt \frac{N-2s}{N}.
\end{equation} Leite and Montenegro [Reference Leite and Montenegro26] showed the existence of positive viscosity solutions to (1.24) under the subcritical condition (1.25). They reduced (1.24) to a single equation and considered energy functional on 
$W_0^{1, 1+\frac1q}(\Omega)\cap W^{2s, 1+\frac1q}(\Omega)$. Choi and Kim [Reference Choi12, Reference Choi and Kim13] studied the related problems with respect to spectral fractional Laplacian.
Different from [Reference Leite and Montenegro26], we work with the functional space 
$\mathcal{W}^{2s, 1+\frac1q}(\Omega)$. Comparing with [Reference Leite and Montenegro26], our setting avoid many regularity problems since every 
$\mathcal{W}^{2s, 1+\frac1q}(\Omega)$ solution is naturally a L 1-weak solution. More precisely, we consider
\begin{equation}
\begin{cases}
\begin{aligned}
(-\Delta)^{s}\left(|(-\Delta)^{s}u|^{\frac1q-1}(-\Delta)^{s}u\right) &=|u|^{p-1}u \;&&\text{in}~\Omega,\\
u &=(-\Delta)^{s}u=0 &&\text{in} ~ \mathbb{R}^N\setminus\Omega.\\
\end{aligned}
\end{cases}
\end{equation}Definition 1.8. We call u an energy solution of (1.26), if 
$u\in \mathcal{W}^{2s, 1+\frac1q}(\Omega)$ and
\begin{equation}
\int_{\Omega}|(-\Delta)^{s}u|^{\frac1q-1}(-\Delta)^{s}u(-\Delta)^{s}\varphi dx=\int_{\Omega}|u|^{p-1}u\varphi dx,\quad \forall\,\varphi\in \mathcal{W}^{2s, 1+\frac1q}(\Omega).
\end{equation} Note that the above terms are well defined since 
$(-\Delta)^{s}u\in L^{1+\frac{1}{q}}(\Omega)$ for any 
$u\in \mathcal{W}^{2s, 1+\frac1q}(\Omega)$. Putting 
$v=|(-\Delta)^{s}u|^{\frac1q-1}(-\Delta)^{s}u$, we can prove that (u, v) is a classical solution to (1.24) if and only if u is an energy solution to (1.26), see proposition 5.1.
Moreover, solutions to (1.26) coincide with critical points of the following C 1 functional in 
$\mathcal{W}^{2s, 1+\frac1q}(\Omega)$,
\begin{equation*}
\mathcal{I}(u) :=\frac{q}{q+1}\left|(-\Delta)^su\right|_{1+\frac1q}^{1+\frac1q}-\frac1{p+1}|u|_{p+1}^{p+1}.
\end{equation*} Consider the Nehari manifold associated to 
$\mathcal{I}$,
\begin{equation}
\mathcal{N}_{\mathcal{I}}:=\{u\in \mathcal{W}^{2s, 1+\frac1q}(\Omega)\backslash \{0\}: \langle \mathcal{I}'(u), u\rangle=0\},
\end{equation}and the ground state level is defined as
\begin{equation}
c_{\mathcal{I}}:=\underset{u\in \mathcal{N}_{\mathcal{I}}}{\inf}\mathcal{I}(u).
\end{equation} We establish the existence of positive solutions to (1.24) by showing that 
$c_{\mathcal{I}}$ can be attained.
Theorem 1.9. Assume that 
$s\in (0, 1)$, 
$N \gt 2s$, 
$p, q$ satisfy (1.25) and 
$pq\neq 1$. Then
When Ω is a ball and pq < 1, the uniqueness of positive classical solution to (1.24) ensures that the solution is radially symmetric. For pq > 1, we can not claim directly as in [Reference Bonheure, Moreira dos Santos and Ramos5] that the ground state solution is radially symmetric, but we can work simply with 
$\mathcal{W}_{{\mathrm{rad}}}^{2s, r}(B_R)$, the subset of radial functions in 
$\mathcal{W}^{2s, r}(B_R)$. In fact, for 
$f\in C(\overline{B_R})$, 
$f\ge 0$, if u, w are respectively solutions of
\begin{equation*}
\begin{cases}
\begin{aligned}
(-\Delta)^{s} u&=f &&\text{in}~ B_R,\\
u&=0~&&\text{in} ~ \mathbb{R}^N\setminus\overline{B_R},\\
\end{aligned}
\end{cases}
\quad\quad
\begin{cases}
\begin{aligned}
(-\Delta)^{s} w&=f^{\#} &&\text{in}~B_R,\\
w&=0~&&\text{in} ~ \mathbb{R}^N\setminus\overline{B_R},\\
\end{aligned}
\end{cases}
\end{equation*} where 
$f^{\#}$ denotes the radial decreasing rearrangement of f, we do not have always 
$u^{\#}\le w$, see [Reference Ferone and Volzone23], while it is true for s = 1.
Remark 1.10. For any classical solution (u, v) of (1.24), we have
\begin{equation*}
\int_{\Omega}|v|^{q+1} dx=\langle (-\Delta)^s u, v \rangle=\langle (-\Delta)^s v, u \rangle=\int_{\Omega}|u|^{p+1} dx.
\end{equation*} Combining with 
$\langle \mathcal{I}'(u), u\rangle=0$, there holds 
$\mathcal{K}(u, v)=\mathcal{I}(u)$. Hence 
$\mathcal{K}(u, v)$ reaches the minimal energy among all classical solutions if u is a ground state of (1.26).
The paper is organized as follows. In §2, we introduce some preliminary results. In §3, we construct an interpolation space setting and prove theorem 1.4. In §4, we construct the framework of dual method to show theorem 1.5. In §5, we revisit the fractional Lane–Emden system (1.24).
2. Notations and preliminary
Throughout this paper, we use the following notations.
• We denote by
$\nabla H$ the gradient of H in
$\mathbb R^2$, and we write
$H_z(u, v) := \nabla H(u(x), v(x))$ for the functional variable case where
$z=(u, v)$. The same convention is used for the Legendre–Fenchel transform
$H^*$.•
$C, C', C_1, C_2,...$ denote always generic positive constants.• For
$r\in [1, \infty]$, we denote by
$|u|_r$ the usual
$L^r(\Omega)$ or
$L^r(\mathbb R^N)$ norm.• For any
$p, q \gt 0$, we set
$X:=L^{1+\frac1p}(\Omega)\times L^{1+\frac1q}(\Omega)$ with norm
$\|(f, g)\|_{X}=|f|_{1+\frac1p}+|g|_{1+\frac1q}$, hence its dual space is
$X^*=L^{p+1}(\Omega)\times L^{q+1}(\Omega)$.• The Hamiltonian functional is denoted by
\begin{equation*}
\mathcal{H}(u, v) =\int_{\Omega}H(u, v) dx,\quad (u, v)\in X^*.
\end{equation*}
2.1. Basic properties for L 1-weak solutions
Here we show some elementary facts for L 1-weak solutions. Many of them were established for solutions in 
$X_0^s(\Omega)$, however it’s worthy to check carefully for L 1-weak solutions.
Let us begin with the comparison principle. We say that 
$(-\Delta)^su \le (-\Delta)^sv$ in the L 1-weak sense, if 
$u, v\in L^1(\Omega)$ and
\begin{equation}
\int_{\Omega}(u-v)(-\Delta)^s\varphi dx\le 0,\quad \forall\,\varphi\in \mathcal{T}_s(\Omega) \mbox{with } (-\Delta)^s\varphi\ge 0.
\end{equation} By the definition of 
$\mathcal{T}_s(\Omega)$, as 
$(-\Delta)^s\varphi$ fulfilled 
$C_c^\infty(\Omega)$, we get immediately
Lemma 2.1. If 
$(-\Delta)^su \le (-\Delta)^sv$ in the L 1-weak sense, then 
$u\le v$ a.e. in Ω.
The next one is a Hopf lemma result for L 1-weak supersolutions.
Lemma 2.2. Assume that 
$f\in L^1(\Omega; \delta^s dx)$ is nonnegative and 
$f\ne 0$. Then there exists C > 0 depending only on f, Ω, N and s such that for any 
$u\in L^1(\Omega)$ satisfying 
$(-\Delta)^s u\ge f$, we have 
$u \gt C\delta^s$ a.e. in Ω.
Proof. Fix 
$f_k:=\min\{k, f\}$ with k large enough such that 
$f_k\ne 0$. Let 
$u_k \subset X_0^s(\Omega)$ solve (ensured by lemma 1.1)
\begin{equation}
(-\Delta)^s u_k = f_k \;\;\text{in}~\Omega,\quad u_k =0 \;\;\text{in}~\mathbb{R}^N\setminus\Omega.
\end{equation} As 
$f_k \in L^\infty(\Omega)$, by maximum principle (see [Reference Silvestre36, Proposition 2.2.8], [Reference Brasco and Franzina8, Theorem A.1]) and Hopf’s lemma (see [Reference Caffarelli, Roquejoffre and Sire11, Proposition 2.7], [Reference Ros-Oton and Serra35, Lemma 3.2] and [Reference Pezzo and Quaas31, Theorem 1.5]) for bounded source, there holds 
$u_k \gt C\delta^s$ for some C > 0. Using the previous lemma, we get 
$u\ge u_k$, and finish the proof.
Inspired by [Reference Brezis, Cazenave, Martel and Ramiandrisoa10], we show the following existence and uniqueness result for L 1-weak solution.
Proposition 2.3. Given any 
$f\in L^1(\Omega; \delta^s dx) $, there exists a unique 
$L^1$-weak solution 
$u$ to (1.11). Moreover, there exists 
$C=C(N, s, \Omega) \gt 0$ such that
\begin{equation}
|u|_1\le C\|f\|_{L^1(\Omega; \delta^s dx) }.
\end{equation}Proof. Assume first 
$f\ge 0$, otherwise we decompose 
$f = f_+ - f_-$. Let 
$\xi\in X_0^s(\Omega)$ be the solution of
\begin{equation}
(-\Delta)^s \xi = 1 \;\;\text{in}~\Omega,\quad
\xi =0 \;\;\text{in}~\mathbb{R}^N\setminus\Omega.
\end{equation} For any 
$k\in {\mathbb N}$, set 
$f_k:=\min\{k, f\}$, and 
$u_k\geq 0$ the solution of (2.2). As 
$u_k \in X_0^s(\Omega)$, we can use it as test function to (2.4). Applying the estimates in lemma 1.1 for ξ,
\begin{equation*}
|u_k|_1 = \langle (-\Delta)^s \xi, u_k\rangle = \langle (-\Delta)^s u_k, \xi\rangle = \int_\Omega f_k\xi dx \leq C\int_\Omega f_k\delta^sdx \leq C\|f\|_{L^1(\Omega; \delta^s dx)}.
\end{equation*} Similarly, considering 
$u_k - u_l$, we have 
$|u_k-u_l|_{L^1} \leq C\|f_k-f_l\|_{L^1(\Omega; \delta^s dx)}$, so 
$\{u_k\}$ is a Cauchy sequence, hence a convergent sequence in 
$L^1(\Omega)$. We get easily a L 1-weak solution to (1.11) by taking u, the limit of uk in 
$L^1(\Omega)$. The uniqueness is ensured by the comparison principle.
Remark 2.4. If 
$\Omega = B_R$ is a ball and 
$f\in L^1(B_R, \delta^s dx)$ is radial, then there exists a unique radially symmetric solution to (1.11).
Next, we state the regularity results for the unique L 1-weak solution when f admits further integrability. Here we summarize the results in [Reference Leonori, Peral, Primo and Soria27] and [Reference Biccari, Warma and Zuazua4, Lemma 2.5].
Proposition 2.5. Assume that 
$s\in (0, 1)$, 
$N \gt 2s$, Ω is a smooth bounded domain. For any 
$f\in L^r(\Omega)$ with 
$r\ge 1$, the unique 
$L^1$-weak solution 
$u$ to (1.11) satisfies
(i) If
$1 \le r \lt \frac{N}{2s}$,
$u\in L^{\gamma}(\Omega)$ and
$|u|_{\gamma}\le C(\Omega, r, s)|f|_r$, where
\begin{equation*}
\gamma={\frac{Nr}{N-2rs}}\; {\mathrm{if} }\; r \gt 1, \quad 1\le \gamma \lt \frac{N}{N-2s}\;{\mathrm{if} }\; r=1;
\end{equation*}(ii) If
$r \ge \frac{2N}{N+2s}$, we have
$u\in X_0^s(\Omega)$ and
$\|u\|_{X_0^s(\Omega)}\le C(\Omega, r, s)|f|_r$;(iii) If
$r \gt \frac{N}{2s}$, there holds
$u\in L^{\infty}(\Omega)$ and
$|u|_{\infty}\le C(\Omega, r, s)|f|_{r}$;(iv) If
$r=\frac{N}{2s}$, there is a constant
$\alpha(\Omega, f, s) \gt 0$ such that
\begin{equation*}
\int_{\Omega} e^{\alpha |u|} dx \lt \infty.
\end{equation*}In particular,
$u\in L^\gamma(\Omega)$ for all
$1\le \gamma \lt \infty$.
2.2. Embedding and density
Here we expose some basic properties of 
$\mathcal{W}^{2s, r}(\Omega)$, including (compact) embedding results and density property, which will be important for constructing the variational framework.
Proposition 2.6. Assume that 
$s\in (0,1)$, 
$N \gt 2s$, 
$r\ge 1$, Ω is a smooth bounded domain. Then 
$\mathcal{W}^{2s, r}(\Omega)$ is a Banach space and 
$\mathcal{T}_s(\Omega)$ is dense in 
$\mathcal{W}^{2s, r}(\Omega)$. On the other hand, we have continuous embedding 
$\mathcal{W}^{2s, r}(\Omega)\subset L^{\gamma}(\Omega)$ for
\begin{equation}
\begin{cases}
\begin{aligned}
&1\le \gamma \lt \frac{N}{N-2s}, \;\; &&\text{if}~ r=1;\\
&1\le \gamma\le\frac{Nr}{N-2sr}, \;\; &&\text{if}~1 \lt r \lt \frac{N}{2s};\\
&1\le \gamma \lt \infty, && \text{if}~2sr = N;\\
& \gamma = \infty, && \text{if}~ 2sr \gt N.\\
\end{aligned}
\end{cases}
\end{equation}Moreover, the above embeddings are compact provided
\begin{equation}
1\le \gamma \lt \frac{Nr}{N-2sr}\;\; \text{if}~N \gt 2sr;\quad \mbox{or} \;\; 1\le \gamma \lt \infty, \;\; \text{if}~N\le 2sr.
\end{equation}Proof. The density of 
$\mathcal{T}_s(\Omega)$ in 
$\mathcal{W}^{2s, r}(\Omega)$ is obvious by definition and lemma 1.1.
The completeness of 
$\mathcal{W}^{2s, r}(\Omega)$ is a simple fact. Let 
$\{u_j\}\subset \mathcal{W}^{2s, r}(\Omega)$ be a Cauchy sequence, by definition, 
$\{(-\Delta)^s u_j\}$ is a Cauchy sequence in 
$L^r({\Omega})$, hence converges to some 
$v\in L^r({\Omega})$. Using proposition 2.3, there is a unique L 1-weak solution u satisfying (1.11) with f = v. Clearly 
$u \in \mathcal{W}^{2s, r}(\Omega)$, and uj converges to u in 
$\mathcal{W}^{2s, r}(\Omega)$.
On the other hand, in virtue of proposition 2.5, the embedding of 
$\mathcal{W}^{2s, r}(\Omega)\subset L^\gamma(\Omega)$ is continuous if (2.5) holds true, we shall only prove that the embedding is compact provided (2.6).
Without loss of generality, we only consider the case 
$N \gt 2sr$, 
$1 \lt r$, 
$1 \le\gamma \lt \frac{Nr}{N-2sr}$. If 
$r\ge \frac{2N}{N+2s}$, then 
$\mathcal{W}^{2s, r}(\Omega)\subset X_0^s(\Omega)$ is continuous by proposition 2.5. By Hölder inequality, fix any 
$1 \lt \gamma \lt \frac{Nr}{N-2sr}$,
\begin{equation}
|u|_{{\gamma}}\le |u|^{1-\theta}_{{1}}|u|^{\theta}_{{\frac{Nr}{N-2sr}}} \quad \mbox{with } \theta=\frac{Nr-\gamma (N-2sr)}{Nr-(N-2sr)}.
\end{equation} Since 
$\mathcal{W}^{2s, r}(\Omega)$ is continuously embedded in 
$L^{\frac{Nr}{N-2sr}}(\Omega)$ and compactly embedded in 
$L^1(\Omega)$, 
$\mathcal{W}^{2s, r}(\Omega)$ is compactly embedded in 
$L^{\gamma}(\Omega)$.
Suppose now 
$1 \lt r \lt \frac{2N}{N+2s}$, we first claim that the embedding 
$\mathcal{W}^{2s, r}(\Omega)\subset L^{1}(\Omega)$ is compact. Indeed, let 
$\{u_j\}$ be a bounded sequence in 
$\mathcal{W}^{2s, r}(\Omega)$, then 
$\{u_j\}$ is bounded in 
$L^{\frac{Nr}{N-2sr}}(\Omega)$ and the sequence 
$f_j = (-\Delta)^{s} u_j$ is bounded in 
$L^r(\Omega)$. Let A > 0 be chosen later, we denote 
$g_j = f_j\chi_{|f_j(x)|\le A}$ and 
$h_j = f_j\chi_{|f_j(x)| \gt A}$. By Hölder and Chebychev inequalities,
\begin{align*}
\|h_j\|_{L^1(\Omega; \delta^s dx)} \leq C|h_j|_1\leq C|f_j|_r\Big(\int_{|f_j(x)| \gt A} dx\Big)^{1-\frac{1}{r}} \le C'A^{\frac{1}{r}-1}.
\end{align*} Let 
$\{v_j\}$ and 
$\{w_j\}$ be respectively L 1-weak solutions of
\begin{equation}
\left\{
\begin{split}
(-\Delta)^{s} v_j &=g_j &&\text{in}~\Omega,\\
v_j&=0 &&\text{in} ~ \mathbb{R}^N\setminus\Omega,\\
\end{split}
\right. \quad \left\{
\begin{split}
(-\Delta)^{s} w_j&=h_j &&\text{in}~\Omega,\\
w_j&=0 &&\text{in} ~ \mathbb{R}^N\setminus\Omega.\\
\end{split}
\right.
\end{equation} Clearly 
$u_j=v_j+w_j$. Given any ɛ > 0, we fix A large enough such that
\begin{equation*}
|w_j|_1\le C\|h_j\|_{L^1(\Omega; \delta^s dx)} \le C'A^{\frac{1}{r}-1} \le \varepsilon.
\end{equation*} Moreover, 
$\{v_j\}$ is bounded in 
$X_0^s(\Omega)$ hence relatively compact in 
$L^1(\Omega)$. Therefore 
$\{u_j\}$ is precompact, or equivalently relatively compact in 
$L^1(\Omega)$, which means 
$\mathcal{W}^{2s, r}(\Omega)$ is compactly embedded in 
$L^{1}(\Omega)$. Applying again the interpolation inequality (2.7), 
$\{u_j\}$ is compact in 
$L^\gamma(\Omega)$ for 
$1\le\gamma \lt \frac{Nr}{N-2sr}$. So we are done.
2.3. A linking theorem
The Linking theorem 2.8 below was given by Felmer [Reference Felmer22, Theorem 3.1], it is useful for the proof of theorem 1.4. For the sake of completeness, we recall also the definition of Palais–Smale condition.
Definition 2.7. For a Banach space 
$X$, 
${\mathcal I}\in C^1(X, \mathbb{R})$ is said satisfying the Palais–Smale condition at level 
$c\in \mathbb{R}$ (for short 
$(PS)_c$) if any sequence 
$\{u_j\}\subset X$ satisfying
\begin{equation*}
{\mathcal I}(u_j)\rightarrow c \quad \mbox{and}\quad {\mathcal I}'(u_j) \rightarrow 0 \; \mbox{in }\; X^*
\end{equation*} admits a convergent subsequence in 
$X$. We say that 
${\mathcal I}$ satisfies the Palais–Smale condition ((PS) for short) if 
$(PS)_c$ is satisfied for all 
$c\in \mathbb R$.
Theorem 2.8. Let 
$(H, \langle \cdot, \cdot \rangle)$ be a Hilbert space such that 
$H=H_1\oplus H_2$. Suppose that 
$\mathcal{I}\in C^1(H)$ satisfies (PS) condition and 
$\mathcal{I}(z)=\frac12 \langle Lz, z\rangle-\mathcal{J}(z)$, where
(i)
$L : H\to H$ is a bounded, self-adjoint linear operator and
$L(H_1) \subset H_1$;
$\mathcal{J}' : H\to H^*$ is compact;(ii) There exist two linear, bounded, invertible operators
$B_1$,
$B_2 : H\to H$ such that
$\widehat{B}_{\tau}=P_2B_1^{-1}e^{\tau L}B_2: H_2\to H_2$ is invertible for any
$\tau\ge 0$. Here
$P_2: H\to H_2$ is the projection along
$H_1$ and
$e^{\tau L}=\sum_{n\in {\mathbb N}}\frac{(\tau L)^n}{n!}$;(iii) Let
$e_1\in H_1$ with
$\|e_1\|=1$. Let ρ > 0,
$R_1 \gt \frac{\rho}{\|B_1^{-1}B_2 e_1\|}$,
$R_2 \gt \rho$ and define
\begin{align*}
& S =\{B_1z_1 : z_1\in H_1, \|z_1\|=\rho\}, \\ &
Q =\{B_2(te_1+z_2) : 0\le t\le R_1, z_2\in H_2, \|z_2\|\le R_2\}.
\end{align*}Suppose that
$\mathcal{I}(z)\ge \sigma \gt 0$ on S and
$\mathcal{I}(z)\le 0$ on
$\partial Q$.
Then 
$\mathcal{I}$ has a critical point 
$z_0 \in H$ such that 
$\mathcal{I}(z_0)\ge \sigma$.
2.4. Legendre–Fenchel transform
Let V be a Banach space and 
$V^*$ be its dual space. For a function 
$G: V\to \mathbb{R}\cup \{+\infty\}$, 
$G\not\equiv +\infty$, the function 
$G^*: V^*\to \mathbb{R}\cup \{+\infty\}$ given by
\begin{equation}
G^*(u^*)=\sup\{\langle u^*, v \rangle-G(v): v\in V\},\quad \forall\, u^*\in V^*,
\end{equation} is named the Legendre–Fenchel transform of 
$G$. The following are some basic properties of 
$G^*$.
Lemma 2.9. Let 
$G^*$ be the Legendre–Fenchel transform of 
${\textit{G}}$.
(a) The Legendre–Fenchel transform reverses the order, that is,
$\widetilde{G}\ge G$ implies
$\widetilde{G}^*\le G^*$.(b) If
$G\in C^1(V, \mathbb{R})$ is convex and
$G^*\in C^1(V^*, \mathbb{R})$, then for
$v\in V$,
$u^*\in V^*$, there holds
\begin{equation*}
G^*(u^*)+G(v)=\langle u^*, v\rangle \iff v=(G^*)'(u^*) \iff u^*=G'(v).
\end{equation*}(c) If
$V=\mathbb{R}^N$,
$G$ is lower semi-continuous, strictly convex, and
$\lim_{|x|\to \infty}\frac{G(x)}{|x|}=+\infty$, then
$G^*\in C^1(\mathbb{R}^N, \mathbb{R})$.(d) If
$V=\mathbb{R}^N$,
$G(x)=\frac1{p+1}|x|^{p+1}$ with p > 0, then
$G^*(x)=\frac{p}{p+1} |x|^{1+\frac1p}$. If
$V=L^{p+1}(\Omega)$ with p > 0,
$G(v)=\frac{1}{p+1}|v|_{p+1}^{p+1}$, then
$G^*(u^*)=\frac{p}{p+1} |u^*|^{1+\frac1p}_{1+\frac1p}$.(e) Let
$a \gt 0$ and
$G_a(u):=aG(u)$, then
$G^*_a(u^*)=aG^*\big(\frac{u^*}{a}\big)$.
Proof. For 
${\rm (a)}$–
${\rm (d)}$, see respectively [Reference Struwe37, Chapter I, Section 6.2], [Reference Struwe37, Chapter I, Lemma 6.3], [Reference Mawhin and Willem28, Proposition 2.4] and [Reference Ekeland and Temam21, Proposition 4.2, Remark 4.1]. For 
${\rm (e)}$, one can check directly by means of (2.9).
3. Proof of theorem 1.4
In this section, we will use the interpolation space method to handle theorem 1.4. Throughout this section, we assume that H satisfies 
$(H1)$–
$(H4)$.
Let us begin with the understanding of interpolation spaces E α (
$\alpha \ge 0$) and interpolation operators A α given in (1.14) and (1.15). First, 
$E^0=L^2(\Omega)$ and 
$A^{0} = {\rm id}_{L^2(\Omega)}$. Clearly for 
$\alpha \ge 0$, A α is an isomorphism and its inverse 
$A^{-\alpha}: L^2(\Omega)\to E^{\alpha}$ can be denoted by
\begin{equation*}
A^{-\alpha}u=\sum_{j\ge 1}\lambda_j^{\frac{-\alpha}{2s}}a_j\varphi_j\quad \forall\, u=\sum_{j\ge 1}a_j\varphi_j\in L^2(\Omega).
\end{equation*}More generally, for convenience, we define formally, when it makes sense,
\begin{align}
A^\beta\Big(\sum_{j\ge 1}a_j\varphi_j\Big) =\sum_{j\ge 1}\lambda_j^{\frac{\beta}{2s}}a_j\varphi_j, \quad \forall\;\beta \in \mathbb R.
\end{align}Proposition 3.1. For any 
$s \in (0, 1)$, there hold 
$E^{2s}=\mathcal{H}^{2s}(\Omega)$, and 
$A^{2s}u=(-\Delta)^s u$ for all 
$u\in E^{2s}$.
Proof. For any 
$u\in \mathcal{H}^{2s}(\Omega)$, let
\begin{equation*}
u=\sum_{j\ge 1} a_j \varphi_j\in L^2(\Omega), \quad (-\Delta)^s u=\sum_{j\ge 1} b_j \varphi_j\in L^2(\Omega).
\end{equation*}Observe that
\begin{equation*}
b_j=\int_{\Omega}(-\Delta)^s u \varphi_j dx=\int_{\Omega} (-\Delta)^s \varphi_j udx=\int_{\Omega} \lambda_j \varphi_j udx=a_j\lambda_j.
\end{equation*} So 
$ (-\Delta)^s u= \sum_{j\ge 1} a_j \lambda_j \varphi_j$ and 
$\sum_{j\ge 1}a^2_j\lambda^2_j \lt \infty$ as 
$u\in \mathcal{H}^{2s}(\Omega)$. It follows that 
$\mathcal{H}^{2s}(\Omega)\subset E^{2s}$ and 
$(-\Delta)^s u = A^{2s}u$ for all 
$u\in \mathcal{H}^{2s}(\Omega)$.
Conversely, let 
$u\in E^{2s}$, then 
$u=\sum_{j\ge 1} a_j\varphi_j$ with 
$\sum_{j\ge 1} \lambda_j^2a_j^2 \lt \infty$. Let 
$\varphi \in \mathcal{H}^{2s}(\Omega)$, 
$\varphi=\sum_{j\ge 1} b_j\varphi_j$, we consider
\begin{equation}
f(\varphi)=\int_{\Omega}u (-\Delta)^s\varphi dx = \sum_{j\ge 1} a_jb_j\lambda_j.
\end{equation}Therefore
\begin{equation*}
|f(\varphi)| \le \Big|\sum_{j\ge 1} a_j^2\lambda_j^2\Big|^{\frac12}\Big|\sum_{j\ge 1} b_j^2\Big|^{\frac12}= \lVert u\rVert_{E^{2s}}| \varphi|_{2}.
\end{equation*} Using Riesz’s representation theorem and the density of 
$\mathcal{H}^{2s}(\Omega)$ in 
$L^2(\Omega)$ (all eigenfunctions belong to 
$\mathcal{H}^{2s}(\Omega)$), we deduce that there exists 
$v\in L^2(\Omega)$ such that
\begin{equation}
f(\varphi)=\int_{\Omega}v\varphi dx.
\end{equation} Combining (3.2) and (3.3), we get 
$(-\Delta)^su=v$, so 
$(-\Delta)^su\in L^2(\Omega)$ and 
$u\in \mathcal{H}^{2s}(\Omega)$.
Remark 3.2. By the definition of E α, if 
$0 \lt \alpha \lt 2s$, E α is the real interpolation space 
$[L^2(\Omega), E^{2s}(\Omega)]_{\alpha / 2s}$. Therefore, applying proposition 2.6 and proposition 3.1, 
$E^\alpha \subset L^{p+1}(\Omega)$ is a continuous embedding when 
$\frac{1}{p+1}\ge \frac12-\frac{\alpha}{N}$ and this embedding is compact provided that the strict inequality is valid, see [Reference Persson30] and [Reference Adams and Fournier1, Sections 7.22, 7.23].
By (1.17), (1.18), and remark 3.2, fix 
$\alpha\in (0, 2s)$ such that
\begin{equation}
N\left(\frac12-\frac{1}{\max\{p,q\}+1}\right) \lt \alpha \lt N\left(\frac{1}{\min\{p, q\}+1}-\frac{N-4s}{2N}\right),
\end{equation}so that
\begin{equation*}\mathbf{E}_{\alpha} :=E^{\alpha}\times E^{2s-\alpha}\subset X^*\end{equation*} and 
$\mathcal{E}: \mathbf{E}_{\alpha}\to \mathbb{R}$ given in (1.16) is well defined. In this section, we always consider 
$\mathbf{E}_{\alpha}$ with 
$\alpha \in (0, 2s)$ satisfying (3.4). 
$\mathbf{E}_{\alpha}$ is a Hilbert space with scalar product
\begin{equation*}
( (u_1, v_1), (u_2, v_2))_{\mathbf{E}_{\alpha}}=(u_1, u_2)_{E^\alpha}+(v_1, v_2)_{E^{2s-\alpha}}
\end{equation*} and 
$\mathbf{E}_{\alpha}=\mathbf{E}^+_{\alpha}\oplus \mathbf{E}^-_{\alpha}$, where
\begin{equation*}
\begin{split}
\mathbf{E}^+_{\alpha}&=\big\{(u, A^{-2s+\alpha}\circ A^{\alpha} u): u\in E^\alpha\big\},\quad \mathbf{E}^-_{\alpha}=\big\{(u, -A^{-2s+\alpha}\circ A^{\alpha} u): u\in E^\alpha\big\}.
\end{split}
\end{equation*} In the spirit of (3.1), we will write 
$A^{-2s+\alpha}\circ A^{\alpha}$ as 
$A^{-2s+2\alpha}$ for briefness, similarly we write 
$A^{2s-2\alpha}$ for 
$A^{-\alpha}\circ A^{2s-\alpha}$. If 
$z=(u,v)\in \mathbf{E}_{\alpha}$, denote
\begin{equation*}
\begin{split}
z^+&=(u^+,v^+)=\left(\frac{u+A^{2s-2\alpha}v}{2},\frac{v+A^{-2s+2\alpha}u}{2} \right)\in \mathbf{E}^+_{\alpha},\\
z^-&=(u^-,v^-)=\left(\frac{u-A^{2s-2\alpha}v}{2},\frac{v-A^{-2s+2\alpha} u}{2}\right)\in \mathbf{E}^-_{\alpha},
\end{split}
\end{equation*} then 
$z=z^++z^-$. Notice that although the quadratic part of 
$\mathcal{E}$ turns out to be strongly indefinite, it is positive definite in 
$\mathbf{E}^+_{\alpha}$ and negative definite in 
$\mathbf{E}^-_{\alpha}$. Consider the bilinear form 
$B: \mathbf{E}_{\alpha} \times \mathbf{E}_{\alpha} \to \mathbb{R}$ defined by
\begin{equation}
B((u,v), (\varphi, \psi))=\int_{\Omega}A^\alpha uA^{2s-\alpha}\psi\, +A^{2s-\alpha} vA^{\alpha}\varphi\, dx.
\end{equation} It is easy to see that B is continuous and symmetric, which induces a self-adjoint bounded linear operator 
$L: \mathbf{E}_{\alpha}\to \mathbf{E}_{\alpha}$ satisfying
\begin{equation}
B((u, v), (\varphi, \psi))=\big\langle L(u, v), (\varphi, \psi)\big\rangle_{\mathbf{E}_{\alpha}},\quad\forall\; (u, v), (\varphi, \psi)\in \mathbf{E}_{\alpha}.
\end{equation}Remark 3.3. By definition, we see easily that
\begin{equation*}
L(u, v)=(A^{2s-2\alpha}v, A^{2\alpha-2s}u ),\quad \forall\; (u, v)\in \mathbf{E}_{\alpha},
\end{equation*} and 
$L|_{\mathbf{E}_{\alpha}^+}={\rm id}_{\mathbf{E}_{\alpha}^+}$, 
$L|_{\mathbf{E}_{\alpha}^-}=-{\rm id}_{\mathbf{E}_{\alpha}^-}$.
Proposition 3.4. Let 
$p$, 
$q$ satisfy (1.17), (1.18), and 
$(u,v)\in \mathbf{E}_{\alpha}$ be a critical point of 
$\mathcal{E}$. Then 
$(u, v)$ is a L 1-weak solution to (1.1).
Proof. Since 
$(u,v)\in \mathbf{E}_{\alpha}$ is a critical point of 
$\mathcal{E}$, 
$\langle \mathcal{E}'(u,v), (\varphi,\psi)\rangle=0$ for all 
$(\varphi,\psi)\in \mathbf{E}_{\alpha}$, i.e.
\begin{equation*}
\int_{\Omega}A^\alpha uA^{2s-\alpha}\psi\, dx+\int_{\Omega}A^{2s-\alpha} vA^{\alpha}\varphi\, dx-\int_{\Omega}H_v(u,v)\psi\, dx-\int_{\Omega}H_u(u, v)\varphi\, dx=0.
\end{equation*}Taking φ = 0, we get
\begin{equation}
\int_{\Omega}A^\alpha uA^{2s-\alpha}\psi\, dx=\int_{\Omega}H_v(u,v)\psi\, dx.
\end{equation} It follows from proposition 3.1 that 
$\mathcal{T}_s(\Omega)\subset E^{2s}\subset E^{2s-\alpha}$. Let 
$\psi\in \mathcal{T}_s(\Omega)$ with 
$\psi=\sum_{j\ge 1} b_j \varphi_j$, and decompose 
$u=\sum_{j\ge 1} a_j \varphi_j$. By proposition 3.1,
\begin{equation}
\int_{\Omega}A^\alpha uA^{2s-\alpha}\psi\, dx=\sum_{j\ge 1} a_jb_j\lambda_j^{s}=\int_{\Omega} uA^{2s}\psi\, dx=\int_{\Omega} u(-\Delta)^s\psi\, dx.
\end{equation}Combining (3.7) with (3.8), we obtain
\begin{equation*}\int_{\Omega} u(-\Delta)^s\psi\, dx=\int_{\Omega}H_v(u,v)\psi\, dx.\end{equation*} Similarly, for any 
$\varphi \in \mathcal{T}_s(\Omega)$, there holds
\begin{equation*}\int_{\Omega} v(-\Delta)^s\varphi\, dx=\int_{\Omega}H_u(u,v)\varphi\, dx.\end{equation*}So (u, v) is a L 1-weak solution to (1.1).
To find a nontrivial critical point of 
$\mathcal{E}$ by theorem 2.8, our next step is to check the (PS) condition for 
$\mathcal{E}$. We start from the following compactness property.
Lemma 3.5. Let p, q satisfy (1.17), (1.18), 
$H$ satisfy 
$(H4)$ and α satisfy (3.4), then 
$H_z: \mathbf{E}_{\alpha}\to X$ is compact.
Proof. By 
$(H4)$ and remark 3.2, Hz is well defined from 
$\mathbf{E}_{\alpha}$ to 
$X^*$. Let 
$\{(u_j, v_j)\}$ be bounded in 
$\mathbf{E}_{\alpha}$, we aim to prove that 
$\{(H_u(u_j, v_j), H_v(u_j, v_j))\}$ has a convergent subsequence in X.
Since 
$\mathbf{E}_{\alpha}$ is a Hilbert space, up to a subsequence, there exists 
$(u, v)\in \mathbf{E}_{\alpha}$ such that 
$(u_j, v_j) \rightharpoonup (u, v)$ weakly in 
$\mathbf{E}_{\alpha}$. By virtue of the compact embeddings 
$E^\alpha\subset L^{p+1}(\Omega)$, 
$E^{2s-\alpha}\subset L^{q+1}(\Omega)$ (see remark 3.2 and (3.4)), there holds, up to a subsequence,
\begin{align*}
& u_j\to u ~\text{in}~L^{p+1}(\Omega), \;\; v_j\to v ~\text{in}~L^{q+1}(\Omega) \quad\mbox{and} \\ & u_j(x)\to u(x),\;\; v_j(x)\to v(x)\;\; {\rm a.e.~in }~\Omega..
\end{align*} Using [Reference Brezis9, Theorem 4.9], there exists 
$(\Phi_1, \Phi_2)\in X^*$ such that 
$|u_j|\le \Phi_1, |v_j|\le \Phi_2$ a.e. for all j. Applying 
$(H4)$ again, we obtain
\begin{equation*}
\begin{aligned}
&|H_u(u_j, v_j)|\le C\big(| \Phi_1|^p+| \Phi_2|^\frac{p(q+1)}{(p+1)}+1\big)\in L^{1+\frac1p}(\Omega),\\
&|H_v(u_j, v_j)|\le C\big(| \Phi_2|^q+| \Phi_1|^\frac{q(p+1)}{(q+1)}+1\big)\in L^{1+\frac1q}(\Omega).
\end{aligned}
\end{equation*}Let
\begin{equation*}
w_j=H_u(u_j, v_j),\quad w=H_u(u, v),\quad y_j=H_v(u_j, v_j),\quad y=H_v(u, v).
\end{equation*} Using the Lebesgue dominated convergence theorem, we get 
$(w_j, y_j)\to (w, y)$ in X. The proof is done.
Remark 3.6. Clearly, under the assumptions of lemma 3.5, 
$\mathcal{H}': \mathbf{E}_{\alpha}\to \mathbf{E}_{\alpha}^*$ is compact, which can be seen as follows
\begin{equation*}
\mathbf{E}_{\alpha}\underset{H_z}{\to} X\underset{\rm id}{\to} \mathbf{E}_{\alpha}^*.
\end{equation*} The first mapping above is compact by lemma 3.5, and the last one is continuous by the embedding 
$\mathbf{E}_{\alpha}\subset X^*$.
Lemma 3.7. Let p, q satisfy (1.17), (1.18), H satisfy 
$(H2)$ and 
$(H4)$. Then the functional 
$\mathcal{E}$ satisfies 
$(PS)_c$ condition for any 
$c\in \mathbb{R}$.
Proof. Let 
$\{(u_j, v_j)\}\subset \mathbf{E}_{\alpha}$ be any 
$(PS)_c$ sequence of 
$\mathcal{E}$ with 
$c\in\mathbb{R}$. By definition, when j goes to infinity,
\begin{equation}
\mathcal{E}(u_j,v_j)=\int_{\Omega}A^{\alpha}u_jA^{2s-\alpha}v_j\, dx-\int_{\Omega}H(u_j, v_j)\, dx=c+o(1),
\end{equation} and for all 
$(\varphi, \psi)\in \mathbf{E}_{\alpha}$,
\begin{equation}
\begin{aligned}
\langle \mathcal{E}'(u_j, v_j), (\varphi, \psi)\rangle& =\int_{\Omega}A^{\alpha}u_jA^{2s-\alpha}\psi \, dx+\int_{\Omega}A^{2s-\alpha}v_jA^{\alpha}\varphi\, dx-\int_{\Omega}H_u(u_j, v_j) \varphi\, dx\\
&\quad -\int_{\Omega}H_v(u_j, v_j) \psi dx\\
& =o(1)\lVert (\varphi, \psi) \lVert_{\mathbf{E}_{\alpha}}.
\end{aligned}
\end{equation} We first claim that 
$\{(u_j, v_j)\}$ is bounded in 
$\mathbf{E}_{\alpha}$. Take first
\begin{equation*}(\varphi, \psi)=\left(\frac{q+1}{p+q+2}u_j, \frac{p+1}{p+q+2}v_j\right).\end{equation*} From (3.10), (3.9), and 
$(H2)$, as j goes to infinity,
\begin{equation*}
\begin{split}
& \quad c+o(1)\|(u_j, v_j)\|_{\mathbf{E}_{\alpha}}+o(1)\\
& = \mathcal{E}(u_j,v_j)-\left\langle \mathcal{E}'(u_j, v_j), \left(\frac{q+1}{p+q+2}u_j, \frac{p+1}{p+q+2}v_j\right)\right\rangle\\
& = \frac{(p+1)(q+1)}{p+q+2}\\
&\quad\left(\frac{1}{p+1}\int_{\Omega} H_u(u_j, v_j) u_j\, dx+ \frac{1}{q+1}\int_{\Omega}H_v(u_j, v_j) v_j\, dx-\int_{\Omega}H(u_j, v_j)\, dx\right)\\
& \quad + \frac{pq-1}{p+q+2}\int_{\Omega}H(u_j, v_j)dx\\
&\quad \ge \frac{pq-1}{p+q+2}\int_{\Omega}H(u_j, v_j)dx -C.
\end{split}
\end{equation*} From 
$(H2)$ (see [Reference Felmer22, Lemma 1.1]), there exist 
$C_1, C_2 \gt 0$ such that
\begin{equation*}
H(u,v)\ge C_1(|u|^{p+1}+|v|^{q+1})-C_2, \quad \forall\; (u, v)\in \mathbf{E}_{\alpha}.
\end{equation*}Thus,
\begin{equation}
|u_j|_{p+1}^{p+1}+|v_j|_{q+1}^{q+1} \leq C+ o(1)\|(u_j, v_j)\|_{\mathbf{E}_{\alpha}}.
\end{equation}On the other hand, it follows from (3.10) with ψ = 0 that
\begin{equation}
\int_{\Omega}A^{2s-\alpha}v_jA^{\alpha}\varphi dx = \int_{\Omega}H_u(u_j, v_j) \varphi dx+o(1)\|\varphi\|_{E^\alpha},
\end{equation}and if φ = 0,
\begin{equation}
\int_{\Omega}A^{\alpha}u_jA^{2s-\alpha}\psi dx = \int_{\Omega}H_v(u_j, v_j) \psi dx+o(1)\|\psi\|_{E^{2s-\alpha}}.
\end{equation} Applying 
$(H4)$, Hölder inequality and 
$E^{\alpha}\subset L^{p+1}(\Omega)$, there holds
\begin{equation}
\int_{\Omega}\left|H_u(u_j, v_j) \varphi \right| dx\le C\left( |u_j|_{p+1}^p+|v_j|_{q+1}^{p(q+1)/(p+1)} + 1\right)\|\varphi\|_{E^\alpha}.
\end{equation} Since 
$A^{2s-\alpha}$ is an isomorphism from 
$E^{2s-\alpha}$ to 
$L^2(\Omega)$ and due to (3.12), (3.14), we have
\begin{align}
& \|v_j\|_{E^{2s-\alpha}} = |A^{2s-\alpha}v_j |_2= \underset{\|\varphi\|_{E^{\alpha}}=1}{\sup} \left| \int_{\Omega}A^{2s-\alpha}v_j A^{\alpha}\varphi dx \right| \nonumber\\ &\quad \le C\left( |u_j|_{p+1}^p+|v_j|_{q+1}^{p(q+1)/(p+1)}+1\right).
\end{align}Similarly, we obtain
\begin{equation}
\begin{aligned}
\|u_j\|_{E^{\alpha}} \le C\left( |v_j|_{q+1}^q+|u_j|_{p+1}^{q(p+1)/(q+1)}+1\right).
\end{aligned}
\end{equation} Combining (3.11) with (3.15) and (3.16), it is clear that 
$\{(u_j, v_j)\}$ is bounded in 
$\mathbf{E}_{\alpha}$. Therefore, up to a subsequence, there is 
$(u, v)\in \mathbf{E}_{\alpha}$ such that 
$(u_j, v_j) \rightharpoonup (u, v)$ weakly in 
$\mathbf{E}_{\alpha}$. Let
\begin{equation*}
w_j=H_u(u_j, v_j),\quad w=H_u(u, v),\quad y_j=H_v(u_j, v_j),\quad y=H_v(u, v).
\end{equation*} By lemma 3.5, up to a subsequence, 
$(w_j, y_j)\to (w, y)$ in X. In view of (3.12) and (3.13), let j tend to 
$\infty$, there holds
\begin{equation}
\begin{aligned}
\int_{\Omega}A^{2s-\alpha}vA^{\alpha}\varphi dx = \int_{\Omega}w \varphi dx,\quad \int_{\Omega}A^{\alpha}v A^{2s-\alpha}\psi dx = \int_{\Omega}y \psi dx.
\end{aligned}
\end{equation} Thanks to (3.12), (3.17), Hölder inequality and 
$E^{\alpha}\subset L^{p+1}(\Omega)$, we arrive at
\begin{equation*}
\begin{aligned}
\|v_j-v\|_{E^{2s-\alpha}}&=\underset{\|\varphi\|_{E^{\alpha}}=1}{\sup}\int_{\Omega}A^{2s-\alpha}(v_j-v)A^{\alpha}\varphi dx \\
&= \underset{\|\varphi\|_{E^{\alpha}}=1}{\sup}\left(\int_{\Omega}(w_j-w) \varphi dx+o(1)\|\varphi\|_{E^\alpha}\right)\\
&\le C\left(|w_j-w|_{1+\frac1p}+o(1)\right)
\end{aligned}
\end{equation*} which implies 
$v_j\to v$ in 
$E^{2s-\alpha}$. In the same way, we obtain 
$u_j\to u$ in E α.
In the sequel, we shall verify that 
$\mathcal{E}$ possesses the linking structure stated in theorem 2.8.
Lemma 3.8. Let 
$p$, 
$q$ satisfy (1.17), (1.18), H satisfy 
$(H1)$–
$(H4)$. Then there exist two linear, bounded and invertible operators 
$B_1$, 
$B_2: \mathbf{E}_{\alpha}\to \mathbf{E}_{\alpha}$ such that for all 
$\tau\ge 0$, 
$\widehat{B}_{\tau}=P_2 B_1^{-1}e^{\tau L}B_2 : \mathbf{E}_{\alpha}^-\to \mathbf{E}_{\alpha}^-$ is invertible, where P 2 is the projection of 
$\mathbf{E}_{\alpha}$ onto 
$\mathbf{E}_{\alpha}^-$. Moreover, let
\begin{equation*}e_1=\Big(\lambda_1^{-\frac{\alpha}{2s}}\varphi_1, \lambda_1^{\frac{-2s+\alpha}{2s}}\varphi_1\Big)\in \mathbf{E}_{\alpha}^+,\end{equation*} then there exist constants 
$R_2 \gt \rho \gt 0$, 
$R_1 \gt \rho/\|B_1^{-1}B_2 e_1\|_{\mathbf{E}_{\alpha}}$ and σ > 0 such that
\begin{equation*}
(G1)\;\; \mathcal{E}(z)\ge \sigma \gt 0 \;\mbox{on } S_\rho^+;\quad\quad \quad (G2)\;\; \mathcal{E}(z)\le 0 \;\mbox{on } \partial Q\end{equation*} where 
$S_\rho^+:=\{B_1z^+: z^+\in \mathbf{E}^+_{\alpha}, \|z^+\|=\rho\}$ and
\begin{equation*}
Q:=\big\{B_2(te_1+z^-): 0\le t\le R_1, \, z^-\in \mathbf{E}_{\alpha}^-,\, \|z^-\|_{\mathbf{E}_{\alpha}}\le R_2\big\}.
\end{equation*}Proof. Since 
$\frac1{p+1}+\frac{1}{q+1} \lt 1$, we can select 
$\mu, \nu\ge 1$ satisfying
\begin{equation}
\frac{1}{p+1} \lt \frac{\mu}{\mu+\nu},\quad \frac{1}{q+1} \lt \frac{\nu}{\mu+\nu}.
\end{equation} Let 
$R_1 \gt 1$, 
$0 \lt \rho \lt 1$ and define
\begin{equation*}
B_1(u, v)=(\rho^{\mu-1}u, \rho^{\nu-1}v), \;\; B_2(u, v)=(R_1^{\mu-1}u, R_1^{\nu-1}v),\quad \forall\; (u, v)\in \mathbf{E}_{\alpha},
\end{equation*} Since 
$L|_{\mathbf{E}_{\alpha}^-}=-id_{\mathbf{E}_{\alpha}^-}$, 
$\widehat{B}_\tau$ is invertible in 
$\mathbf{E}_{\alpha}^-$. According to the definitions of B 1 and B 2, we set
\begin{equation*}
S_\rho^+ =\big\{(\rho^{\mu-1}u^+, \rho^{\nu-1}v^+): \|(u^+, v^+)\|_{\mathbf{E}_{\alpha}}=\rho,\, z^+=(u^+, v^+)\in \mathbf{E}^+_{\alpha}\big\},
\end{equation*}and
\begin{equation*}
\begin{aligned}
Q=\Big\{t(&R_1^{\mu-1}\widetilde{\varphi}_1, R_1^{\nu-1}\lambda_1^{\frac{-s+\alpha}{s}}\widetilde{\varphi}_1)+(R_1^{\mu-1}u^-, R_1^{\nu-1}v^-): 0\le t\le R_1, \\
&\|(u^-, v^-)\|_{\mathbf{E}_{\alpha}}\le R_2,\, z^-=(u^-, v^-)\in \mathbf{E}^-_{\alpha}\Big\},
\end{aligned}
\end{equation*} where 
$\widetilde{\varphi}_1=\lambda_1^{-\frac{\alpha}{2s}}\varphi_1$.
Verification of 
$(G1)$. For any 
$(\rho^{\mu-1}u^+, \rho^{\nu-1}v^+)\in S_\rho^+$, by 
$(H3)$ and 
$\mathbf{E}_{\alpha}\subset L^{p+1}(\Omega)\times L^{q+1}(\Omega)$, when ρ is small, by (3.18), there exist 
$C_1, C_2 \gt 0$ such that
\begin{align*}
\begin{aligned}
\mathcal{E}(\rho^{\mu-1}u^+, \rho^{\nu-1}v^+)\ge \frac12 &\rho^{\mu+\nu-2}\|z^+\|^2_{\mathbf{E}_{\alpha}}-C_1\rho^{(\mu-1)(p+1)} \\ & \int_{\Omega}|u^+|^{p+1} dx -C_1\rho^{(\nu-1)(q+1)} \int_{\Omega}|v^+|^{q+1} dx\\
\ge \frac12 &\rho^{\mu+\nu-2}\|z^+\|^2_{\mathbf{E}_{\alpha}}-C_2\rho^{(\mu-1)(p+1)}\|z^+\|^{p+1}_{\mathbf{E}_{\alpha}} \\ & -C_2\rho^{(\nu-1)(q+1)} \|z^+\|_{\mathbf{E}_{\alpha}}^{q+1}\\
\ge \frac12 &\rho^{\mu+\nu}-C_2\rho^{\mu (p+1)}-C_2 \rho^{\nu (q+1)},
\end{aligned}
\end{align*} which yields that there exist 
$\sigma,\, \rho \gt 0$ satisfying 
$\mathcal{E}(\rho^{\mu-1}u^+, \rho^{\nu-1}v^+)\ge \sigma \gt 0$.
Verification of 
$(G2)$. We proceed by the following three steps.
Step 1. For any 
$(R_1^{\mu-1}u^-, R_1^{\nu-1}v^-)\in Q\cap \{t=0\}$, it follows from 
$(H1)$ that
\begin{equation*}
\begin{aligned}
\mathcal{E}(R_1^{\mu-1}u^-, R_1^{\nu-1}v^-)\le R_1^{\mu+\nu-2}\int_{\Omega}A^{\alpha}u^- A^{2s-\alpha}v^-\, dx \le -R_1^{\mu+\nu-2}\int_{\Omega}|A^{\alpha}u^-|^2 dx \le 0.
\end{aligned}
\end{equation*} Step 2. For any 
$\widetilde{z}_{R_1}=(R_1^{\mu}\widetilde{\varphi}_1 +R_1^{\mu-1}u^-,\, R_1^{\nu}\lambda_1^{\frac{-s+\alpha}{s}}\widetilde{\varphi}_1 +R_1^{\nu-1}v^-)\in Q\cap \{t=R_1\}$, we write 
$u^-=r\varphi_1+w$ where 
$w\in E^{\alpha}$ is orthogonal to φ 1 in 
$L^2(\Omega)$.
Let 
$r\ge 0$. By direct computations,
\begin{equation}
r+\lambda_1^{-\frac{\alpha}{2s}}t=\int_{\Omega}(t\widetilde{\varphi}_1+ u^-){\varphi}_1 dx\le |t\widetilde{\varphi}_1+ u^-|_{p+1}|\varphi_1|_\frac{p+1}{p}\le C|t\widetilde{\varphi}_1+ u^-|_{p+1}.
\end{equation} Set 
$\widetilde{z}_t=t(R_1^{\mu-1}\widetilde{\varphi}_1, R_1^{\nu-1}\lambda_1^{\frac{-s+\alpha}{s}}\widetilde{\varphi}_1)+(R_1^{\mu-1}u^-, R_1^{\nu-1}v^-)$ with 
$t\ge 0$. Using 
$(H2)$ (see [Reference Felmer22, Lemma 1.1]) and (3.19), we have
\begin{equation}
\begin{aligned}
\mathcal{H}(\widetilde{z}_t)& \ge C_1R_1^{(p+1)(\mu-1)}\int_{\Omega}|t\widetilde{\varphi}_1+u^-|^{p+1} dx + C_1R_1^{(q+1)(\nu-1)}\\
&\int_{\Omega}|t\lambda_1^{\frac{-s+\alpha}{s}}\widetilde{\varphi}_1 +v^-|^{q+1} dx-C_2\\
& \ge C_3 R_1^{(p+1)(\mu-1)}(r+\lambda_1^{-\frac{\alpha}{2s}}t)^{p+1}-C_2\\
& \ge C_4 R_1^{(p+1)(\mu-1)}t^{p+1}-C_2.
\end{aligned}
\end{equation} If r < 0, let 
$v^-=-A^{-2s+\alpha}\circ A^{\alpha}u^- = -A^{-2s+2\alpha} u^-$,
\begin{equation*}
\begin{aligned}
\int_{\Omega} v^- \varphi_1dx = \int_{\Omega}(-A^{-2s+2\alpha}u^-) \varphi_1 dx
& = \int_{\Omega}-A^{-2s+2\alpha} (r\varphi_1+w) \varphi_1 dx\\
&= -\lambda_1^{\frac{-s+\alpha}{s}}r-\int_{\Omega}\varphi_1A^{-2s+2\alpha} w dx\\
&= -\lambda_1^{\frac{-s+\alpha}{s}}r.
\end{aligned}
\end{equation*}Thus,
\begin{equation*}
\lambda_1^{\frac{-s+\alpha}{s}}(\lambda_1^{-\frac{\alpha}{2s}}t-r)=\int_{\Omega}(t\lambda_1^{\frac{-s+\alpha}{s}}\widetilde{\varphi}_1 +v^-)\varphi_1\, dx\le |t\lambda_1^{\frac{-s+\alpha}{s}}\widetilde{\varphi}_1 +v^-|_{q+1}|\varphi_1|_\frac{q+1}{q}.
\end{equation*}Hence, similar to (3.20), one has
\begin{equation}
\mathcal{H}(\widetilde{z}_t)\ge C_5 R_1^{(q+1)(\nu-1)}t^{q+1}-C_6.
\end{equation}Therefore, either (3.20) or (3.21) holds, which yields
\begin{equation*}
\mathcal{E}(\widetilde{z}_{R_1})\le {R_1^{\mu+\nu}}- \frac{R_1^{\mu+\nu-2}}{2}\|z^-\|_{\mathbf{E}_{\alpha}}^2-C_4 R_1^{(p+1)\mu}+C_2,
\end{equation*}or
\begin{equation*}
\mathcal{E}(\widetilde{z}_{R_1})\le {R_1^{\mu+\nu}}- \frac{R_1^{\mu+\nu-2}}{2}\|z^-\|_{\mathbf{E}_{\alpha}}^2-C_5 R_1^{(q+1)\nu}+C_6.
\end{equation*} In both cases, we can choose 
$R_2=2{R_1}$ large so that 
$\mathcal{E}(\widetilde{z}_{R_1})\le 0$.
Step 3. For any 
$\|(u^-, v^-)\|_{\mathbf{E}_{\alpha}}=R_2$, there holds
\begin{equation*}\widetilde{z}_t=t(R_1^{\mu-1}\widetilde{\varphi}_1, R_1^{\nu-1}\lambda_1^{\frac{-s+\alpha}{s}}\widetilde{\varphi}_1)+(R_1^{\mu-1}u^-, R_1^{\nu-1}v^-)\in Q\cap \{\|z^-\|_{\mathbf{E}_{\alpha}}=R_2\}.\end{equation*}Since either (3.20) or (3.21) holds, we get either
\begin{equation*}
\mathcal{E}(\widetilde{z}_t)\le {R_1^{\mu+\nu-2}t^2}- \frac{R_1^{\mu+\nu-2}}{2}R_2^2-C_4 R_1^{(p+1)(\mu-1)}t^{p+1}+C_2,
\end{equation*}or
\begin{equation*}
\mathcal{E}(\widetilde{z}_t)\le {R_1^{\mu+\nu-2}t^2}- \frac{R_1^{\mu+\nu-2}}{2}R_2^2-C_5 R_1^{(q+1)(\nu-1)}t^{q+1}+C_6.
\end{equation*} Choosing 
$R_2=2{R_1}$ large, we get 
$\mathcal{E}(\widetilde{z}_t)\le 0$ for 
$0\le t \le R_1$.
Lemma 3.9. Under the assumption 
$(H4)$ with subcritical 
$(p, q)$, any weak solution 
$(u, v)\in X^*$ of (1.1) is a classical solution.
Proof. Since 
$(u, v)\in X^*$ and 
$(H4)$, one has 
$H_z(u, v)\in X$, so
\begin{equation}
(u, v)\in\mathcal{W}^{2s, 1+\frac{1}{q}}(\Omega)\times \mathcal{W}^{2s, 1+\frac{1}{p}}(\Omega).
\end{equation} Consider first 
$p, q \gt \frac{2s}{N-2s}$, proposition 2.6 yields then that 
$u\in L^{\frac{N(q+1)}{Nq-2s(q+1)}}(\Omega)$ and 
$v\in L^{\frac{N(p+1)}{Np-2s(p+1)}}(\Omega)$. The assumption 
$(H4)$ deduces that 
$(u, v)\in \mathcal{W}^{2s, r_1}(\Omega)\times \mathcal{W}^{2s, t_1}(\Omega)$ where
\begin{equation*}
r_1=\min\left\{\frac{N(p+1)}{[Nq-2s(p+1)]q}, \frac{N(q+1)^2}{[Nq-2s(q+1)](p+1)q}\right\},
\end{equation*}and
\begin{equation*}
t_1=\min\left\{\frac{N(q+1)}{[Nq-2s(q+1)]p}, \frac{N(p+1)^2}{[Np-2s(p+1)](q+1)p}\right\}.
\end{equation*} Since 
$(p, q)$ is subcritical, we have
\begin{equation}
r_1 \gt r_0:=1+\frac{1}{q}, \quad t_1 \gt t_0:=1+\frac{1}{p}.
\end{equation} In the same way, whenever 
$r_n, t_n \lt \frac{N}{2s}$, let
\begin{equation*}
r_{n+1}=\min\left\{\frac{Nt_n}{(N-2st_n)q}, \frac{Nr_n(q+1)}{(N-2sr_n)(p+1)q}\right\}
\end{equation*}and
\begin{equation*}
t_{n+1}=\min\left\{\frac{Nr_n}{(N-2sr_n)p}, \frac{Nt_n(p+1)}{(N-2st_n)(q+1)p}\right\},
\end{equation*} then 
$(u, v)\in \mathcal{W}^{2s, r_{n+1}}(\Omega)\times \mathcal{W}^{2s, t_{n+1}}(\Omega)$, and 
$r_n \lt r_{n+1}, t_n \lt t_{n+1}$ due to (3.23). We claim that
\begin{align}
\mbox{there exists } n \text{ such that } (u, v)\in \mathcal{W}^{2s, r_n}(\Omega)\times \mathcal{W}^{2s, t_n}(\Omega) \text{ with } \min\{r_n, t_n\} \gt \frac{N}{2s}.
\end{align} This claim can be proved by contradiction. Assume for example 
$r_n\to r \le \frac{N}{2s}$, 
$t_n\to t \le \frac{N}{2s}$, other cases can be ruled out similarly. We can check that
\begin{equation*}r=\frac{N}{2s}\frac{pq-1}{pq+q}, \quad t=\frac{N}{2s}\frac{pq-1}{pq+p},\end{equation*} which contradicts 
$r \gt r_0$ since (p, q) is subcritical.
Therefore, after finite number of iterations, (3.24) holds true. Then 
$u, v\in C^s(\overline{\Omega})$, thanks again to proposition 2.6 and lemma 1.1.
Other situations with p or 
$q\le \frac{2s}{N-2s}$ can be handled easily using proposition 2.6 and lemma 1.1, so we omit the details.
Proof of theorem 1.4 completed
By lemmas 3.7–3.8, remarks 3.3, 3.6, applying theorem 2.8 and proposition 3.4, there exists a L 1-weak solution 
$(u, v)\in X^*$. By lemma 3.9, (u, v) is a classical solution, so the proof is completed.
4. Dual method and proof of theorem 1.5
This section is devoted to the proof of theorem 1.5. We state first some basic properties for H, 
$H^*$, 
$\mathcal{H}$, and 
$\mathcal{H}^*$. As
\begin{equation}
H(u, v)= \int_0^1 H_u(tu, tv)u + H_v(tu, tv)v dt,
\end{equation}there holds readily
Lemma 4.1. For 
$H$ satisfying 
$(H5)$, there are positive constants 
$A_1$, 
$A_2$ such that
\begin{equation*}
A_1 \left(|u|^{p+1}+ |v|^{q+1}\right)\le {H}(u, v)\le A_2\big(|u|^{p+1}+ |v|^{q+1}\big), \quad \forall\; (u, v) \in \mathbb R^2.
\end{equation*} Let 
$H^*$ be the Legendre–Fenchel transform of convex functional H, and 
$\mathcal{H}^*$ be the Legendre–Fenchel transform of 
$\mathcal{H}: X^*\to \mathbb{R}$, i.e.
\begin{equation*}
H^*(u, v)=\underset{(t, s)\in \mathbb{R}^2}{\sup}\left\{tu+sv-H(t, s)\right\},\quad \forall\; (u,v)\in \mathbb{R}^2,
\end{equation*}and
\begin{equation}
\mathcal{H}^*(f, g)=\underset{(u, v)\in X^*}{\sup} \left\{\int_{\Omega} (fu+gv) dx-\mathcal{H}(u, v)\right\},\quad \forall\; (f, g)\in X.
\end{equation}Using lemma 4.1 and lemma 2.9 (d)(e), recalling that the Legendre–Fenchel transform reverses the order by lemma 2.9 (a), we have
Lemma 4.2. Assume that 
$H$ is convex satisfying 
$(H5)$, then there are positive constants A 3, A 4 such that
\begin{equation*}
A_3\left(|u|^{1+\frac1p}+|v|^{1+\frac1q}\right)\le {H}^*(u, v)\le A_4\left(|u|^{1+\frac1p}+|v|^{1+\frac1q}\right), \quad \forall\; (u,v)\in \mathbb R^2
\end{equation*}and
\begin{equation*}
A_3 \Big(|f|_{1+\frac1p}^{1+\frac1p}+ |g|_{1+\frac1q}^{1+\frac1q}\Big)\le \mathcal{H}^*(f, g)\le A_4 \Big(|f|_{1+\frac1p}^{1+\frac1p}+ |g|_{1+\frac1q}^{1+\frac1q}\Big), \quad \forall\; (f, g) \in X.
\end{equation*} Consider 
$\nabla H$ and 
$\nabla H^*$ as mappings of 
$\mathbb R^2$. If 
$\nabla H$ is invertible, we denote by 
$(\nabla H)^{-1}$ its inverse.
Lemma 4.3. Assume that 
$H$ satisfies 
$(H5)$–
$(H6)$, then 
$H^*\in C^1(\mathbb{R}^2, \mathbb{R})$. Moreover 
$\nabla H : \mathbb{R}^2\to \mathbb{R}^2$ is a homeomorphism, and 
$(\nabla H)^{-1} = \nabla H^*$.
Proof. By lemma 4.1, there holds 
$\lim_{|z|\to \infty}{\frac{H(z)}{|z|}}=\infty$. Applying lemma 2.9 (c), 
$H^*\in C^1(\mathbb{R}^2)$.
As 
$\nabla H: \mathbb{R}^2\to \mathbb{R}^2$ is strictly monotone, 
$\nabla H$ is injective. We claim that 
$\nabla H$ is surjective. This fact should be known, we give a proof here for the sake of completeness. Indeed, for any 
$R_1 \gt 0$, by 
$(H5)$, there exists 
$R_2 \gt 0$ such that if 
$|z|\ge R_2$, then 
$|\nabla H(z)| \gt R_1$. Obviously, for any 
$w \in B_{R_1}$, 
$\text{deg} (\nabla H, B_{R_2}, w)$ is well defined and
\begin{equation}
\text{deg} (\nabla H, B_{R_2}, w)=\text{deg} (\nabla H, B_{R_2}, 0),
\end{equation} where ‘deg’ denotes topological degree, see [Reference Willem38, Appendix]. We define a homotopy 
$F: B_{R_2}\times [0,1]\to \mathbb{R}^2$,
\begin{equation}
F(z, t)=\nabla H\left(\frac{z}{1+t}\right)-\nabla H\left(\frac{-tz}{1+t}\right).
\end{equation} For any 
$z\in \partial B_{R_2}$ and 
$t\in [0, 1]$, we have 
$\frac{z}{1+t}\neq \frac{-tz}{1+t}$, hence 
$F(z, t)\neq 0$ as 
$\nabla H$ is injective. Therefore,
\begin{equation*}
\text{deg} (F(\cdot, 1), B_{R_2}, 0)=\text{deg} (F(\cdot, 0), B_{R_2}, 0)=\text{deg} (\nabla H, B_{R_2}, 0).
\end{equation*} According to Borsuk theorem [Reference Willem38, Theorem D.17], 
$\text{deg} (F(\cdot, 1), B_{R_2}, 0)$ is odd, which implies
\begin{equation*}\text{deg} (\nabla H, B_{R_2}, 0)\neq 0.\end{equation*}So by (4.3) we get
\begin{equation*}\text{deg} (\nabla H, B_{R_2}, w)\neq 0,\quad \forall \; w\in B_{R_1},\end{equation*} which deduces that there exists some 
$z\in B_{R_2}$ such that 
$\nabla H(z)=w$. Due to the arbitrariness of R 1, we complete the proof of this claim.
By Brouwer’s invariance of domain theorem, it’s known that any continuous bijection of Euclidean space is indeed a homeomorphism; we can get this fact also by noting that 
$\nabla H$ is proper, i.e. the preimage of any bounded set is bounded seeing 
$(H5)$, so 
$\nabla H$ is a closed mapping, hence a homeomorphism. Finally, it follows from lemma 2.9 (b) that 
$(\nabla H)^{-1}=\nabla H^*$.
Now we consider the functional version of the previous lemma. Recall that 
$H_z(u, v) = \nabla H(u, v)$ and 
$H_z^*(u, v) = \nabla H^*(u, v)$.
Lemma 4.4. 
$H_z$ is a homeomorphism from 
$X^*$ onto 
$X$, and 
$H_z^{-1}= H_z^*$.
Proof. We claim that Hz is continuous. Consider 
$(u_n, v_n)\to (u, v)$ in 
$X^*$, as before, up to a subsequence, there exists 
$(\Phi_1, \Phi_2)\in X^*$ such that 
$|u_n|\le \Phi_1$, 
$|v_n|\le \Phi_2$ a.e. for all n; 
$u_n \to u$ and 
$v_n \to v$ a.e. in Ω. Using 
$(H5)$ once more, we obtain
\begin{equation*}
|H_u(u_n, v_n)|\le C\Big(|\Phi_1|^p+ |\Phi_1|^{\alpha-1}|\Phi_2|^{\beta}\Big)\in L^{1+\frac1p}(\Omega)
\end{equation*}and
\begin{equation*}
|H_v(u_n, v_n)| \le C\Big(|\Phi_2|^q+ |\Phi_1|^{\alpha}|\Phi_2|^{\beta-1}\Big) \in L^{1+\frac1q}(\Omega).
\end{equation*} By Lebesgue’s theorem, we deduce that 
$H_z(u_n, v_n) \to H_z(u, v)$ in X, which means Hz is continuous.
Given any 
$(f, g)\in X$, 
$(u, v)= (\nabla H)^{-1}(f, g) = \nabla H^*(f, g) \in X^*$ thanks to 
$(H5)$, since 
$|u(x)|^p \leq C|f(x)|$ and 
$|v(x)|^q \leq C|g(x)|$. The bijectivity of Hz is an easy consequence by definition. The continuity of 
$H_z^{-1} = H^*_z$ can be proved as for Hz, we omit the details.
Finally, we state properties of the Legendre–Fenchel transform 
$\mathcal{H}^*$ for the Hamiltonian functional 
$\mathcal{H}$.
Lemma 4.5. 
$\mathcal{H}^{*}\in C^1(X, \mathbb{R})$. More precisely,
\begin{equation}
\mathcal{H}^{*}(f, g)=\int_{\Omega}H^*(f, g) dx, \quad
\left\langle (\mathcal{H}^{*})'(f, g), (\widetilde{f}, \widetilde{g}) \right\rangle=\int_{\Omega}H^*_u(f, g) \widetilde{f}+H^*_v(f, g)\widetilde{g} dx.
\end{equation}Proof. For any 
$(f, g)\in X$, by lemma 2.9 (b) and lemma 4.3, there exist functions 
$u, v$ such that
\begin{equation*}
H^*(f, g)=fu + gv -H(u, v),
\end{equation*} and 
$(u, v) = H_z^*(f, g) = \nabla H^*(f, g)$. By lemma 4.4, we have 
$(u, v)\in X^*$. In view of (4.2), it follows that
\begin{align*}
\begin{aligned}
\mathcal{H}^*(f, g) = & \underset{(\widetilde{u}, \widetilde{v})\in X^*}{\sup} \int_{\Omega}f\widetilde{u}+g\widetilde{v}-H(\widetilde{u}, \widetilde{v}) dx \\
& \le \int_{\Omega}\underset{(t, s)\in \mathbb{R}^2}{\sup} \left\{f(x)t+g(x)s-H(t, s) \right\} dx\\
& = \int_{\Omega}H^*(f, g) dx\\
& = \int_{\Omega} fu+gv-H(u,v) dx\\
& \le \mathcal{H}^*(f, g),
\end{aligned}
\end{align*} which implies the expression of 
$\mathcal{H}^*$. By 
$H^*\in C^1(\mathbb{R}^2, \mathbb{R})$ and lemma 4.4, we get the expression of 
$(\mathcal{H}^{*})'$ and 
$\mathcal{H}^{*}\in C^1(X^*, \mathbb{R})$.
Next, we explain how to find weak solutions to (1.1). Clearly, 
$(-\Delta)^s$ is an isomorphism of 
$\mathcal{W}^{2s, r}(\Omega)$ onto 
$L^r(\Omega)$. Therefore we can denote its inverse by 
$\mathcal{A}: L^r(\Omega)\to \mathcal{W}^{2s, r}(\Omega)$. Define 
$\mathcal{J}: X\to \mathbb{R}$ by
\begin{equation*}
\mathcal{J}(f, g)=\mathcal{H}^*(f, g)-\int_{\Omega}g\mathcal{A}f dx.
\end{equation*} Using proposition 2.6, 
$\mathcal{J}$ is well defined.
Lemma 4.6. 
$\mathcal{A}$ is self-adjoint in the following sense:
\begin{equation}
\int_{\Omega}g\mathcal{A}f dx=\int_{\Omega}f\mathcal{A}g dx, \quad \forall\; (f, g)\in X.
\end{equation}Proof. In fact, this is a direct consequence of the self-adjointness of 
$(-\Delta)^s$ over 
$X_0^s(\Omega)$ and density argument. Let 
$\{f_n\}, \{g_n\}\subset C_c^\infty (\Omega)$ satisfying 
$(f_n, g_n)\to (f, g)$ in X. Set 
$(u_n, v_n) = (\mathcal{A}f_n, \mathcal{A}g_n)\in X_0^{s}(\Omega)^2$. Then
\begin{equation*}
\begin{split}
\int_{\Omega}g_n\mathcal{A}f_n dx=\int_{\Omega} u_n (-\Delta)^sv_n dx=\int_{\Omega}v_n (-\Delta)^s u_n dx=\int_{\Omega}f_n\mathcal{A}g_n dx.
\end{split}
\end{equation*} Taking 
$n\to \infty$, we get (4.6).
We are now in position to present the dual method setting.
Proposition 4.7. If 
$(f, g)\in X$ is a critical point of 
$\mathcal{J}$, then 
${H}_z^{*}(f, g)$ is a 
$L^1$-weak solution to (1.1).
Proof. Since 
$(f, g)\in X$ is a critical point of 
$\mathcal{J}$,
\begin{align}
0=&\langle \mathcal{J}'(f, g), (f_1, g_1) \rangle= \nonumber\\ & \langle (\mathcal{H}^*)'(f, g), (f_1, g_1 )\rangle-\int_{\Omega}g\mathcal{A}f_1 dx-\int_{\Omega}g_1\mathcal{A}f dx,\quad \forall\; (f_1, g_1)\in X.
\end{align} Let 
$(u, v)=H_z^*(f, g)$, then by lemma 4.5, we have
\begin{equation}
\langle (\mathcal{H}^*)'(f, g), (f_1, g_1 )\rangle=\int_{\Omega}uf_1 dx +\int_{\Omega}v g_1 dx.
\end{equation} As 
$H_z^*=H_z^{-1}$, we obtain 
$(f, g)=(H_u(u, v), H_v(u,v))$. If we choose 
$g_1=0$, then using (4.7), (4.8), and lemma 4.6, it holds that
\begin{equation*}
\int_{\Omega}uf_1 dx=\langle (\mathcal{H}^*)'(f, g), (f_1, 0) \rangle=\int_{\Omega}g\mathcal{A}f_1 dx=\int_{\Omega}f_1\mathcal{A}g dx, \quad \forall\; f_1 \in L^{1+\frac{1}{p}}(\Omega).
\end{equation*} Hence 
$u=\mathcal{A}g$, i.e. 
$(-\Delta)^s u=g=H_v(u, v)$ in the weak sense. Similarly, 
$(-\Delta)^s v=H_u(u, v)$.
4.1. Proof of theorem 1.5 completed
Now we establish a mountain pass structure to get existence of nontrivial critical points of 
$\mathcal{J}$. Set
\begin{equation*}S_\rho:=\{(\rho^{k-1}f, \rho^{l-1}g): (f, g)\in X, \,\|(f, g)\|_X=\rho\}\end{equation*} where 
$k, l \gt 1$ satisfy
\begin{equation}
\frac{p}{p+1} \gt \frac{k}{k+l}, \quad \frac{q}{q+1} \gt \frac{l}{k+l}.
\end{equation} Assume 
$0 \lt \rho \lt 1$. For any 
$(\rho^{k-1}f, \rho^{l-1}g)\in S_\rho$, by lemma 4.2 and proposition 2.6,
\begin{equation*}
\begin{split}
\mathcal{J}(\rho^{k-1}f, \rho^{l-1}g)& =\mathcal{H}^*(\rho^{k-1}f, \rho^{l-1}g)-\rho^{k+l-2}\int_{\Omega}g\mathcal{A}f dx\\
& \ge C_1 \rho^{(k-1)(1+1/p)}|f|_{1+\frac1p}^{1+\frac1p}+C_1\rho^{(l-1)(1+1/q)}|g|_{1+\frac1q}^{1+\frac1q}-C_2\rho^{k+l}\\
& \ge C_3\rho^{\max\{k(1+1/p),\, l(1+1/q)\}}-C_2\rho^{k+l}.
\end{split}
\end{equation*} According to (4.9), there exist 
$\rho_0 \gt 0, \beta \gt 0$ such that 
$\mathcal{J}(\rho_0^{k-1}f, \rho_0^{l-1}g) \gt \beta$ if 
$(\rho_0^{k-1}f, \rho_0^{l-1}g)\in S_{\rho_0}$. On the other hand, we fix some 
$(f_0, g_0)\in X$ with
\begin{equation*}\int_{\Omega}f_0\mathcal{A}g_0 dx \gt 0.\end{equation*}By lemma 4.2,
\begin{equation*}
\mathcal{J}(\rho^{k}f_0, \rho^{l}g_0)\le C_4 \rho^{k(1+1/p)}|f_0|_{1+\frac1p}^{1+\frac1p}+C_4\rho^{l(1+1/q)} |g_0|_{1+\frac1q}^{1+\frac1q}-\rho^{k+l}\int_{\Omega}f_0\mathcal{A}g_0 dx.
\end{equation*} Consequently, due to (4.9), there exists 
$\rho_1 \gt \rho_0$ such that 
$\mathcal{J}(\rho_1^{k}f_0, \rho_1^{l}g_0) \lt 0$. Let
\begin{equation*}
\Gamma:=\{\gamma\in C([0, 1], X): \gamma(0)=0,\, \gamma(1)=(\rho_1^{k}f_0, \rho_1^{l}g_0) \}.
\end{equation*} It is clear that 
$\gamma([0, 1])\cap S_{\rho_0}\neq \varnothing$ for any 
$\gamma\in \Gamma$, from which we obtain a mountain pass structure of 
$\mathcal{J}$ around 0, and define the mountain pass level by
\begin{equation*}
c:=\underset{\gamma\in \Gamma}{\inf}\,\underset{t\in [0, 1]}{\max}\, \mathcal{J}(\gamma(t))\ge \beta \gt 0.
\end{equation*}Lemma 4.8. 
$\mathcal J$ satisfies the Palais–Smale condition.
Proof. Let 
$\{(f_n, g_n)\}\subset X$ be a Palais–Smale consequence of 
$\mathcal{J}$ at level c, that is,
\begin{equation}
\mathcal{J}(f_n, g_n)\to c\;\; \mbox{and }\; \mathcal{J}'(f_n, g_n)\to 0\;\; \text{in}~X^*, \quad \mbox{as } n\to \infty.
\end{equation}From (4.10) and (4.6), it follows that
\begin{equation}
\mathcal{H}^*(f_n, g_n)-\langle (\mathcal{H}^*)'(f_n, g_n), ((1-\theta) f_n, \theta g_n)\rangle=c+o(1)\|(f_n, g_n)\|_{X}+o(1)
\end{equation} where θ is given in 
$(H7)$. Let 
$(u_n, v_n)=H^*_z(f_n, g_n)$. According to lemma 2.9 (b) and lemma 4.5,
\begin{equation}
\mathcal{H}^*(f_n, g_n)=\int_{\Omega}f_nu_n+g_nv_ndx - \mathcal{H}(u_n, v_n),
\end{equation} which together with (4.11) and 
$(H7)$ implies that
\begin{equation*}
\begin{aligned}
&\quad c+o(1)\|(f_n, g_n)\|_{X}+o(1)\\
& = \int_{\Omega}f_nu_n+g_nv_n dx - \mathcal{H}(u_n, v_n)-\langle (\mathcal{H}^*)'(f_n, g_n), ((1-\theta) f_n, \theta g_n)\rangle\\
& = \theta\int_{\Omega} H_u(u_n, v_n)u_ndx+(1-\theta)\int_{\Omega} H_v(u_n, v_n)v_ndx- \mathcal{H}(u_n, v_n)\\
& \ge C_1\int_{\Omega}|u_n|^{p+1}dx+C_1\int_{\Omega}|v_n|^{q+1}dx-C_2|\Omega|.
\end{aligned}
\end{equation*} Combining the above inequality, (4.12), lemma 4.1, lemma 4.2, and Young’s inequality, 
$\{(f_n, g_n)\}$ is bounded in X.
Furthermore, by proposition 2.6, 
$\mathcal{A}: X\to X^*$ is compact, hence 
$\{(\mathcal{A}f_n, \mathcal{A}g_n)\}$ is compact in 
$X^*$. By (4.10),
\begin{equation}
\mathcal{J}'(f_n, g_n)=H^*_z(f_n, g_n)-(\mathcal{A}g_n, \mathcal{A}f_n)=o(1)\quad \text{in}~X^*.
\end{equation} Since 
$H^*_z$ is a homeomorphism from X onto 
$X^*$, 
$\{(f_n, g_n)\}$ is compact in X. The proof is done.
As 
$\mathcal J$ satisfies the (PS) condition, applying mountain pass theorem [Reference Ambrosetti and Rabinowitz2] and proposition 4.7, we get a nontrivial solution 
$(u, v)\in X^*$ of (1.1). Since 
$(H5)$ implies 
$(H4)$, by lemma 3.9, (u, v) is a classical solution.
5. Positive solutions for fractional Lane–Emden system
As a very special case of (1.1), we consider the system (1.24) under the subcritical assumption (1.25). Using lemma 3.9, we know that any energy solution of (1.26) belongs to 
$C^s(\overline{\Omega})$.
Proposition 5.1. Let 
$p,\, q$ satisfy (1.25). Then 
$u$ is an energy solution to (1.26) if and only if 
$(u, v)$ is a classical solution to (1.24), where 
$v=|(-\Delta)^{s}u|^{\frac1q-1}(-\Delta)^{s}u$.
Proof. First, assume that u is an energy solution to (1.26). It follows from (1.27) that
\begin{equation}
\int_{\Omega}v(-\Delta)^{s}\varphi dx=\int_{\Omega}|u|^{p-1}u\varphi dx\quad \forall\; \varphi\in\mathcal{T}_s(\Omega).
\end{equation} Next, since 
$v=|(-\Delta)^{s}u|^{\frac1q-1}(-\Delta)^{s}u$ and 
$u\in\mathcal{W}^{2s, 1+\frac1q}(\Omega)$, we have
\begin{equation}
\int_{\Omega}u(-\Delta)^{s}\varphi dx=\int_{\Omega}|v|^{q-1}v\varphi dx\quad \forall\; \varphi\in\mathcal{T}_s(\Omega).
\end{equation} By (5.1) and (5.2), (u, v) is a weak solution to (1.24). Since u is an energy solution to (1.26), using lemma 3.9, 
$u,\, v\in C^s(\overline{\Omega})$.
Conversely, if (u, v) is a classical solution to (1.24), (5.1), and (5.2) will hold, which implies
\begin{equation}
\int_{\Omega}|(-\Delta)^{s}u|^{\frac1q-1}(-\Delta)^{s}u(-\Delta)^{s}\varphi dx=\int_{\Omega}|u|^{p-1}u\varphi dx\quad \forall\; \varphi\in \mathcal{W}^{2s, 1+\frac1q}(\Omega).
\end{equation}Thus, u is an energy solution to (1.26).
We denote
\begin{equation}
\mathcal{R}(u):=\langle \mathcal{I}'(u), u\rangle=\left|(-\Delta)^su\right|_{1+\frac1q}^{1+\frac1q}-|u|_{p+1}^{p+1},\quad \forall\; u\in \mathcal{W}^{2s, 1+\frac1q}(\Omega).
\end{equation} For any 
$u\in\mathcal{N}_{\mathcal{I}}$ (see (1.28)), there holds
\begin{equation*}
\begin{aligned}
\langle \mathcal{R}'(u), u\rangle =\Big(1+\frac{1}{q}\Big)\left|(-\Delta)^su\right|_{1+\frac1q}^{1+\frac1q}-(p+1)|u|_{p+1}^{p+1} &=\Big(\frac{1}{q}-p\Big)\left|(-\Delta)^su\right|_{1+\frac1q}^{1+\frac1q},
\end{aligned}
\end{equation*}which implies
\begin{equation}
\begin{cases}
\langle \mathcal{R}'(u), u\rangle \lt 0,\quad \text{if}~~pq \gt 1,\\
\langle \mathcal{R}'(u), u\rangle \gt 0,\quad \text{if}~~pq \lt 1.
\end{cases}
\end{equation} Using the implicit function theorem, if pq ≠ 1, 
$\mathcal{N}_{\mathcal{I}}$ is a C 1-submanifold of 
$\mathcal{W}^{2s, 1+\frac1q}(\Omega)$ with codimension 1.
Lemma 5.2. If 
$pq\neq 1$, then 
$\mathcal{N}_{\mathcal{I}}$ is non-empty. Moreover, when 
$pq \gt 1$, 
$\mathcal{N}_{\mathcal{I}}$ is far away from zero and 
$\mathcal{I}$ constrained on 
$\mathcal{N}_{\mathcal{I}}$ has a positive lower bound; when 
$pq \lt 1$, 
$\mathcal{N}_{\mathcal{I}}$ is bounded in 
$\mathcal{W}^{2s, 1+\frac1q}(\Omega)$.
Proof. For any 
$u\in \mathcal{W}^{2s, 1+\frac1q}(\Omega)\backslash\{0\}$, when pq > 1 (resp. pq < 1), let tu be the maximum (resp. minimum) point of 
$\mathcal{I}(tu)$ for t > 0. Then 
$t_uu\in \mathcal{N}_{\mathcal{I}},$ so 
$\mathcal{N}_{\mathcal{I}}$ is non empty. For any 
$u\in \mathcal{N}_{\mathcal{I}}$, using (5.4) and 
$\langle \mathcal{I}'(u), u\rangle=0$, we have
\begin{equation}
\begin{aligned}
\left|(-\Delta)^su\right|_{1+\frac1q}^{1+\frac1q}=|u|_{p+1}^{p+1}\le C\left|(-\Delta)^su\right|_{1+\frac1q}^{p+1}.
\end{aligned}
\end{equation} Therefore if pq > 1, there exists some 
$C_0 \gt 0$ such that
\begin{equation}
\left|(-\Delta)^su\right|_{1+\frac1q} \gt C_0,\quad \forall\; u\in \mathcal{N}_{\mathcal{I}}.
\end{equation} Next, if 
$u\in \mathcal{N}_{\mathcal{I}}$, by (5.4),
\begin{equation*}
\begin{aligned}
\mathcal{I}(u)=\left(\frac{q}{q+1}-\frac1{p+1}\right)\left|(-\Delta)^su\right|_{1+\frac1q}^{1+\frac1q}.
\end{aligned}
\end{equation*} From (5.7) and pq > 1, it follows that 
$\mathcal{I}$ has a positive lower bound on 
$\mathcal{N}_{\mathcal{I}}$. When pq < 1, we know from (5.6) that 
$\mathcal{N}_{\mathcal{I}}$ is bounded in 
$\mathcal{W}^{2s, 1+\frac1q}(\Omega)$.
In the sequel, we consider the functional 
$\mathcal{I}$ constrained on 
$\mathcal{N}_\mathcal{I}$. A constrained critical point u of 
$\mathcal{I}|_{\mathcal{N}_\mathcal{I}}$, means that there exists a Lagrange multiplier 
$\lambda\in \mathbb{R}$ such that
\begin{equation}
\mathcal{I}'(u)=\lambda \mathcal{R}'(u) \quad \text{in}~\mathcal{W}^{2s, 1+\frac1q}(\Omega)^*,
\end{equation} where 
$\mathcal{W}^{2s, 1+\frac1q}(\Omega)^*$ is the dual space of 
$\mathcal{W}^{2s, 1+\frac1q}(\Omega)$. In particular, u is a critical point of 
$\mathcal{I}$ whenever λ = 0.
Lemma 5.3. Assume that 
$pq\neq 1$. Then any constrained critical point of 
$\mathcal{I}|_{\mathcal{N}_\mathcal{I}}$ is a critical point of 
$\mathcal{I}$.
Proof. Suppose that u is a constrained critical point, there is 
$\lambda\in\mathbb{R}$ such that (5.8) holds. Thus,
\begin{equation}
0=\langle \mathcal{I}'(u), u \rangle=\lambda \langle \mathcal{R}'(u), u\rangle=\lambda \left(\frac{1}{q}-p\right)\left|(-\Delta)^su\right|_{1+\frac1q}^{1+\frac1q},
\end{equation} which deduces λ = 0. Consequently, u is a critical point of 
$\mathcal{I}$
Lemma 5.4. For 
$pq\neq 1$, 
$\mathcal{I}|_{\mathcal{N}_\mathcal{I}}$ satisfies the (PS) condition.
Proof. Let 
$\{u_n\}$ be a (PS) sequence at level 
$c\in \mathbb{R}$ for 
$\mathcal{I}|_{\mathcal{N}_\mathcal{I}}$, that is 
$\{u_n\} \subset \mathcal{N}_\mathcal{I}$, 
$\lambda_n \in \mathbb R$ such that
\begin{equation}
\mathcal{I}(u_n)=c+o(1) \quad \mbox{and} \quad \mathcal{I}'(u_n)-\lambda_n\mathcal{R}'(u_n)=o(1)\;\; \text{in}~\mathcal{W}^{2s, 1+\frac1q}(\Omega)^*.
\end{equation}Therefore
\begin{equation*}
\begin{aligned}
c+o(1) =\mathcal{I}(u_n)-\frac{1}{p+1}\mathcal{R}(u_n) &=\left(\frac{q}{q+1}-\frac{1}{p+1}\right)\left|(-\Delta)^su_n\right|_{1+\frac1q}^{1+\frac1q}.
\end{aligned}
\end{equation*} By pq ≠ 1, one can conclude the boundedness of 
$\{u_n\}$ in 
$\mathcal{W}^{2s, 1+\frac1q}(\Omega)$, which yields also that 
$\{\langle \mathcal{R}'(u_n), u_n\rangle\}$ is bounded. Up to a subsequence, assume that
\begin{equation}
\langle \mathcal{R}'(u_n), u_n\rangle\to m \in \mathbb R.
\end{equation} Let m ≠ 0, (5.10) yields 
$o(1)=\langle\mathcal{I}'(u_n)-\lambda_n\mathcal{R}'(u_n), u_n\rangle=\lambda_nm+o(1),$ which yields 
$\lambda_n\to 0$. Up to a new subsequence, there exists 
$u\in \mathcal{W}^{2s, 1+\frac1q}(\Omega)$ such that 
$(-\Delta)^su_n\rightharpoonup (-\Delta)^su$ in 
$L^{\frac{q+1}{q}}(\Omega)$, and 
$u_n\to u$ in 
$L^{p+1}(\Omega)$, so 
$u\in \mathcal{N}_{\mathcal{I}}$. By direct computations,
\begin{equation*}
o(1)=|u_n|_{p+1}^{p+1}-|u|_{p+1}^{p+1}=\left|(-\Delta)^su_n\right|_{1+\frac1q}^{1+\frac1q}-\left|(-\Delta)^su\right|_{1+\frac1q}^{1+\frac1q},
\end{equation*} which together with the weak convergence of 
$(-\Delta)^s u_n$ in 
$L^{1 +\frac{1}{q}}(\Omega)$ yields that 
$u_n\to u$ in 
$\mathcal{W}^{2s, 1+\frac1q}(\Omega)$.
Let now m = 0. It is clear from (5.11) and 
$u_n\in \mathcal{N}_{\mathcal{I}}$ that
\begin{equation*}
\begin{aligned}
o(1)=\langle \mathcal{R}'(u_n), u_n\rangle=\Big(\frac{1}{q}-p\Big)\left|(-\Delta)^su_n\right|_{1+\frac1q}^{1+\frac1q},
\end{aligned}
\end{equation*} which deduces that 
$\left|(-\Delta)^su_n\right|_{1+\frac1q}\to 0.$ So we are done.
Proof of theorem 1.9 completed
(i) We apply a deformation lemma on C 1 manifold (see [Reference Bonnet7]) to 
$\mathcal{I}|_{\mathcal{N}_\mathcal{I}}$. By lemma 5.2, it is easy to see that 
$c_{\mathcal{I}} \lt 0$ if pq < 1 and 
$c_{\mathcal{I}} \gt 0$ if pq < 1. By lemma 5.4, 
$c_{\mathcal{I}}$ is attained by some 
$u\in \mathcal{N}_\mathcal{I}$. Hence, u is a constrained critical point of 
$\mathcal{I}|_{\mathcal{N}_\mathcal{I}}$. Using lemma 5.3, u is also a critical point of 
$\mathcal{I}$, hence an energy solution of (1.26). Let w be the weak solution of
\begin{equation*}
(-\Delta)^{s} w =|(-\Delta)^{s}u| \; \text{in }\;\Omega,\quad u=0~\;\text{in} \;\;\mathbb{R}^N\setminus\Omega.
\end{equation*} By lemma 2.1, 
$w\ge \pm u$, so 
$w\ge |u|$. Consequently, 
$\mathcal{I}(tw)\le \mathcal{I}(tu)$ for all t > 0. Thus there exists a unique 
$t_w \gt 0$ such that 
$t_ww\in \mathcal{N}_\mathcal{I}$ and 
$\mathcal{I}(t_ww)\le \mathcal{I}(u)=c_{\mathcal{I}}.$ Then tww is a minimizer for 
$c_{\mathcal{I}}$. By means of lemma 2.2, one has 
$t_ww \gt 0$. The proof can be concluded by proposition 5.1.
(ii) Assume that 
$u_1, u_2\in \mathcal{W}^{2s, 1+\frac1q}(\Omega)$ are two distinct positive energy solutions of (1.26). We know that 
$u_1, u_2\in C^s(\overline{\Omega})$. Denote
\begin{equation*}
\beta_0:=\sup\{l\in\mathbb{R}: u_2\ge lu_1 ~\text{a.e. in} ~\Omega \} \quad\text{and}\quad
\beta_1:=\sup\{l\in\mathbb{R}: u_1\ge lu_2 ~\text{a.e. in} ~\Omega \}.
\end{equation*} From lemma 2.2 and lemma 1.1, there exist 
$C_1, C_2 \gt 0$ such that 
$C_1\delta^s\le u_1, u_2\le C_2\delta^s$. Therefore 
$0 \lt \beta_0$, 
$\beta_1 \lt \infty$. In order to prove 
$u_1=u_2$, it suffices to prove that 
$\beta_0 \ge1$ and 
$\beta_1\ge 1$. Notice that
\begin{equation*}
\begin{aligned}
(-\Delta)^{s}\left(|(-\Delta)^{s}u_1|^{\frac1q-1}(-\Delta)^{s}u_1\right)&=u_1^{p}\ge \beta_1^{p} u_2^{p}\\
&=\beta_1^{p}(-\Delta)^{s}\left(|(-\Delta)^{s}u_2|^{\frac1q-1}(-\Delta)^{s}u_2\right)\\
&=(-\Delta)^{s}\left(|(-\Delta)^{s}(\beta_1^{pq}u_2)|^{\frac1q-1}(-\Delta)^{s}(\beta_1^{pq}u_2)\right).\\
\end{aligned}
\end{equation*} By comparison principle, we get 
$u_1\ge \beta_1^{pq}u_2$. So 
$\beta_1\ge\beta_1^{pq}$, hence 
$\beta_1\ge 1$. A similar argument gives 
$\beta_0\ge 1$. By proposition 5.1, the classical solution of (1.24) is unique.
(iii) Replacing 
$\mathcal{W}^{2s, 1+\frac1q}(B_R)$ by 
$\mathcal{W}_{{\rm rad}}^{2s, 1+\frac1q}(B_R)$, we can check that all the above steps work, so we have a positive classical radially symmetric solution.
Acknowledgements
The author completed this paper during his Ph.D. in East China Normal University. The author would like to thank Prof. Dong Ye for careful reading of the manuscript.
