  We consider the pricing of VIX options in the rough Bergomi model. In this
setting, the VIX random variable is defined by the one-dimensional integral of
the exponential of a Gaussian process with correlated increments, hence
approximate samples of the VIX can be constructed via discretization of the
integral and simulation of a correlated Gaussian vector. A Monte-Carlo
estimator of VIX options based on a rectangle discretization scheme and exact
Gaussian sampling via the Cholesky method has a computational complexity of
order $\mathcal{O}(\varepsilon^{-4})$ when the mean-squared error is set to
$\varepsilon^2$. We demonstrate that this cost can be reduced to
$\mathcal{O}(\varepsilon^{-2} \log^2(\varepsilon))$ combining the scheme above
with the multilevel method, and further reduced to the asymptotically optimal
cost $\mathcal{O}(\varepsilon^{-2})$ when using a trapezoidal discretization.
We provide numerical experiments highlighting the efficiency of the multilevel
approach in the pricing of VIX options in such a rough forward variance
setting.
