Entanglement entropy dynamics of non-Gaussian states in free boson systems: random sampling approach

We develop a random sampling method for calculating the time evolution of the Rényi entanglement entropy after a quantum quench from an insulating state in free boson systems. Because of the non-Gaussian nature of the initial state, calculating the Rényi entanglement entropy calls for the exponential cost of computing a matrix permanent. We numerically demonstrate that a simple random sampling method reduces the computational cost of a permanent; for an $N_{\mathrm{s}}\times N_{\mathrm{s}}$ matrix corresponding to $N_{\mathrm{s}}$ sites at half filling, the sampling cost becomes $\mathcal{O}(2^{\alpha N_{\mathrm{s}}})$ with a constant $\alpha\ll 1$, in contrast to the conventional algorithm with the $\mathcal{O}(2^{N_{\mathrm{s}}})$ number of summations requiring the exponential-time cost. Although the computational cost is still exponential, this improvement allows us to obtain the entanglement entropy dynamics in free boson systems for more than $100$ sites. We present several examples of the entanglement entropy dynamics in low-dimensional free boson systems.