Pure Neumann problem for slightly supercritical Lane-Emden system

We are concerned with the following Lane-Emden system with Neumann boundary conditions: {-Delta u(1)=|u(2)|(p & varepsilon;-1)u(2), in Omega, -Delta u(2)=|u(1)|q(& varepsilon;-1)u(1), in Omega, partial derivative(nu)u(1)=partial derivative(nu)u(2)=0, on partial derivative Omega where Omega=B-1(0) is the unit ball in R-n (n >= 4) centered at the origin, p(& varepsilon;)=p+alpha & varepsilon;,q(& varepsilon;)=q+beta & varepsilon; with alpha,beta>0 and 1/p+1+1/q+1=n-2/n. We construct solutions based on the Lyapunov-Schmidt reduction argument incorporating the zero-average condition by certain symmetries. It's important to emphasize that we are simultaneously examining two scenarios: when p is greater than n/n-2 and when p is less than n/n-2. The coupling mechanisms within the Lane-Emden system differ significantly between these two cases, resulting in notable variations in solution behaviors. This distinction also leads to varying research challenges. Presently, there is few literature that takes both of these ranges into consideration when addressing solution construction. Consequently, this aspect constitutes the main feature and new ingredient of our work.