We study the boundary value problem $\begin{array}{rcl} {\rm div}(|\n u|^{m-2}\n u) + u^av^b &=& 0\quad \mbox{ in } {\Omega}, \vspace{\jot}\\ {\rm div}(|\n v|^{m-2}\n v) + u^cv^d &=& 0 \quad \mbox{ in } {\Omega}, \vspace{\jot}\\ u =v &= & 0 \quad \mbox{ on } {\partial}{\Omega},\vspace{\jot}\\ \end{array}$ where ${\Omega}\subset\mathbb{R}^n$ ($n\ge2$) is a bounded connected smooth domain, and the exponents $m>1$ and $a,b,c,d\ge0$ are non-negative numbers. Under appropriate conditions on the exponents $m$, $a$, $b$, $c$ and $d$, a variety of results on a priori estimates and existence of positive solutions has been established.