Why quantum states are classified using only conserved quantities? While studying quantum mechanics from standard textbooks I always felt some conceptual gap that was never mentioned or explained. In what follow I tried to formulate my question, please be patient with me.
For a quantum particle in an infinite potential well the stationary states are labelled by the quantum number $n$ which labels the eigenenergies. An eigenenergy, that corresponds to a stationary state, does not change with time, hence is a conserved quantity.
For a spinless electron in Coulomb potential, to model the hydrogen atom, again we have the same story, the stationary states are labeled by the quantum numbers $n$, $l$, $m$ which corresponds to conserved quantities.
My question is rather general since I am trying to understand conceptually why only conserved quantities are used to label the quantum states.
I mean how would someone think in advance that he has to look for conserved quantities, and then use such conserved quantities to label the states ? 
 A: It is not just a matter of symmetries, though the latter play an important role.
The rationale for the labeling is that one wants to have a simple basis for the Hilbet space of wavefunctions in which the energy is diagonal. Such a basis has two advantages: (i) The time-dependent solution is easy to compute. (ii) One has a good starting point for perturbation theory if one adds additional interactions. Indeed, the labeling is usually done for a highly symmetric exactly solvable prototype Hamiltonian, and realistic systems are obtained by adding appropriate interactions.
That the energy is diagonal implies that one of the labels must be related to the energy. (It is often not the energy itself, but some quantity simply related to it). All other labels serve to distinguish states with the same energy. These form a small Hilbert space; it is 1-dimensional when the eigenvalue is nondegenerate, but if the system has a nontrivial symmetry (an operator commuting with the Hamiltonian) then the eigenspaces are usually higher-dimensional. 
The main point now is that by the spectral theorem, any maximal set of commuting selfadjoint operators has 1-dimensional common eigenspaces, hence their joint eigenvectors form a basis completely labelled with the corresponding eigenvalues. Therefore the goal must be to extend the Hamiltonian $H$ to such a complete set by adding enough additional operators $Q$. of course, any such operator must commute with $H$ and with all operators already in the set. But $[H,Q]=0$ implies that $Q$ is a conserved quantity. Hence one has to look for conserved quantities.
Now if a system has a symmetry group preserving the Hamiltonians, its generators commute with the Hamiltonian and are primary candidates for choice.
But often they do not commute with each other, so that one has to be somewhat selective. For example when the group is $SO(3)$, one generally uses one of the generators, typically $J_3$, but this does not yet give a complete set of labels. One also needs to add Casimir operators of the Lie algebra (i.e., central elements in the universal enveloping algebra, in case of $SO(3)$ it is $J^2$. For small groups, this sort of heuristics is sufficient to get a complete set.
For bigger symmetry groups, the problem of finding a complete set is more involved, amounting to find Casimirs of a chain of subgroups of the group. 
If this is not enough, one speaks about the problem of ''missing labels''.
A search for this keyword (together woth Lie algebras) in http://scholar.google.com turns up relevant papers for specific instances.
A: I believe the whole point is in symmetries; for every symmetry of a Hamiltonian at hand you identify, you can construct an operator which commutes with H. When you collect a complete set of operators, then it can be proven that every eigenstate of H can be uniquely labeled by set of eigenvalues of all operators from complete set of your operators. Therefore, it is much easier to label wavefunctions through small set of real eigenvalues, than writing down infinite numbers of projections of a L^2 set member. 
A: This is an accident of the fact that when you are given a quantum system and its solution, the classical limit must have been integrable, because you could solve it. In the integrable case, there are action variables $J_n$ which fill up the whole phase space with invariant tori, defined by fixing the action variables, and letting the canonically conjugate angle variables vary over the full range The action variables $J_n$'s are conserved quantities by definition. Integrability means that any point in the space is labeled by the $J$'s and the conjugate $\theta$'s.
In the semiclassical limit of large number labels for the energy levels, the energy is given by setting the angle variables to integers times $2\pi$ (in units where $\hbar=1$)
$$ J_k = 2\pi n_k $$
and near large J_k, 
$$ {\partial H\over \partial J_k} = {1 \over T_k} $$
Where $T_k$ is the classical period of the action variable. This makes the energy levels spacing a multiperiodic projected lattice near any one energy $E_0$
$$ E = E_0 + \sum_k {2\pi \Delta n_k\over T_k} = \sum_k \omega_k \Delta n_k $$
You can extract the classical periods of the angle variables from the energy level diagram. The energy level diagram near any (large) energy will look periodic with multiple incommensurate periods.
If you want a quantum system where the energy states are not labelled by classical conserved quantities, at least not by useful ones, take any classically chaotic system and consider the quantum version. You need at least two degrees of freedom, so the simplest example is the quantum double pendulum, where there are two angles $\theta$ and $\phi$, and the Lagrangian is given by
$$L = {1\over 2} (\dot{\theta}^2 + \dot{\phi}^2) - \cos(\theta) - A cos(\phi)$$
If you choose an A so that the system is classically chaotic, the eignestates of this are going to have discrete energy levels, but the labels cannot correspond to action variables, because there are none. The energy levels in a chaotic quantum system obey random matrix statistics in the semiclassical limit, they do not form multiperiodic sequences locally, like the integrable case.
A: Also, another way of phrasing Stipe Galic's answer is:  even if you choose to adopt the Schroedinger picture, in which the states are changing, provided you've labelled your states with the values of conserved quantities, then the labels are maintained through the time evolution.
A: A point that is being missed here completely is that an eigenstate (a state described by conserved quantities) won't mix when you evolve it in time. It will disperse (broaden) spatially, but it's quantum numbers will remain the same (in the absence of an external potential).
A: Are you familiar with Rydberg formula (1888 year)?
http://en.wikipedia.org/wiki/Rydberg_formula
Which is extension from older Balmer formula (1885 year)?
http://en.wikipedia.org/wiki/Balmer_series
They both have $n$ parameter, full integer positive.
Data for Balmer series are gathered while ionization of Hydrogen gas by high voltage. Equation is describing observation made in laboratory.
Schrödinger equation is just extension made 36 years later.
