  We look for least energy solutions to the cooperative systems of coupled
Schr\"odinger equations
  \begin{equation*}
  \begin{cases}
  -\Delta u_i + \lambda_i u_i = \partial_iG(u)\quad \mathrm{in} \ \mathbb{R}^N,
\ N \geq 3,
  u_i \in H^1(\mathbb{R}^N),
  \int_{\mathbb{R}^N} |u_i|^2 \, dx \leq \rho_i^2
  \end{cases}
  i\in\{1,\dots,K\}
  \end{equation*} with $G\geq 0$, where $\rho_i>0$ is prescribed and
$(\lambda_i, u_i) \in \mathbb{R} \times H^1 (\mathbb{R}^N)$ is to be
determined, $i\in\{1,\dots,K\}$.
  Our approach is based on the minimization of the energy functional over a
linear combination of the Nehari and Poho\v{z}aev constraints intersected with
the product of the closed balls in $L^2(\mathbb{R}^N)$ of radii $\rho_i$, which
allows to provide general growth assumptions about $G$ and to know in advance
the sign of the corresponding Lagrange multipliers. We assume that $G$ has at
least $L^2$-critical growth at $0$ and admits Sobolev critical growth. The more
assumptions we make about $G$, $N$, and $K$, the more can be said about the
minimizers of the corresponding energy functional. In particular, if $K=2$,
$N\in\{3,4\}$, and $G$ satisfies further assumptions, then $u=(u_1,u_2)$ is
normalized, i.e., $\int_{\mathbb{R}^N} |u_i|^2 \, dx=\rho_i^2$ for
$i\in\{1,2\}$.
