We study the $L^2$-gradient flow of functionals $\mathcal F$ depending on the eigenvalues of Schrödinger potentials $V$ for a wide class of differential operators associated to closed, symmetric, and coercive bilinear forms, including the case of all the Dirichlet forms (as for second order elliptic operators in Euclidean domains or Riemannian manifolds). We suppose that $\mathcal F$ arises as the sum of a $-θ$-convex functional $\mathcal K$ with proper domain $\mathbb{K}\subset L^2$ forcing the admissible potentials to stay above a constant $V_{\rm min}$ and a term $\mathcal H(V)=φ(λ_1(V),\cdots,λ_J(V))$ which depends on the first $J$ eigenvalues associated to $V$ through a $C^1$ function $φ$. Even if $\mathcal H$ is not a smooth perturbation of a convex functional (and it is in fact concave in simple important cases as the sum of the first $J$ eigenvalues) and we do not assume any compactness of the sublevels of $\mathcal K$, we prove the convergence of the Minimizing Movement method to a solution $V\in H^1(0,T;L^2)$ of the differential inclusion $V'(t)\in -\partial_L^-\mathcal F(V(t))$, which under suitable compatibility conditions on $φ$ can be written as \[ V'(t)+\sum_{i=1}^J\partial_iφ(λ_1(V(t)),\dots, λ_J(V(t)))u_i^2(t)\in -\partial_F^-\mathcal K(V(t)) \] where $(u_1(t),\dots, u_J(t))$ is an orthonormal system of eigenfunctions associated to the eigenvalues $(λ_1(V(t)), ,\dots,λ_J(V(t)))$.