In quantum Brownian motion, bosonic/harmonic oscillator bath and interaction described by Hamiltonian
$$ H_B = \sum_{n}\hbar\omega_n(b_n^\dagger b_n) \\ H_I = -\sigma_x \otimes B $$
and
$$ B = \sum_n \kappa_n\sqrt{\frac{\hbar}{2m_n\omega_n}}(b_n+b_n^\dagger)$$
We have two correlation function to describe time evolution of system (See Breuer p.174)
$$ D(\tau) \equiv iTr_R([B,B(-\tau)]\rho_B) \\ D_1(\tau) \equiv Tr_R(\{B,B(-\tau)\}\rho_B) $$
And using Spectral Density (approximating infinite mode) we get (setting $\hbar=1$)
$$ D(\tau)=2\int_0^\infty d\omega \, J(\omega)\sin(\omega\tau) \\ D_1(\tau)=2\int_0^\infty d\omega \, J(\omega)\coth \left(\frac{\hbar\omega}{2k_BT} \right)\cos(\omega\tau) $$
In many references (including Breuer) ohmic spectral density with Drude-Lorentz cutoff is used $$ J(\omega)=\alpha^2\omega\frac{\lambda^2}{\omega^2+\lambda^2}$$
And correlation function can be easily analytically calculated.
My question is, if we define more general spectral density
$$J(\omega)=\alpha^2\omega^s\omega_{ph}^{1-s}\frac{\lambda^2}{\omega^2+\lambda^2}$$
With $S = 1$ (Ohmic), $S > 1$ (Super-Ohmic), $S < 1$ (Sub-Ohmic). In Sub-ohmic region, working the integral for $D(\tau)$ in Mathematica gives Hypergeometric function which is more complicated than Ohmic case but still a well behaved function. The problem arise when choosing spectral density as a Super-Ohmic
$$ D(\tau)=2\alpha^2\lambda^2\omega_ph^{1-s}\int_0^\infty \frac{\omega^s \sin(\omega\tau)}{\omega^2+\lambda^2}$$
For $S > 1$ the integral is divergent. From this result, is the boson bath "limited" to the ohmic and sub-ohmic spectral density(is there any microscopical model to suggest this) ? How the divergence of correlation function can be interpret in relation to super-ohmic ?