Arbitrary number of thermally induced phase transitions in different universality classes in $XY$ models with higher-order terms

We propose generalized variants of the $XY$ model capable of exhibiting an arbitrary number of phase transitions only by varying temperature. They are constructed by supplementing the magnetic coupling with $n_t-1$ nematic terms of exponentially increasing order with the base $q=2,3,4$ and $5$, and increasing interaction strength. It is found that for $q=2,3$ and $4$ with sufficiently large coupling strength of the final term, the models exhibit a number of phase transitions equal to the number of the terms in the generalized Hamiltonian. Starting from the paramagnetic phase, the system transitions through the cascade of $n_t-1$ nematic phases of the orders $q^{k}$, $k=n_t-1,n_t-2,\hdots,1$, that are characterized by $q^{k}$ preferential spin directions symmetrically disposed around the circle, to the ferromagnetic (FM) phase at the lowest temperatures. Besides the BKT transition from the paramagnetic phase, all the remaining transitions have a non-BKT nature: depending on the value of $q$ they belong to either the Ising ($q=2$ and $4$) or the three-states Potts ($q=3$) universality class. For $q=5$, due to the interplay between different terms, the phase transitions between the ordered phases observed for $q<5$ split into two and the number of the ordered phases increases to $2n_t-1$. These phases are characterized by a domain structure with the gradually increasing short-range FM ordering within domains that extends to different kinds of FM ordering in the last two low-temperature phases. The respective transitions do not seem to obey any universality.