Summary: In this paper, we investigate the following nonlinear and non-homogeneous elliptic system \[ \begin{cases} -\operatorname{div}(a_1(|\nabla u|) \nabla u) = \lambda_1 F_u(x,u,v) - \lambda_2 G_u(x,u,v) - \lambda_3 H_u(x,u,v)\quad &\text{in}\; \Omega, \\ -\operatorname{div}(a_2(|\nabla v|) \nabla v) = \lambda_1 F_v(x,u,v) - \lambda_2 G_v(x,u,v) - \lambda_3 H_v(x,u,v) &\text{in}\; \Omega, \\ u = v = 0 &\text{on}\; \partial\Omega, \end{cases} \] where \(\Omega\) is a bounded domain in \(\mathbb{R}^N(N \geq 1)\) with smooth boundary \(\partial \Omega\), \(\lambda_1\), \(\lambda_2\), \(\lambda_3\) are three parameters, \(\phi_i(t) = a_i(|t|)t\) (\(i = 1,2\)) are two increasing homeomorphisms from \(\mathbb{R}\) onto \(\mathbb{R}\), and functions \(F\), \(G\), \(H\) are of class \(C^1(\Omega \times \mathbb{R}^2, \mathbb{R})\) and satisfy some reasonable growth conditions. By using a three critical points theorem due to B. Ricceri, we obtain that system has at least three solutions. With some additional conditions, by using a four critical points theorem due to G. Anello, we obtain that system has at least four solutions.