This article deals with a two-predator and one-prey system with prey-taxis \(\begin{aligned} \left\{ \begin{array}{l} \begin{aligned} & u_t =d_1\Delta u -\chi _1\nabla \cdot (u\nabla w)+\gamma _1uF(w)-uh(u)-\beta _1uv,& x\in \Omega ,\ t>0& ,\\ & v_t =d_2\Delta v -\chi _2\nabla \cdot (v\nabla w)+\gamma _2vF(w)-vg(v)-\beta _2uv,& x\in \Omega ,\ t>0& ,\\ & w_t = d_3\Delta w-(u+v)F(w)+f(w),& x\in \Omega ,\ t>0& \\ \end{aligned} \end{array} \right. \end{aligned}\) in a bounded and smooth domain \(\Omega \subset \mathbb R^n\) \((n\ge 2)\) under homogeneous Neumann boundary conditions. In the two-dimensional case, we confirm the global boundedness of the classical solution for all suitably regular initial data. Moreover, in the higher-dimensional counterpart, the uniform boundedness of the global classical solution is also established provided that the sensitivity of the functional response function F(w) at high prey densities is sufficiently large. Finally, we study the large-time behavior of the prey-only, semi-coexistence and coexistence steady states by constructing several suitable Lyapunov functionals.
