  We introduce Rayleigh functional for nonlinear systems. It is defined using
the energy functional and the normalization properties of the variables of
variation. The key property of the Rayleigh quotient for linear systems is
preserved in our definition: the extremals of the Rayleigh functional coincide
with the stationary solutions of the Euler-Lagrange equation. Moreover, the
second variation of the Rayleigh functional defines stability of the solution.
This gives rise to a powerful numerical optimization method in the search for
the energy minimizers. It is shown that the well-known imaginary time
relaxation is a special case of our method. To illustrate the method we find
the stationary states of Bose-Einstein condensates in various geometries.
Finally, we show that the Rayleigh functional also provides a simple way to
derive analytical identities satisfied by the stationary solutions of the
critical nonlinear equations. We also show that the functional used in Phys.
Rev. A 66, 036612 (2002) is erroneous.
