\documentclass[reqno]{amsart}
\usepackage{hyperref}

\AtBeginDocument{{\noindent\small
\emph{Electronic Journal of Differential Equations},
Vol. 2017 (2017), No. 113, pp. 1--10.\newline
ISSN: 1072-6691. URL: http://ejde.math.txstate.edu or http://ejde.math.unt.edu}
\thanks{\copyright 2017 Texas State University.}
\vspace{8mm}}

\begin{document}
\title[\hfilneg EJDE-2017/113\hfil Nonlinear elliptic systems]
{Existence of positive solutions to nonlinear elliptic systems 
involving gradient term and reaction potential}

\author[A. Attar, R. Bentifour \hfil EJDE-2017/113\hfilneg]
{Ahmed Attar, Rachid Bentifour}

\address{Ahmed Attar \newline
Laboratoire d'Analyse Nonlin\'eaire et Math\'ematiques Appliqu\'ees,
D\'epartement de \newline Math\'ematiques,
Universit\'e Abou Bakr Belka\"{\i}d, Tlemcen,
Tlemcen 13000, Algeria}
\email{ahm.attar@yahoo.fr}


\address{Rachid Bentifour \newline
Laboratoire d'Analyse Nonlin\'eaire et Math\'ematiques Appliqu\'ees,
D\'epartement de \newline Math\'ematiques,
Universit\'e Abou Bakr Belka\"{\i}d, Tlemcen,
Tlemcen 13000, Algeria}
\email{rachidbentifour@gmail.com}

\dedicatory{Communicated by Jesus Ildefonso Diaz}

\thanks{Submitted July 6, 2016. Published April 26, 2017.}
\subjclass[2010]{35K15, 35K55, 35K65, 35B05, 35B40}
\keywords{Elliptic system; Schauder fixed point theorem; gradient dependance}

\begin{abstract}
 In this note we study  the elliptic system
 \begin{gather*}
 -\Delta u  =  z^p+f(x) \quad \text{in }\Omega , \\
 -\Delta z = |\nabla u|^{q}+g(x) \quad \text{in }\Omega , \\
 z,u > 0 \quad \text{in }\Omega ,\\
 z=u= 0 \quad \text{on }\partial \Omega,
 \end{gather*}
 where $\Omega \subset \mathbb{R}^{N}$ is a bounded domain, $p>0$, $0<q\le 2$
 with $pq<1$ and $f,g$ are two nonnegative measurable functions.
 The main result of this work is to analyze the interaction between the potential
 and the gradient terms in order to get the existence of a positive solution.
\end{abstract}

\maketitle
\numberwithin{equation}{section}
\newtheorem{theorem}{Theorem}[section]
\newtheorem{lemma}[theorem]{Lemma}
\newtheorem{remark}[theorem]{Remark}
\newtheorem{definition}[theorem]{Definition}
\allowdisplaybreaks

\section{Introduction}

In this work we study the elliptic system
\begin{equation}
\begin{gathered}
-\Delta u  =  z^p+f(x) \quad \text{in }\Omega , \\
-\Delta z  =  |\nabla u|^{q}+g(x) \quad \text{in }\Omega , \\
z,u  > 0 \quad \text{in }\Omega , \\
z=u  =  0 \quad \text{on }\partial \Omega,
\end{gathered}  \label{P}
\end{equation}
where $p>0$, $0<q\le 2$ and $f,g$ are non negative measurable functions. 
Our goal is to prove the existence of a positive solution under some suitable
 hypotheses on the data.

Elliptic systems with gradient appear when dealing of the modeling of an 
electrochemical engineering problem, see \cite{F}.
We refer also to \cite{C} and \cite{Diez} for other applications of these 
class of systems.

Recently, in \cite{BOP}, the authors consider the system
\begin{gather*}
-\operatorname{div}(b(x,z)\nabla u)= f(x) \quad \text{in }\Omega , \\
-\operatorname{div}(a(x,z)\nabla z)=b(x,z)|\nabla u|^{2}\quad \text{in }\Omega , \\
z=u=0 \quad \text{on }\partial \Omega ,
\end{gather*}
where $a(x,s), b(x,s)$ are positive and coercive Caratheodory functions. 
Under the hypothesis that $f\in L^m(\Omega)$ with $m\ge \frac{2N}{N+2} $,
 they proved the existence and regularity of a positive solution.

When the gradient appears as an absorption term, the system becomes
\begin{gather*}
-\operatorname{div}(a(x,z)\nabla u)= f\quad \text{in }\Omega , \\
-\operatorname{div}(b(x,z)\nabla z)+K(x,z)|\nabla u|^{2}=g \quad \text{in }\Omega , \\
z=u=0 \quad  \text{on }\partial \Omega.
\end{gather*}
This system was studied in \cite{BOPu}. 
It is clear that in this case a priori estimate can be obtained easily 
and existence is allowed for $L^1$ data.

In \cite{AS} the authors deal with the so-quoted 
``elliptic system with triangular structure'', namely they consider the system
\begin{equation}\label{tri}
\begin{gathered}
-\Delta u_i=f_i(x,u,\nabla u)+F_i(x)  \quad\text{in }\Omega,\\
u_i=0  \quad \text{on }\partial \Omega ,
\end{gathered}
\end{equation}
where $\sum_{1\le j \le i}f_j\le 0$ and $1\le i\le m$.  
It is clear that under the above condition on $\{f_i\}_i$, the gradient terms 
in \eqref{tri} have an absorption effect and then a priori estimates can be 
inferred directly. We refer also to \cite{Gur} where a variation of the 
system \eqref{tri} is studied in a radial domain with blow-up boundary conditions.

In the case where $p=1$, $q=2$, $g=f=0$ and $\Omega=B_R(0)$, the system is reduced 
to 
\begin{equation}\label{diez}
-\Delta u =z, -\Delta z =|\nabla u|^{2} \quad \text{ in }B_R(0).
\end{equation}
Using the radial structure of the previous system, the authors in \cite{Diez} 
were able to reduce \eqref{diez} to the study of first-order ODEs and then 
they proved existence and uniqueness of a nonnegative large radial solution 
to \eqref{diez}.

The parabolic version of problem \eqref{diez} is studied as a modification 
of the classical Boussinesq approximation for buoyancy-driven flows of 
viscous incompressible fluids, see \cite{Diez2,Diez3} for more details in this 
direction.

In our case the situation is quite different and we need to  analyze 
the approximated system to get a priori estimates.
This note is organized as follows, in Section \ref{s2} we introduce 
some preliminaries results, like the functional setting and some other useful tools.
Section \ref{s3} is dedicated to prove our main existence result.
 Notice that, as a consequence of the existence results we will be able 
to show an existence result for the Bi-Laplacian operator with gradient 
term $|\nabla u|^q$ under suitable hypothesis on $q$.

In Section \ref{s4} we give some optimal conditions and we collet some open 
problems. In the first subsection we prove non existence results, that, 
in some sense, justify the conditions imposed on $p$ and $q$  to get the 
existence of positive solution for all $f,g\in L^2(\Omega)$. Some interesting 
open problems related to \eqref{P} are given in the last subsection.

\section{Preliminaries}\label{s2}

In this section, we begin by recalling some useful results. 
Since we are considering problems with general datum, we will use the concept
 of weak solution.

\begin{definition}\label{def:faible} \rm
Let $f, g \in L^1(\Omega)$ be nonnegative functions. 
Assume that $p>0$ and $0<q\le 2$, we say that $(u,z)\in L^1(\Omega)\times L^1(\Omega)$ 
is a weak solution of \eqref{P}, if $|\nabla u|^q\in L^1(\Omega)$, $z^p \in L^1(\Omega)$ 
and for all $\varphi\in \mathcal{C}^\infty_0(\Omega)$, we have
\begin{equation}\label{eq1f}
 \int_\Omega (-\Delta \varphi) u
= \int_\Omega z^p\varphi+ \int_\Omega f \varphi,\quad\text{and} \int_\Omega (-\Delta \varphi)z
= \int_\Omega| \nabla u| ^q\varphi+ \int_\Omega g\varphi.
\end{equation}
\end{definition}

Notice that, since $(z^p+f)\in L^1(\Omega)$, then we can see $u$ as a weak solution 
of the problem 
$$
-\Delta u  =  z^p+f(x)   \text{ in  }\Omega, \quad u=0  \text{  on   }  \partial\Omega.
$$
Therefore, by a result in \cite{BZ} we know that $u\in W^{1,\sigma}_0(\Omega)$ for all 
$ \sigma<\frac{N}{N-1}$, more precisely
we will use the following result proved in the appendix of \cite{BP}.

\begin{lemma}\label{BP}
 Assume that $u\in L^1_{\rm loc}(\Omega)$ is such that $\Delta u\in
 L^1_{\rm loc}(\Omega)$, then for all $p\in [0,\frac{N}{N-1})$, 
and for any open sets
 $\Omega_1\subset \Omega_2\subset \overline{\Omega}_2\subset\Omega$, there exists a
 positive constant $C\equiv C(p,\Omega_1,\Omega_2,N)$ such that
 \begin{equation}\label{tt1}
 \|u\|_{W^{1,p}(\Omega_1)}\le C\int_{\Omega_2} (|u|+|\Delta u|)\,dx.
 \end{equation}
Moreover if $u\in L^1(\Omega)$ and $\Delta u\in L^1(\Omega)$, then the above
 estimate holds globally in the domain $\Omega$.
\end{lemma}

To prove the main existence result, we use the next Schauder fixed point Theorem.

\begin{theorem}\label{pointfixe}
Let $T$ be a continuous and compact mapping of a Banach space  into itself, 
such that the set
$$
\{x\in X: x=\lambda Tx  \text{ for some  } 0\le \lambda \le 1\}
$$
is bounded. Then $T$ has a fixed point.
\end{theorem}

\section{Existence results}\label{s3}

We begin by considering an approximating problem with regular data. 
More precisely we have the next existence result.

\begin{theorem}\label{exist1}
Let $\Omega\subset \mathbb{R}^N$ be a bounded domain and suppose that $f,g \in L^{\infty}(\Omega)$ 
are nonnegative functions. Then for all $p>0$, $0<q\le 2$ and for all $\varepsilon>0$, 
the system
\begin{equation}\label{P1}
\begin{gathered}
-\Delta u=\frac{z^p}{1+\varepsilon z^p}+f(x)  \quad \text{in }\Omega , \\ 
-\Delta z=\frac{|\nabla u|^{q}}{1+\varepsilon|\nabla u|^{q}}+g(x)\quad\text{ in }\Omega , \\
z=u=0 \quad \text{on }\partial \Omega ,
\end{gathered}
\end{equation}
has a positive solution $(u, z)\in (W_0^{1,2}(\Omega))^2\cap (L^\infty(\Omega))^2$.
\end{theorem}

\begin{proof} 
We will use a fixed point argument. Let $u\in L^1(\Omega)$ be fixed and define 
$(\varphi, z)$ to be the unique solution of the system
\begin{equation}\label{P100}
\begin{gathered}
-\Delta \varphi =  h_\varepsilon(x,u)=\frac{ u^{p}_+}{1+\varepsilon u^{p}_+}+ f(x)
\quad \text{in }\Omega, \\
-\Delta z  =  \frac{|\nabla \varphi|^{q}}{1+\varepsilon|\nabla \varphi|^{q}}+g(x)
\quad \text{in }\Omega,\\
\varphi=z = 0  \quad  \text{on }\partial \Omega.
\end{gathered}
\end{equation}
It is clear that $h_\varepsilon\in L^{\infty}(\Omega)$, thus 
$\varphi\in \mathbb{X}(\Omega)\equiv \mathcal{C}^{1,\sigma}(\Omega)
\cap L^\infty(\Omega)\cap W^{1,2}_0(\Omega)$. Thus $z$ is well defined and 
$z\in \mathbb{X}(\Omega)$.
Hence we can define the operator 
$T:  L^1(\Omega) \to L^1(\Omega)$, $T(u)=z$.
We claim that $T$ satisfies the conditions of Schauder fixed point Theorem. 
The proof of the claim will be done in several steps.
\smallskip

\noindent\textbf{Step I:} $T$ is continuous. 
Let $\{u_n\}_n \subset L^1(\Omega)$ be such that $u_n\to u$ strongly in 
$L^1(\Omega)$. We set $z_n=T(u_n)$ and $z=T(u)$, then $(\varphi_n,z_n)$  and 
$(\varphi,z)$ satisfy
\begin{equation}
\begin{gathered}
-\Delta \varphi_n =  h_\varepsilon(x,u_n) \quad \text{ in }\Omega, \\
-\Delta z_n  =  \frac{|\nabla \varphi_n|^{q}}{1+\varepsilon|\nabla \varphi_n|^{q}}+g(x)
\quad \text{in }\Omega,\\
\varphi_n=z_n = 0 \quad \text{on }\partial \Omega,
 \end{gathered} \label{Pn1}
\end{equation}
 and
\begin{equation}
\begin{gathered}
-\Delta \varphi =  h_\varepsilon(x,u) \quad \text{in }\Omega, \\
-\Delta z  =  \frac{|\nabla \varphi|^{q}}{1+\varepsilon|\nabla \varphi|^{q}}+g(x)
\quad \text{in }\Omega,\\
\varphi=z = 0 \quad \text{on }\partial \Omega,
 \end{gathered}
\end{equation}
Notice that
$$
-\Delta (\varphi_n-\varphi)=h_\varepsilon(x,u_n)-h_\varepsilon(x,u).
$$
Taking into account that $|h_\varepsilon(x,s)|\le C(\varepsilon)$ and since $u_n\to u$ strongly 
in $L^1(\Omega)$, it holds $h_\varepsilon(.,u_n)\to h_\varepsilon(.,u)$ strongly in $L^a(\Omega)$ 
for all $a>0$. Moreover, using H\"{o}lder and Poincar\'e inequalities, 
it follows that
$$
\int_\Omega |\nabla (\varphi_n-\varphi)|^2\,dx
\le \int_\Omega (h_\varepsilon(x,u_n)-h_\varepsilon(x,u))^2\,dx \to 0\quad \text{as   }n\to \infty.
$$
Hence $\varphi_n\to \varphi$   strongly in  $W_0^{1,2}(\Omega)$.  
Now going back to the problems of $z_n$ and $z$ and since
$$
\frac{|\nabla \varphi_n|^{q}}{1+\varepsilon|\nabla \varphi_n|^{q}} 
\to \frac{|\nabla \varphi|^{q}}{1+\varepsilon|\nabla \varphi|^{q}}
\quad \text{strongly in $L^a(\Omega)$ for all }a>1, 
$$
it follows that $z_n \to z$   strongly in  $W_0^{1,2}(\Omega)$.
 Hence $z_n\to z$   strongly in    $L^1(\Omega)$.  Then $T$ is continuous.
\smallskip

\noindent\textbf{Step II:} $T$ is compact.
Consider now a sequence $\{u_n\}_n$ such that $\|u_n\|_{L^1(\Omega)}\le C$. 
As above we set $z_n=T(u_n)$ and define $\varphi_n$ as the unique solution 
of the first problem in \eqref{Pn1}. It is clear that $\{\varphi_n\}_n$ 
is bounded in $L^\infty(\Omega)\cap W_0^{1,2}(\Omega)$. Then up to a subsequence not relabeled, 
$\varphi_n\rightharpoonup \varphi$ weakly in $W_0^{1,2}(\Omega)$ and strongly in $L^a(\Omega)$ 
for all $a<2^*$. Thus $\varphi\in W_0^{1,2}(\Omega)\cap L^\infty(\Omega)$. 
Using $(\varphi_n-\varphi)$ as a test function in the equation of $\varphi_n$, 
there results that
$$
\int_\Omega |\nabla (\varphi_n-\varphi)|^2\,dx
\le \int_\Omega \nabla\varphi \nabla(\varphi-\varphi_n)\,dx+o(1).
$$
Since $\varphi_n\rightharpoonup \varphi$ weakly in $W_0^{1,2}(\Omega)$, it follows that
 $\varphi_n\to \varphi$ strongly in $W_0^{1,2}(\Omega)$.

Hence up to a subsequence, we reach that $z_n\to z$ strongly in $W_0^{1,2}(\Omega)$ 
and in particular in $L^1(\Omega)$. Hence $T$ is a compact operator.


To complete the proof of the claim we just have to show that 
$T(B_R(0))\subset B_R(0)$ for some ball $B_R(0)\subset L^1(\Omega)$.
Notice that by using $\varphi_n$ as test function in the first equation 
of \eqref{Pn1}, it follows that $\|\varphi_n\|_{W_0^{1,2}(\Omega)}<C(\varepsilon, \Omega)$. 
On the other hand, using $z_n$ as a test function in the second equation 
of \eqref{Pn1}, we obtain $\|z_n\|_{W_0^{1,2}(\Omega)}<C_1(\varepsilon, \Omega)$.
Thus $\|z_n\|_{L^1(\Omega)}<C_2(\varepsilon, \Omega)$. Hence choosing $R>C_2(\varepsilon, \Omega)$, 
we conclude that if $\|u\|_{L^1(\Omega)}\le R$, then $\|z\|_{L^1(\Omega)}\le R$. 
Hence the claim follows.

Therefore, by Schauder fixed point Theorem, we obtain the existence of $u$
such that $T(u)=u$. It is clear that $u>0$ in $\Omega$, hence $(u,z)$ solves 
the system \eqref{P1}. Now, by classical regularity results and the 
previous a priori estimates we obtain easily that
$(u, z)\in (W_0^{1,2}(\Omega))^2\cap (L^\infty(\Omega))^2$.
\end{proof}

Now, we are able to state the main result in this note.

\begin{theorem}\label{existence}
Let $\Omega \subset \mathbb{R}^N$  be a bounded domain. Suppose that $p>0$, $0<q<2$ 
with $pq<1$, then for all $f, g \in L^2(\Omega)$,
then system \eqref{P} has a positive solution $(u,z)$ such that
$(u, z^{\frac{\alpha+1}{2}})\in W_0^{1,2}(\Omega)\\times W_0^{1,2}(\Omega)$ where $\alpha>0$ 
satisfies $p<\frac{\alpha+1}{2}<\frac{1}{q}$.
\end{theorem}

\begin{proof}
We proceed by approximation. Let $\{f_n\}_n, \{g_n\}_n\subset
L^{\infty}(\Omega)$ be such that $f_n\uparrow f$ and $g_n\uparrow g$ strongly 
in $L^2(\Omega)$. Let $(u_n, z_n)\in [W_0^{1,2}(\Omega)\cap L^\infty(\Omega)]^2$ be the unique
positive solution to the approximate system
\begin{equation}\label{Pn-zn}
\begin{gathered}
-\Delta z_n=\frac{|\nabla u_n|^{q}}{1+\frac{1}{n}|\nabla u_n|^{q}}+g_n(x)
\quad \text{in }\Omega , 
\end{gathered} 
\end{equation}
\begin{equation} \label{Pn-un}
\begin{gathered}
z_n=0 \quad  \text{on  }\partial \Omega; \\[4pt]
-\Delta u_n=\frac{z_n^p}{1+z^p_n}+ f_n(x) \quad \text{in  }\Omega , \\
u_n=0 \quad \text{on  }\partial \Omega.
\end{gathered} 
\end{equation}
Notice that the existence of $(u_n,z_n)$ follows by using Theorem \ref{exist1}.

Fix $\alpha>0$ such that the above condition on $\alpha$ holds. Using $z_n^\alpha$ 
as a test function in \eqref{Pn-zn}, it follows that
$$
-\int_{\Omega } \Delta z_nz_n^{\alpha }\,dx
=\int_\Omega \frac{|\nabla u_n|^{q}}{1+\frac{1}{n}|\nabla u_n|^{q}}
z_n^{\alpha }\,dx+\int_\Omega g_n z_n^{\alpha }\,dx.
$$
Thus by Young and H\"{o}lder inequalities we obtain
\begin{align*}
&\frac{4\alpha }{(\alpha +1)^{2}}\int_\Omega 
\big| \nabla z_n^{\frac{\alpha +1}{2}}\big| ^{2}\,dx \\
&\leq \frac{q}{2} \int_\Omega | \nabla u_n| ^{2}\,dx
+ \frac{2-q}{2}{\int_\Omega }z_n^{\alpha \frac{2}{2-q}}\,dx
+\| g_n\|_{L^2(\Omega)}\Big(\int_\Omega z_n^{2\alpha}\ \,dx\Big)^{1/2}.
\end{align*}
Let us estimate each term in the left hand side of the previous inequality.

Using Sobolev and Young inequalities we easily reach
$$
\Big(\int_\Omega z_n^{2\alpha}\ \,dx\Big)^{1/2}
\le \varepsilon\int_\Omega |\nabla z_n^{\frac{\alpha+1}{2}}|^2\ \,dx+ c(\varepsilon).
$$
Furthermore, by the fact that $\alpha<\frac{2-q}{q}$, it follows that 
$\frac{2\alpha }{2-q}<2^{\ast }\frac{\alpha +1}{2}$. Hence
$$
\int_\Omega z_n^{\alpha \frac{2}{2-q}}\,dx
\le C(\Omega) \Big(\int_\Omega z_n^{2^{\ast }\frac{\alpha +1}{2}}\,dx
 \Big)^{1/\beta}
\le C(\Omega) \Big(\int_\Omega |\nabla z_n^{\frac{\alpha +1}{2}}|^2\,dx
\Big)^{\frac{2^*}{2\beta}},
$$
where $\beta= \frac{2^{\ast }(\alpha +1)(2-q)}{4\alpha }$. 
It is clear that $ \frac{2^*}{2\beta}<1$. 
Therefore, combining the above estimates we have
\begin{equation}\label{ahmed}
\int_\Omega | \nabla z_n^{\frac{\alpha +1}{2}}| ^{2}\,dx
\leq C_1\int_\Omega | \nabla u_n| ^{2}\,dx+C_2.
\end{equation}
Thus
\begin{equation}\label{estimation1}
\|z^{\frac{\alpha +1}{2}}_n\|_{L^{2^{\ast }}(\Omega )}^{2}
\leq C_1\| u_n\|_{W_0^{1,2}(\Omega)}^{2}+C_3.
\end{equation}

Let us choose $u_n$ as a test function in \eqref{Pn-un}, we obtain
\begin{equation}\label{rachid}
\int_\Omega| \nabla u_n| ^{2}\,dx=\int_\Omega (z_n^p+f_n)u_n\,dx.
\end{equation}
It is clear that
$$
\int_\Omega f_nu_n\,dx\le C\|f\|_{L^2(\Omega)}\|u_n\|_{W_0^{1,2}(\Omega)}.
$$
Now, using H\"older inequality and taking into consideration the estimate 
\eqref{estimation1}, we obtain
\begin{align*}
\int_\Omega z_n^pu_n\,dx 
&\le  \Big(\int_\Omega z_n^{2^*\frac{\alpha +1}{2}}\,dx\Big)^{\frac{2p}{2^*(\alpha+1)}}
\Big(\int_\Omega u_n^{\frac{2^*(\alpha +1)}{2^*(\alpha+1)-2p}}\,dx
 \Big)^{\frac{2^*(\alpha+1)-2p}{2^*(\alpha+1)}}\\
&\le  C_2\Big(\int_\Omega |\nabla u_n|^2 \,dx\Big)^{\frac{p}{\alpha+1}} 
\Big(\int_\Omega u_n^{2^*}\,dx\Big)^{\frac{1}{2^*}},
\end{align*}
this is true because  $\frac{2^*(\alpha +1)}{2^*(\alpha+1)-2p}\le 2^*$. Hence
$$
  \int_\Omega z_n^pu_n\,dx
\le C\Big(\int_\Omega |\nabla u_n|^2 \,dx\Big)^{\frac{p}{\alpha+1}+\frac 12}+C_1(\Omega).
$$
Going back to \eqref{rachid} and taking into consideration that 
$\frac{p}{\alpha+1}+\frac 12<1$.
 We conclude that
 $$
\int_\Omega |\nabla u_n|^2 \,dx\le C \quad \text{for all } n .
$$
 Hence we obtain the existence of $u\in W_0^{1,2}(\Omega)$ such that, up to 
subsequences not relabeled, $u_n\rightharpoonup u$ weakly in $W_0^{1,2}(\Omega)$ 
and $u_n\to u$ strongly in $L^\sigma(\Omega)$ for all $\sigma<2^*$.

Now, by \eqref{ahmed} we conclude that 
$$
\|z_n^{\frac{\alpha +1}{2}}\|_{W_0^{1,2}(\Omega)}\leq C_1 \quad \text{for all } n.
$$ 
Hence we obtain the existence of a measurable function $z$ such that 
$z^{\frac{\alpha +1}{2}}\in W_0^{1,2}(\Omega)$ and, up to subsequences not relabeled,
 $z_n^{\frac{\alpha +1}{2}}\rightharpoonup z^{\frac{\alpha +1}{2}}$ 
weakly in $W_0^{1,2}(\Omega)$ and $z_n\to z$ strongly in $L^\sigma(\Omega)$ for all 
$\sigma<\frac{2^*(\alpha +1)}{2}$.

Since $p<\frac{\alpha+1}{2}$, then $\frac{ 2^*p}{2^*-1}<\frac{2^*(\alpha +1)}{2}$.
 Thus $z^{p}_n\to z^p$ strongly in $L^{\frac{2^*p}{2^*-1}}(\Omega)$. Hence classical 
results for elliptic problem allows us to conclude that
$$
u_n\to u \quad  \text{strongly in } W_0^{1,2}(\Omega).
$$
As a conclusion we obtain that $(u,z)$ is a solution to system \eqref{P} 
in the sense of Definition \ref{def:faible} with 
$(u,z^{\frac{\alpha+1}{2}})\in (W_0^{1,2}(\Omega))^2$.
\end{proof}

As a direct application of the Theorem \ref{existence}, we obtain the next
 existence result for the Bi-Laplacian problem with gradient term.

\begin{theorem}\label{existence00}
Let $\Omega \subset \mathbb{R}^N$ be a bounded domain. Suppose that $q<1$ and 
$g\in L^2(\Omega)$, then the problem
\begin{equation}\label{bii}
\begin{gathered}
\Delta^2 u  =  |\nabla u|^q+ g(x) \quad \text{in }\Omega , \\
\Delta u =  u=0 \quad \text{on }\partial \Omega,
\end{gathered}
\end{equation}
has a positive solution $u$ such that
$u\in W_0^{1,2}(\Omega)$ and $|\Delta u|^{\frac{\alpha+1}{2}}\in W_0^{1,2}(\Omega)$ where
 $\alpha$ satisfies $1<\frac{\alpha+1}{2}<\frac{1}{q}$.
\end{theorem}

\begin{proof}
Taking into consideration the result of Theorem \ref{existence} with 
$f\equiv 0$ and $p=1$, it follows that the system
\begin{equation}\label{P1009}
\begin{gathered}
-\Delta u= z \quad \text{in }\Omega , \\
-\Delta z= |\nabla u|^{q}+g(x)\quad \text{in }\Omega , \\
z=u=0 \quad \text{on }\partial \Omega ,
\end{gathered}
\end{equation}
has a solution $(u,z)$ with
$(u, z^{\frac{\alpha+1}{2}})\in W_0^{1,2}(\Omega)\\times W_0^{1,2}(\Omega)$ and 
$1<\frac{\alpha+1}{2}<\frac{1}{q}$.  Hence
$$
\Delta^2 u = |\nabla u|^q+ g(x) \text{ in }\Omega,
$$
and the result follows.
\end{proof}


\begin{remark}\rm 
(1) Following closely the above arguments, we can prove that the existence
 result holds for all $f\in L^{\frac{2^*}{2^*-1}}(\Omega)$ and 
$g\in L^{\frac{2^*}{2^*-(2-q)}}(\Omega)$.

(2) The same arguments can be used to treat the quasi-linear system
 \begin{equation}
 \begin{gathered}
 -\Delta_p u=v^r+f(x) \quad  \text{in }\Omega , \\
 -\Delta_p v=|\nabla u|^{q}+g(x) \quad  \text{in }\Omega , \\
 v,u>0 \quad  \text{in }\Omega , \\
 v=u=0 \quad  \text{on }\partial \Omega .
 \end{gathered}  \label{PP}
 \end{equation}
\end{remark}

In this case,  we have the next existence result.

\begin{theorem}   Assume that
$r>0, 0<q<p$ with $rq<(p-1)^2$ then for all $f,g \in L^{p'}(\Omega)$, 
then system \eqref{PP} has a positive solution $(u,v)$ such that
$(u,v^{\frac{\gamma +p-1}{p}})\in W^{1,p}_0(\Omega)\times W^{1,p}_0(\Omega) $ where 
$\gamma>0$ satisfies\, $r< \frac{p-1}{p}(p+\gamma-1)<\frac{(p-1)^2}{q}$.
 \end{theorem}

\section{Optimal results and open problems}\label{s4}

\subsection{Optimality of the obtained results}
 
\begin{theorem}\label{non0}
Assume that $N>4$ and that $q>\frac{2N}{N-2}=2^*$, then there exist 
$f, g \in L^2(\Omega)$ such that the system \eqref{P} has no positive solution.
\end{theorem}

\begin{proof}
We set $f(x)=\frac{1}{|x-x_0|^{2+\sigma}}$ where $x_0\in \Omega$ and $\sigma>0$ to be 
chosen later. Since $N>4$ and $q>\frac{2N}{N-2}$, then the interval 
$(\frac{N-q}{q}, \frac{N-4}{2})$ is not empty. Hence we choose 
$\sigma\in (\frac{N-q}{q}, \frac{N-4}{2})$. It is clear that $f\in L^2(\Omega)$. 
Now, we argue by contradiction. Assume that the system  \eqref{P} has a 
positive solution $(u,z)$ such that $|\nabla u|^q\in L^1(\Omega)$ and 
$(z^{p}+f)\in L^1(\Omega)$. Then $u\in W^{1,q}_0(\Omega)$. Recall that
$$
-\Delta u=z^p+ f\ge \frac{1}{|x-x_0|^{2+\sigma}}\quad \text{in  }\Omega.
$$
Using a simple comparison argument it holds that $u(x)\ge \frac{1}{|x-x_0|^{\sigma}}$ 
in a small ball $B_r(x_0)\subset \subset \Omega$. Since  $u\in W^{1,q}_0(\Omega)$,
 using Sobolev inequality we conclude that $u\in L^{q^*}(\Omega)$. 
Thus $u\in L^{q^*}(B_r(x_0))$. As a consequence we reach that 
$\frac{1}{|x-x_0|^{\sigma}}\in L^{q^*}(B_r(x_0))$. Hence $\sigma q^*<N$; 
which is a contradiction with the choice of $\sigma$.
\end{proof}

Let us begin by showing the optimality of the condition $pq<1$. 
More precisely we have the next non existence result.

 \begin{theorem}\label{non1}
 Assume that $q=2$, then for all $p>1$, there exist $f, g \in L^2(\Omega)$ such 
that the system \eqref{P} has no positive solution.
\end{theorem}

\begin{proof}
Without loss of generality we can assume that $f=\lambda f_1$ and $g=\mu g_1$
 with $f_1, g_1 \in L^\infty(\Omega)$. We argue by contradiction. 
Suppose that the system \eqref{P} has a positive solution $(u,z)$ such that 
$|\nabla u|^q\in L^1(\Omega)$ and $z^{p}\in L^1(\Omega)$. 
Let $\phi\in \mathcal{C}^\infty_0(\Omega)$, using $\phi^2$ as a test function 
in the equation of $u$ in the system \eqref{P}, it follows that
$$
\int_\Omega z^p\phi^2\,dx+\lambda\int_\Omega f_1\phi^2 \,dx=2\int_\Omega \phi\nabla \phi\nabla u \,dx.
$$
Now by Young inequality, it holds
\begin{equation}\label{ttt}
\int_\Omega z^p\phi^2 \,dx+\lambda\int_\Omega f_1\,\phi^2 \,dx
\le \int_\Omega \phi^2|\nabla u|^2 \,dx+\int_\Omega |\nabla \phi|^2\,dx.
\end{equation}
From the second equation in the system \eqref{P} we reach that 
$|\nabla u|^2\le -\Delta z$, thus,
$$
\int_\Omega \phi^2|\nabla u|^2 \,dx\le \int_\Omega \phi^2(-\Delta z)\,dx
=\int_\Omega z(-\Delta \phi)\,dx\le 2\int_\Omega z\phi(-\Delta\phi) \,dx,
$$
where the last estimate follows using Kato inequality.

Since $p>1$, using Young inequality, we conclude that
$$
\int_\Omega \phi^2|\nabla u|^2 \,dx\le \varepsilon\int_\Omega\phi^2 z^p \,dx 
+C(\varepsilon)\int_\Omega|\phi|^{\frac{p-2}{p-1}}|\Delta\phi|^{p'} \,dx.
$$
Choosing $\varepsilon$ small and going back to \eqref{ttt}, we obtain that
$$
\lambda\int_\Omega f_1\,\phi^2 \,dx
\le C(\varepsilon)\int_\Omega|\phi|^{\frac{p-2}{p-1}}|\Delta\phi|^{p'} \,dx+\int_\Omega |\nabla \phi|^2\,dx.
$$
Thus
$$
\lambda\le \frac{C(\varepsilon) \int_\Omega|\phi|^{\frac{p-2}{p-1}}|\Delta\phi|^{p'} \,dx
+\int_\Omega |\nabla \phi|^2\,dx}{ \int_\Omega f_1\,\phi^2 \,dx}.
$$
Setting
$$
M\equiv \inf_{\phi\in \mathcal{C}^\infty_0(\Omega)}
\frac{C(\varepsilon) \int_\Omega|\phi|^{\frac{p-2}{p-1}}|\Delta\phi|^{p'} \,dx
+\int_\Omega |\nabla \phi|^2\,dx}{ \int_\Omega f_1\,\phi^2 \,dx},
$$
then if $\lambda>M$, then system \eqref{P} has no positive solution and  we 
have the conclution.
\end{proof}

\subsection{Some open problems}

In this subsection we collect some interesting open problems.

(1) The case  $pq\ge 1$ and $q\le 2$: the arguments used to  treat the 
case $pq<1$ can not be adapted to the new situation  $pq\ge 1$ and $q\le 2$.  
Hence new arguments are needed to deal with this last case.

(2) If $p=1$, problem \eqref{P} takes the  form
 \begin{gather*}
 -\Delta u=z+f(x) \quad  \text{in }\Omega , \\
 -\Delta z=|\nabla u|^{q}+g(x) \quad  \text{in }\Omega , \\
 v,u>0 \quad \text{in }\Omega , \\
 v=u=0 \quad \text{on }\partial \Omega ,
 \end{gather*}
Now, by computing   $\Delta^2 u$, we reach that
\begin{equation}
 \begin{gathered}
 \Delta^2 u=|\nabla u|^{q}+\lambda h(x) \quad  \text{in }\Omega , \\
 u=\Delta u=0 \quad \text{on }\partial \Omega ,
 \end{gathered} \label{blap}
 \end{equation}
where $\lambda h=-\Delta f+g $. The existence of solution for \eqref{blap} 
is interesting for itself since, in the case where Bi-Laplacian operator 
is substituted by the Laplacian operator, an approach based on the classical 
elliptic capacity $\text{Cap}_{1,q'}$ gives a necessary and sufficient 
condition  to obtain the existence of a positive solution, 
see for instance the nice paper \cite{HMV}. For Bi-Laplacian operator, 
some particular cases were studied in \cite{Diez} with radial structure. 
It seems to be very interesting to get some similar approach in 
the case of Bi-Laplacian operator with gradient term if $q>1$.

\subsection*{Acknowledgements}
The authors would like to express their gratitude to the anonymous referees
 for their comments and suggestions that improve the last version of 
the manuscript.

\begin{thebibliography}{10}

\bibitem{AS} N. Alaa, M. Salim; 
\emph{Existence Result for Triangular Reaction-Diffusion Systems in $L^1$
 Data and Critical Growth with respect to the Gradient}.
Mediterr. J. Math.,  \textbf{10} (2013), 255-275

\bibitem{BP} P. Baras, M. Pierre; 
\emph{Singularit\'es \'eliminables pour des \'equations semi-lin\'eaires}.
Ann. Inst. Fourier, \textbf{34}, no. 1 (1984), 185-206.

\bibitem{BZ} H. Brezis, W. Strauss; 
\emph{Semi-linear second-order elliptic equations in $L^1$}. J. Math. Soc. Japan,
\textbf{25} (1973), no. 4, 565-590.

\bibitem{Bo} L. Boccardo; 
\emph{Dirichlet problems with singular and gradient quadratic lower order terms}. 
ESAIM - Control, Optimisation and Calculus of Variations, \textbf{14}, no. 3 (2008),
 411-426.

\bibitem{BOP} L. Boccardo, L. Orsina, A. Porretta; 
\emph{Existence of Finite Energy Solutions For Elliptic Systems With $L^1$-Value 
Nonlinearities}. Mathematical Models and Methods in Applied Sciences \textbf{ 18}, 
No. 5 (2008), 669-687.

\bibitem{BOPu} L. Boccardo, L. Orsina, I. J. Puel; 
\emph{A quasilinear elliptic system with natural growth terms}. 
Annali di Matematica,  \textbf{194}, no. 3 (2015), 1733-1750.

\bibitem{Diez} J. I. Diaz, M. Lazzo, P. G. Schmidt; 
\emph{Large Solutions for a System of elliptic equation arising from fluid dynamics}. 
Siam Journal on Mathematical Analysis, \textbf{37} (2005), 490-513.

\bibitem{Diez2} J. I. Diaz, J.-M. Rakotoson, P. G. Schmidt; 
\emph{A parabolic system involving a quadratic gradient term
related to the Boussinesq approximation}. Rev. R. Acad. Cien. Serie A. Mat. 
\textbf{ 101}(1), 2007, pp. 113–118

\bibitem{Diez3} J. I. Diaz, J. M. Rakotoson, P. G. Schmidt; 
\emph{Local strong solutions of a parabolic system related to the Boussinesq 
approximation for buoyancy-driven flow with viscous heating}. 
Adv. Differential Equations, \textbf{13},  no. 9-10 (2008), 977–1000.

\bibitem{C} S. Clain, J. Rappaz, M. Swierkosz, R. Touzani; 
\emph{Numerical modeling of induction heating for two dimentional geometrie},
 Math. Models Methods Appl. Sci., \textbf{3}, no. 6 (1993),  805-822.

\bibitem{F} J. R. Ferguson, J. M. Fiard, R. Herbin; 
\emph{a two dimensional simulation of a solid oxide fuel cell}. 
International energy agency Worshop: Fundamental barries of SOFC performence,
 Lausane, Switzerland, Augest 1992.

\bibitem{GH} T. Gallo\"uet, R. Herbin;
\emph{Existence of solution to a coupled elliptic system},
 Applied. Math. Lett., \textbf{7}, no 2 (1994), 49-55.

\bibitem{Gur} G. Singh;
\emph{Classification of radial solutions for semilinear elliptic systems
with nonlinear gradient terms Gurpreet Singh},
 Nonlinear Analysis, \textbf{129} (2015), 77–103.

\bibitem{HMV} K. Hansson, V. G. Maz'ya, I. E. Verbitsky; 
\emph{Criteria of solvability for multidimensional
Riccati equations}. Ark. Mat., \textbf{37} (1999), 87-120.

 \end{thebibliography}


\end{document}