We study positive solutions to the steady-state reaction diffusion systems of the form: \(\begin{aligned} \left\{ \begin{array}{ll} -\Delta u = \lambda f(v)+\mu h(u), & \Omega ,\\ -\Delta v = \lambda g(u)+\mu q(v),& \Omega ,\\ \frac{\partial u}{\partial \eta }+\root \of {\lambda +\mu }\, u=0,& \partial \Omega ,\\ \frac{\partial v}{\partial \eta }+\root \of {\lambda +\mu }\, v=0, & \partial \Omega ,\\ \end{array}\right. \end{aligned}\) where \({\lambda ,\mu }\) are positive parameters, \({\Omega }\) is a bounded domain in \({\mathbb {R}}^{N}(N>1)\) with smooth boundary \({\partial \Omega }\) or \({\Omega =(0,1)}\) , \({ \frac{\partial z}{\partial \eta } }\) is the outward normal derivative of z. Here, we assume that f, g, h, and q are increasing continuous functions such that \(f(0) = g(0) = h(0) = {q}(0) = 0\) and \(f^\prime (0), g^\prime (0), h^\prime (0), q^\prime (0) > 0\) . We further assume that f and g are combined sublinear at infinity (i.e., \(\lim \limits _{s\rightarrow \infty }\frac{f(M g(s))}{s}=0\) for all \(M>0\) ). Under certain additional assumptions on f, g, h, and q, we establish the existence and multiplicity results for the above system. Our existence and multiplicity results are proved using sub-super solution methods.
