Solar activity has a cyclic nature with the ≈11-year Schwabe cycle dominating its variability on the interannual timescale. However, solar cycles are significantly modulated in length, shape and magnitude, from near-spotless grand minima to very active grand maxima. The ≈400year-long direct sunspot-number series is inhomogeneous in quality and too short to study robust parameters of long-term solar variability. The cosmogenic-isotope proxy extends the timescale to twelve millennia and provides crucial observational constraints of the long-term solar dynamo modulation. Here, we present a brief up-to-date overview of the long-term variability of solar activity at centennial -millennial timescales. The occurrence of grand minima and maxima is discussed as well as the existing quasi-periodicities such as centennial Gleissberg, 201-year Suess/de Vries and 2400-year Hallstatt cycles. It is shown that the solar cycles contain an important random component and have no clock-like phase locking implying a lack of long-term memory. A brief yet comprehensive review of the theoretical perspectives to explain the observed features in the framework of the dynamo models is presented, including the nonlinearity and stochastic fluctuations in the dynamo. We keep gaining knowledge of the processes driving solar variability with the new data acquainted and new models developed.