Polarization response including ergodicity breaking and the divergence of relaxation time is reproduced for model Hamiltonians of growing complexity. Systematic derivation of the dynamical equations and its solutions is based on the Fokker-Planck and imaginary time Schrodinger equation techniques with subsequent symplectic integration. Test solutions are addressed to finite size and spatially extended problems with microscopically interpretation of the model parameters as a challenge.