I recently read in a paper that Monte Carlo integration is very slow near phase transitions. I couldn't understand the concept properly. The paper is Phys. Rev. B 95, 035105 (2017).
Near a first-order phase transition, converting one phase into the other typically involves crossing a free energy barrier. Even though the two phases may have the same free energy, the intermediate configurations will involve both phases present at the same time, which implies an interface between them. The contribution of the interfacial free energy will be significant, in general: proportional to the surface tension and the area of the interface (which approaches $\sim L^2$ in a simulation in a box of length $L$). So, crossing this barrier can be an activated process, with a rate $\propto \exp(-\Delta E^\dagger/kT)$, $T$ being the temperature and $k$ Boltzmann's constant, with $\Delta E^\dagger\propto L^2$. This will happen infrequently for large $L$.
Near a continuous phase transition one observes critical slowing down. Long-wavelength fluctuations in the appropriate order parameter (for example, the density in the case of the liquid-vapour critical point) occur. In real experiments, these can be observed, since they scatter light and the system becomes turbid. Very close to the phase transition, in a finite-sized system of length $L$, the wave length of the fluctuations approaches $L$. There is generally a relation between the time scale and the length scale involving a dynamical critical exponent; the bottom line is that the fluctuations become very slow.
In basic Monte Carlo simulations, the elementary moves consist of displacing a single particle or flipping a single spin. So, although they don't exactly follow the dynamics of a real system, they are local in nature: a move doesn't magically overcome the physical timescale challenges described above. However (as described in several references in that paper) attempts have been made to devise smarter Monte Carlo algorithms which do move many atoms, or flip many spins, in a correlated fashion. If you are interested, you should look at "cluster updating" or "Swendsen-Wang", which can be very successful near continuous transitions. Essentially, a cluster of particles is defined, by a probabilistic method, in such a way that the movement or flipping of the cluster as a single entity will be accepted with reasonably high probability, while having a dramatic non-local effect on the configuration of atoms or spins. The moves are still intended to sample the Boltzmann distribution: they are constructed in a way that satisfies detailed balance. However, such approaches don't work in every situation.
I'll also mention in passing that special "barrier-flattening" Monte Carlo methods have been proposed to attempt to tackle the problem of free energy barriers which occur in first-order phase transitions. The Wang-Landau method is the best known, but there are many variants. A common feature of these is that the sampling is not from the Boltzmann distribution, but from one in which the free energy barrier has effectively been cancelled out; then the desired averages over the Boltzmann distribution are reconstructed afterwards. In this case the key point is the modification of the sampled distribution, rather than the devising of smarter moves.
The paper that you cite seems to be an attempt to systematically approach this problem and devise smarter moves with the aid of machine learning.