  The systematic approach to study bound states in gluodynamics is presented.
The method utilizes flow equations together with low-energy phenomenology, that
provides the perturbative renormalization scaling in conjuction with the change
of the basis to constituent gluon states. The renormalized effective
Hamiltonian of gluodynamics up to the second order is obtained at low energies,
which provides a kind of constituent gluon model for glueball bound states. The
approach allows to include perturbative QCD corrections into nonperturbative
calculations of many-body techniques. The performed numerical calculations
support the constituent picture of hadronic observables.
