Finite-Size Thermodynamics of the Two-Dimensional Dipolar Q-Clock Model

1. Introduction

Low-dimensional magnetic systems with long-range interactions provide a fertile setting for exploring the interplay between geometry, anisotropy, and cooperative effects in mesoscopic materials. Among the paradigmatic models used to describe such systems, the Q-state clock model has emerged as a minimal yet versatile framework for two-dimensional spins with discrete rotational symmetry, interpolating smoothly between the Ising ($Q=2$) [1,2] and XY [3,4,5,6,7] ($Q\to \infty$) limits. Even in the absence of long-range couplings, its thermodynamic behavior is remarkably rich, exhibiting discrete-symmetry ordering for all Q and Berezinskii–Kosterlitz–Thouless (BKT) physics for $Q\ge 4$ [8].

The inclusion of dipolar interactions adds a further layer of complexity. Dipolar couplings decay slowly ($\propto r^{-3}$), placing them within the broader class of long-range interactions characterized by nontrivial critical behavior [9,10,11]. Their strong anisotropy promotes frustrated magnetic textures and competing domain structures, playing a central role in modern two-dimensional magnets such as van der Waals materials (e.g., CrI₃, CrGeTe₃, Fe₃GeTe₂) [12,13], as well as in artificial spin-ice arrays and dipolar-coupled nanomagnets, where stripe phases, domain patterns, and topological defects emerge from the competition between short- and long-range interactions [14,15,16,17].

Despite this relevance, the combined effects of discrete clock symmetry, long-range dipolar anisotropy, and finite cluster geometry remain comparatively unexplored, particularly in small systems where the microscopic energy spectrum can be resolved exactly. Finite clusters are of special interest for two complementary reasons. First, finite-size effects and boundary conditions strongly modify thermodynamic responses: they round critical singularities, shift specific-heat maxima, suppress vortex proliferation, and enhance degeneracy-driven features absent in the thermodynamic limit. Second, for sufficiently small lattices, the full spectrum of the Hamiltonian can be obtained exactly, enabling a noise-free evaluation of thermodynamic observables and a direct correspondence between microscopic level crossings and macroscopic response functions [18,19,20].

An important novelty of the present work lies in its systematic application of an exact enumeration method to a dipolar clock model in a regime where long-range anisotropic interactions, discrete symmetry, and finite-size effects act on equal footing. In particular, we show that dipolar interactions alone suffice to induce sharp ground-state rearrangements in clusters as small as $3\times 3$, which manifest as zeros of the specific heat in the limit $T\to 0$. These zero-temperature anomalies provide a direct thermodynamic fingerprint of level crossings in the microscopic spectrum, allowing us to connect caloric features to changes in ground-state degeneracy unambiguously.

Moreover, by comparing different clock symmetries and lattice sizes, we show that the symmetry or asymmetry of the resulting Schottky-like anomalies encodes detailed information about the multiplicity structure of the low-energy levels. In particular, we identify lattice parity as a key organizing principle that controls whether the post-critical caloric response is symmetric or strongly skewed. This establishes a direct link between finite-size geometry, degeneracy patterns, and measurable thermodynamic signatures, which previous studies of dipolar clocks or Ising models did not explicitly address.

In this work, we investigate the dipolar Q-state clock model on small square lattices with free boundaries, combining exact enumeration of all spin configurations with analytical evaluation of canonical thermodynamic observables. We resolve the full sets of energies and degeneracies $\{E_n,c_n\}$ for clock symmetries $Q=2,4,6,8$ on a $3\times 3$ lattice, as well as for the Ising case ($Q=2$) on lattices of size $L=3,4,5$. By varying the dipolar ratio $\alpha =D/J$, we track how the competition between exchange and dipolar interactions reshapes the low-energy spectrum and manifests in the specific heat, entropy, and ground-state multiplicity. This controlled finite-size setting allows us to identify robust signatures of dipolar-induced spectral reorganizations across lattice sizes and discrete symmetries.

Our analysis reveals several key phenomena. For $Q=2$ and $Q=4$, the system exhibits ground-state crossings at well-defined values of $\alpha$, which appear as zeros of the specific heat in the limit $C(T\to 0,\alpha )$. The Schottky-like anomalies that emerge around these points encode the degeneracy structure of the low-lying spectrum, with odd lattices ($L=3,5$) displaying strongly asymmetric peaks due to unbalanced multiplicities, and the even lattice ($L=4$) showing nearly symmetric anomalies associated with paired degeneracies. In contrast, for $Q\ge 6$ the energy landscape becomes progressively smoother, ground-state crossings disappear within the explored range of $\alpha$, and the specific heat exhibits only broad, purely thermal maxima, reflecting the increasing angular continuity of the spin space.

Beyond the strictly low-temperature regime, the Ising model ($Q=2$) develops a double-peak structure in the specific heat over a finite interval of $\alpha$, signaling the coexistence of bulk-like thermal excitations and additional low-energy modes generated by the competition between exchange and dipolar interactions. For higher clock symmetries ($Q\ge 6$), increasing $\alpha$ suppresses this double-peak structure as the dipolar interaction reorganizes the low-energy spectrum and drives the ground state away from an exchange-stabilized ferromagnetic configuration toward textures governed primarily by long-range anisotropy.

Although the finite-size systems considered here do not host a genuine Berezinskii–Kosterlitz–Thouless phase transition, the caloric features observed for $Q\ge 6$ bear a qualitative resemblance to the crossover phenomena associated with vortex–antivortex physics in the two-dimensional XY model [4,5]. In this sense, dipolar interactions effectively smooth the discrete angular structure of the clock model, promoting gradual thermal activation rather than sharp degeneracy-driven anomalies.

Taken together, these results demonstrate that dipolar interactions alone are sufficient to generate nontrivial caloric behavior, including critical-like features associated with ground-state level crossings, in clusters as small as $3\times 3$. The present study therefore provides a controlled finite-size benchmark for interpreting experiments on van der Waals nanomagnets, patterned dipolar arrays, and mesoscopic magnetic textures, where discrete symmetry, long-range anisotropy, and cluster geometry jointly determine the thermodynamic response.

The remainder of this article is organized as follows. Section 2 introduces the Hamiltonian, discusses the role of dipolar anisotropy, and outlines the exact enumeration approach. Section 3 presents the thermodynamic properties of the model, including internal energy, entropy, specific heat, ground-state multiplicity, and the effects of lattice size, parity, and clock symmetry. Section 4 summarizes the main findings and discusses perspectives for extending this exact approach to larger clusters and dipolar-coupled nanomagnetic systems.

2. Model and Methods

2.1. General Definitions

We consider the two-dimensional Q-state clock model on a square lattice of size $L\times L$ (number of sites $N=L^2$) with free boundaries, where each local magnetic moment (“spin”) $S_i$ at site i points along one of Q equally spaced directions in the plane. Writing $S_i=(sin{\theta }_i,cos{\theta }_i)$ with ${\theta }_i=2\pi k/Q$, $k=0,1,\dots ,Q-1$, the model provides a natural interpolation between the Ising limit ($Q=2$) and the XY model ($Q\to \infty$). Throughout, spins are dimensionless unit vectors, and we set $k_B=1$. The working Hamiltonian includes nearest-neighbor exchange, an in-plane Zeeman field, and long-range dipolar couplings,

$$\mathcal{H}=-J\sum _{\langle i,j\rangle }{S}_i·S_j-B\sum _i{n}_B·S_i+D\sum _{i>j}{\displaystyle \frac{S_i·S_j-3(S_i·n_r)(S_j·n_r)}{r_{ij}^3}},$$

(1)

which can be rewritten as:

$$H=-J\sum _{\langle i,j\rangle }cos({\theta }_i-{\theta }_j)-B\sum _{i}cos{\theta }_i-{\displaystyle \frac{D}{2}}\sum _{i>j}{\displaystyle \frac{cos({\theta }_i-{\theta }_j)+3cos({\theta }_i+{\theta }_j-2{\varphi }_{ij})}{r_{ij}^3}}.$$

(2)

where $J>0$ is the ferromagnetic exchange between nearest neighbors $\langle i,j\rangle$, B is an external field (chosen along the y-axis, $n_B=\{0,1\}$), D sets the dipolar strength, $r_{ij}=|{\mathbf{r}}_i-{\mathbf{r}}_j|$, and ${\varphi }_{ij}$ is the polar angle of the vector ${\mathbf{r}}_{ij}$ in the lattice plane (unit vector $n_r$). Unless stated otherwise, energies and temperatures are expressed in units of J, and we use the dimensionless ratio $\alpha =D/J$ to quantify the relative dipolar coupling.

Figure 1. Schematic 2D configuration (Ising limit $Q=2$) on a square lattice with free boundaries for $L=3$.

Figure 1. Schematic 2D configuration (Ising limit $Q=2$) on a square lattice with free boundaries for $L=3$.

Long-range dipolar terms frustrate local alignment and reshape the energy landscape, competing with exchange and field contributions. In two dimensions, this competition intertwines with topology: for $Q\ge 4$ the clock model supports a Berezinskii–Kosterlitz–Thouless (BKT) regime associated with vortex–antivortex physics. In contrast, the underlying discrete Q-fold rotational symmetry ($Z_Q$) of the clock variables leads to an additional discrete-symmetry ordering at lower temperatures [8].

Finite geometry and boundary conditions play a decisive role in this context. In small lattices, where the linear system size L (defined as the number of spins per lattice edge) becomes comparable to the correlation length $\xi \left(T\right)$ measured in lattice units, cooperative fluctuations are strongly suppressed by finite-size effects. In this mesoscopic regime, specific-heat maxima are rounded and shifted, and the proliferation of topological defects is impeded, as correlations cannot fully develop across the system. As a result, intuition based on the thermodynamic limit is not directly transferable to finite samples. Free boundaries further enhance edge contributions and modify bulk-like signatures relative to periodic boundary conditions.

2.2. Computational Strategy and Observables

To capture thermodynamics across different sizes and symmetries, we perform exact spin-configuration enumeration on small lattices. The calculation of thermal averages, a classic problem in statistical mechanics, has been approached through a spectrum of methods: transfer-matrix solutions in one and two dimensions [21,22,23,24], mean-field theory [25,26], the Metropolis algorithm in Monte Carlo simulations [27,28,29], and more recently, through exotic statistical frameworks and diverse machine-learning–based selection techniques [30,31,32,33]. Without loss of generality, we explore exact state counting in order to lay the groundwork for implementing long-range interactions in simulations of this kind.

For different clock models on a $3\times 3$ lattice and for the Ising model with $L=3,4,5$, the complete energy spectrum $\left\{E_n\right\}$ and corresponding degeneracies $\left\{c_n\right\}$ are obtained through exhaustive enumeration of all $Q^N$ spin configurations. By construction, the degeneracies satisfy the closure relation ${\sum }_n{c}_n=Q^N$, ensuring that the full configuration space is accounted for. The energetic levels $E_n$ depend explicitly on the Hamiltonian parameters $\{J,B,D\}$; throughout this work, their relative influence is conveniently characterized by the dimensionless dipolar ratio $\alpha =D/J$. This allows us to construct the canonical partition function

$$Z(T,\{J,B,D\};Q,L)=\sum _n{c}_n{e}^{-E_n\left(\{J,B,D\}\right)/T}.$$

(3)

The complete set of configurations explored, determined by clock symmetry and lattice size, is summarized in Table 1. Here, T denotes the temperature of the thermal bath, while the energetic states incorporate exchange, Zeeman, and dipolar couplings. Thermodynamic observables then follow exactly:

$$\mathrm{U}(T,\alpha )=-{\displaystyle \frac{\partial }{\partial \beta }}lnZ,C(T,\alpha )={\displaystyle \frac{\partial U}{\partial T}},S(T,\alpha )={\displaystyle \frac{U-F}{T}},$$

(4)

with

$$F(T,\alpha )=-TlnZ.$$

(5)

and, when relevant, the magnetization $M(T,\alpha )=-{\displaystyle \frac{\partial F}{\partial B}}$. This enumeration provides noise-free benchmarks for the effects of $\alpha$, Q, and L, and resolves how each Hamiltonian term contributes to the thermodynamic density of states.

When assessing BKT-like behavior for $Q\ge 4$, we track finite-size trends of response functions and, where informative, stiffness-related proxies consistent with the discrete symmetry. Across all calculations, we explore representative sets in the $(J,D)$ plane at fixed B, reporting results as functions of $\alpha$ for Ising and for $Q=4,6,8$ and size L. This mixed exact-numerical protocol isolates finite-size from interaction-induced effects in a controlled manner, enabling a coherent thermodynamic description of the dipolar clock model from the microscopic, enumeration-accessible regime to mesoscopic lattices where collective behavior and long-range couplings jointly dictate the observed caloric and magnetic responses.

2.3. Remark on Computational Feasibility and Model Truncations

The exponential growth of the configuration space with both L and Q is illustrated in Table 1. Even for moderate sizes, the number of spin configurations quickly exceeds any realistic enumeration capacity, reaching ${10}^{32}$–${10}^{44}$ states for $L\ge 6$ and $Q\ge 8$. This scaling motivates the transition from exact enumeration to Monte Carlo sampling. For dipolar interactions, further care must be taken: a naive truncation of the long-range $1/r^3$ potential introduces nonphysical artifacts in the thermodynamic limit [34,35]. For our study, we will limit ourselves to ${10}^8$ configurations.

2.4. Partition Function for Small Systems

We first analyze the thermodynamics of finite lattices by explicitly enumerating microstates. For a given lattice size $L\times L$ ($N=L^2$), clock symmetry Q, external field B (contributing to the lifting of ground-state degeneracies), and dipolar ratio $\alpha =D/J$ (for simplicity $J=1$ and $\alpha$ contains the dipolar information), the canonical partition function reads

$$Z(T,\alpha ;B,Q,L)=\sum _{n=1}^{\lambda }{c}_n{e}^{-E_n/T},$$

(6)

where ${\left\{E_n\right\}}_{n=1}^{\lambda }$ are the distinct energy levels of the Hamiltonian in Equation (2) and $c_n$ their corresponding degeneracies (number of spin configurations yielding $E_n$). Exact enumeration provides $\{E_n,c_n\}$ without statistical error and thus a noise-free benchmark for all thermodynamic observables derived from Z. To gain intuition on how each Hamiltonian contribution shapes the spectrum at small sizes, we also resolve the energy into exchange, Zeeman, and dipolar parts for every configuration,

$$\begin{array}{cc}\hfill E& =E^{\left(\mathrm{J}\right)}+E^{\left(\mathrm{B}\right)}+E^{\left(\mathrm{D}\right)},\hfill \\ \hfill E^{\left(\mathrm{J}\right)}& =-J\sum _{\langle i,j\rangle }cos({\theta }_i-{\theta }_j),\hfill \\ \hfill E^{\left(\mathrm{B}\right)}& =-B\sum _{i}cos{\theta }_i,\hfill \\ \hfill E^{\left(\mathrm{D}\right)}& =-{\displaystyle \frac{D}{2}}\sum _{i>j}{\displaystyle \frac{cos({\theta }_i-{\theta }_j)+3cos({\theta }_i+{\theta }_j-2{\varphi }_{ij})}{r_{ij}^3}}.\hfill \end{array}$$

(7)

From the enumerated set we construct discrete histograms (degeneracy densities) ${\mathcal{G}}^{\left(\nu \right)}\left(E\right)={\sum }_{\mathrm{configs}}{\delta }_{E,E^{\left(\nu \right)}}$ for $\nu \in \{\mathrm{J},\mathrm{B},\mathrm{D}\}$, which compactly encode how each term contributes to the distribution of energies in a finite sample (fixed thermal bath for enumeration).

2.5. Energy Histograms for $Q=2$ and $L=3$

As a concrete reference, Figure 2(a) displays the histograms of $E^{\left(\mathrm{J}\right)}$, $E^{\left(\mathrm{B}\right)}$ (with $B=1$ in units of J), and $E^{\left(\mathrm{D}\right)}$ for the Ising limit $Q=2$ on a $3\times 3$ lattice with free boundaries [36]. The Zeeman contribution exhibits the narrowest spread, consistent with an effectively noninteracting alignment cost relative to the field direction. The exchange term shows a broader distribution reflecting the combinatorial variety of nearest-neighbor bonds. The dipolar contribution is the broadest, owing to its long-range character and angular dependence through $3cos({\theta }_i+{\theta }_j-2{\varphi }_{ij})$, which significantly increases the spread of pairwise energies across the lattice.

2.6. Dependence on the Dipolar Ratio $\alpha =D/J$

Figure 2(b) shows how the spectrum evolves as $\alpha$ is varied for the same $Q=2$, $L=3$ system. Increasing $\alpha$ broadens and skews the dipolar distribution more strongly than the exchange or Zeeman terms, causing the total energy spectrum $\left\{E_n\right\}$ to spread accordingly. At large $\alpha$, finite-size discretization of $\{r_{ij},{\varphi }_{ij}\}$ under free boundaries generates small multi-peak structures intrinsic to mesoscopic clusters but absent in the thermodynamic limit.

2.7. Dependence on the Clock Symmetry Q

Increasing the number of angles from $Q=2$ to $Q=4,6,8$ enriches the set of local orientations $\left\{{\theta }_i\right\}$ and therefore the multiplicity of distinct pair energies, i.e., the interaction energies associated with distinct spin pairs $(i,j)$ for given angular configurations $({\theta }_i,{\theta }_j)$. Figure 3(a) compares histograms at fixed $L=3$ for several Q. As Q increases, the exchange term interpolates toward a quasi-continuous XY-like distribution, and the Zeeman term develops intermediate values between extreme alignments. The dipolar contribution remains the most dispersive and the most sensitive to increasing Q, reflecting the interplay between angular anisotropy and long-range coupling. .

2.8. Finite-Size Effects at Fixed Q

Figure 3(b) shows how the histograms evolve with system size for the Ising case $Q=2$, comparing $L=3,4,5$. Increasing L increments the number of distinct levels and changes their degeneracy patterns. Free boundaries accentuate edge-dominated splittings, which play a key role in understanding finite-size effects on the rounding and shifting of thermodynamic anomalies such as the specific-heat peaks.

The enumerated sets $\{E_n,c_n\}$ obtained in this section directly feed Equation (6), enabling exact evaluation of $\mathrm{U}(T,\alpha )$, $C(T,\alpha )$, $S(T,\alpha )$, and $M(T,\alpha )$ in Sec. 3. These results provide baseline, noise-free thermodynamics for small lattices and establish controlled trends versus $\alpha$, Q, and L for the analysis presented below.

3. Thermal Averages

3.1. Thermodynamic Behavior Under Isolated Interaction Terms

At low temperatures, the Ising model exhibits a ferromagnetic regime driven by the exchange interaction, which favors alignment between nearest neighbors. As temperature increases, thermal fluctuations progressively destroy the ordered state, leading to a paramagnetic phase with vanishing total magnetization. A pure Zeeman term, in contrast, promotes alignment along the external field direction and therefore requires higher temperatures to fully disorder the system. In Q-clock models, increasing Q enhances the number of accessible spin orientations and allows for vortex-like excitations.

Figure 4 shows the temperature dependence of the internal energy and specific heat when each interaction term is activated separately on an $L=3$ lattice. Panel (a) displays $\mathrm{U}\left(T\right)$: the Zeeman term yields the highest internal energy for all Q values, while the exchange interaction produces the lowest energies due to cooperative alignment. The dipolar term stabilizes low-energy configurations, reflecting the strong ordering tendency induced by long-range couplings.

Panel (b) presents the specific heat $C\left(T\right)$. For pure exchange ($J=1$), the cases $Q=6$ and $Q=8$ display a characteristic double-peak structure, associated with a BKT-like crossover followed by discrete $Z_Q$ ordering in finite systems. The Zeeman term ($B=1$) suppresses the low-temperature peak, merging it with the main response into a single dominant feature. Dipolar couplings ($D=1$) broaden the thermal response and shift spectral weight toward higher temperatures, partially merging the two-peak structure for $Q\ge 6$.

Figure 5 shows the entropy $S\left(T\right)$ for the same system under isolated interactions. In the pure exchange case, the system approaches a finite residual entropy $S_0=ln\left(Q\right)$ as $T\to 0$, reflecting the ground-state degeneracy inherent to the Q-clock symmetry. When the Zeeman term is present, this degeneracy is lifted and $S\to 0$ as $T\to 0$. For dipolar interactions, the system retains a finite residual entropy, as long-range couplings favor frustrated configurations rather than a single ferromagnetic ground state.

At high temperatures, the entropy approaches the expected limit $S=Nln\left(Q\right)$, with $N=9$, yielding approximate values of $6.2$, $12.5$, $16.1$, and $18.7$ for $Q=2,4,6,$ and 8, respectively. This asymptotic regime is reached more slowly in the dipolar case, reflecting the reduced statistical weight of nearly ferromagnetic configurations in the presence of competing long-range interactions.

3.2. Finite-Size Effects Under Isolated Interaction Terms

The role of lattice size is examined in terms of the system’s extensivity and edge effects in real space. When thermodynamic quantities such as internal energy, specific heat, and entropy are evaluated per site, one expects them to converge as the number of spins increases, approaching their thermodynamic-limit values. The internal energy per site exhibits distinct behaviors depending on the interaction. In the presence of an external field ($B=1$), the internal energy remains essentially independent of lattice size, reflecting the uniform alignment imposed by the field. In the exchange and dipolar cases, however, the energy per site decreases with increasing lattice size, with the dipolar interaction yielding the lowest energies due to its long-range, partially antiferromagnetic character.

The specific heat per site shows minimal size dependence for the Zeeman interaction but displays a clear trend for the exchange and dipolar terms: the temperature at which the specific heat reaches its maximum increases with system size, signaling enhanced cooperative effects and a shift toward bulk-like behavior. These features are summarized in Figure 6(a)–(b).

The entropy per site, shown separately in Figure 7, reveals complementary finite-size trends. In the exchange case, the residual entropy at low temperatures decreases as the lattice grows, consistent with reduced degeneracy and stronger ordering. The Zeeman field drives the entropy to zero as $T\to 0$, while the dipolar interaction retains a small but finite residual entropy due to frustration. At very high temperatures, all entropy curves converge to their extensive limits. In the zoom of Figure 7, we note that the degeneracy number of energy-minimizing configurations is the same regardless of whether the Hamiltonian contains an exchange or a dipolar term.

3.3. Ising Case: Dependence on the Dipolar Ratio $\alpha =D/J$

We now focus on the Ising limit ($Q=2$) on a $3\times 3$ lattice with $B=0$, varying the dipolar ratio $\alpha \in [0,2]$ in steps of $0.1$. Figure 8(a)–(b) summarize the thermal behavior of the internal energy and specific heat. In Figure 8(a), the internal energy $\mathrm{U}(T,\alpha )$ decreases monotonically with increasing $\alpha$ over the entire temperature range. The inflection point marking the crossover between the low-T ordered regime and the high-T paramagnetic regime shifts to higher temperatures as $\alpha$ increases, reflecting the strengthening of the dipolar interaction.

Figure 8(b) shows the corresponding specific heat. As $\alpha$ increases, the main peak broadens and shifts toward higher temperatures, while its height decreases. In an intermediate window around $\alpha \sim 1$, a low-T shoulder appears, indicative of competing low-energy configurations. For large $\alpha$, this shoulder merges into a single broad maximum. We observe a tendency towards non-zero specific heat values as $T\to 0$ for certain $\alpha$, revealing a thermal signature of changes in the number of energy-minimizing configurations in the ground state of the system.

Figures 9(a)–(b) provide entropy-based and microscopic insight into the same dependence on $\alpha$. Figure 9(a) displays the entropy $S\left(T\right)$ for all values of $\alpha$. As $T\to 0$, the residual entropy plateaus reflect the ground-state multiplicity: near the edges of the explored interval, $S_0=ln2$, while around $\alpha \approx 1$ the plateau increases, consistent with a higher multiplicity such as a fourfold ground state. At high temperature, all curves converge to a disordered state with a resistance to disorder as the dipolar term increases.

Figure 9(b) confirms this interpretation by showing explicitly the ground-state multiplicity as a function of $\alpha$. Distinct plateaus (e.g., two, four, or six configurations) appear in separate $\alpha$ intervals, and the inset highlights the evolution of the lowest-energy levels with $\alpha$. Level crossings between nearly degenerate states account for the observed jumps in the multiplicity and for the entropy plateaus shown in panel (a).

3.4. Q-Clock Models: Dependence on the Dipolar Ratio

Figure 10 summarizes the internal energy for $Q=2,4,6,8$ on a $3\times 3$ lattice (no magnetic field) when the dipolar ratio is varied as $\alpha =D/J\in \{0,0.5,1,1.5,2\}$. At low temperatures, increasing $\alpha$ lowers the energy for all Q, and a clear separation emerges between two groups: $Q=\{2,4\}$ and $Q=\{6,8\}$. This split is consistent with the larger angular freedom of the higher-Q clocks, which enhances the impact of long-range interactions on the low-energy manifold. At high temperatures, each family approaches its interaction-dependent saturation, and for sufficiently strong dipolar coupling (e.g., $\alpha \approx 1.5$) the $Q\ge 4$ curves cluster near the Ising baseline at $\alpha =0$, indicating a partial compensation between exchange and dipolar contributions.

The specific heat reveals a systematic evolution with the dipolar ratio $\alpha =D/J$ across clock symmetries. For $Q=2$ and $Q=4$ (top panel of Figure 11), increasing $\alpha$ shifts the main maximum of $C(T,\alpha )$ to higher temperatures and broadens the curve; in the Ising case ($Q=2$) a low-T shoulder emerges at intermediate ratios, whereas for $Q=4$ the response remains single-peaked with a reduced height. For $Q=6$ and $Q=8$ (bottom panel), the double-peak structure evident at small $\alpha$ is progressively suppressed and merges into a single broad maximum as $\alpha$ grows, while the high-T tails become nearly Q-indistinguishable. These trends indicate that long-range dipolar couplings smoothen exchange-driven anomalies and push the dominant energy fluctuations to higher temperatures.

The entropy exhibits the expected high-temperature ordering with respect to clock symmetry: for fixed $\alpha$, $S/N$ increases monotonically with Q and approaches its asymptotic limit $lnQ$. Increasing the dipolar ratio $\alpha$ systematically lowers $S\left(T\right)$ at intermediate and high temperatures for $Q\ge 4$, indicating that long-range couplings reduce the number of effectively accessible configurations. At low temperatures, the residual entropy $S_0(\alpha ,Q)$ reflects the $\alpha$-dependent ground-state multiplicity: plateaus change with $\alpha$, signaling zero-temperature configuration phases whose degeneracies differ from the pure-exchange expectation. A complementary view is provided by the residual entropy $S_0(\alpha ,Q)$ extracted at low temperature (here evaluated at $T=0.1$): as shown in Figure 13, the dashed baselines at $lnQ$ represent the pure-exchange expectation, while the solid curves display $\alpha$-dependent enhancements and depressions that signal reorganizations of the ground-state manifold. Peaks and step-like features indicate accidental degeneracies created by level crossings as the dipolar term competes with exchange; the effect is broader and more pronounced for larger Q (greater angular freedom), whereas for small and large $\alpha$ the curves tend back toward their $lnQ$ baselines. These trends corroborate the plateaus observed in $S\left(T\right)$ at $T\to 0$ and quantify how long-range couplings reshape the zero-temperature configuration space.

Figure 12. Entropy $S\left(T\right)$ for a $3\times 3$ lattice and clock symmetries $Q=2,4,6,8$ at selected dipolar ratios $\alpha =D/J\in \{0,0.5,1,1.5,2\}$ (color bar from blue to red). For fixed $\alpha$, $S/N$ increases with Q and tends to the high-T limit $lnQ$. Larger $\alpha$ lowers the curves at intermediate and high temperatures for $Q\ge 4$, evidencing the ordering effect of long-range interactions. At low temperatures, $\alpha$-dependent plateaus reveal changes in ground-state multiplicity (residual entropy $S_0$), consistent with $\alpha$-driven configuration phases.

Figure 12. Entropy $S\left(T\right)$ for a $3\times 3$ lattice and clock symmetries $Q=2,4,6,8$ at selected dipolar ratios $\alpha =D/J\in \{0,0.5,1,1.5,2\}$ (color bar from blue to red). For fixed $\alpha$, $S/N$ increases with Q and tends to the high-T limit $lnQ$. Larger $\alpha$ lowers the curves at intermediate and high temperatures for $Q\ge 4$, evidencing the ordering effect of long-range interactions. At low temperatures, $\alpha$-dependent plateaus reveal changes in ground-state multiplicity (residual entropy $S_0$), consistent with $\alpha$-driven configuration phases.

Figure 13. Residual entropy $S_0(\alpha ,Q)$ at $T=0.1$ for a $3\times 3$ lattice and clock symmetries $Q=2,4,6,8$ as a function of the dipolar ratio $\alpha =D/J$. Dashed lines mark the pure-exchange baselines $lnQ$. Deviations from these baselines (peaks/steps) reflect $\alpha$-induced changes in ground-state multiplicity due to level crossings, with more pronounced effects for larger Q.

Figure 13. Residual entropy $S_0(\alpha ,Q)$ at $T=0.1$ for a $3\times 3$ lattice and clock symmetries $Q=2,4,6,8$ as a function of the dipolar ratio $\alpha =D/J$. Dashed lines mark the pure-exchange baselines $lnQ$. Deviations from these baselines (peaks/steps) reflect $\alpha$-induced changes in ground-state multiplicity due to level crossings, with more pronounced effects for larger Q.

3.5. Approach to Magnetization

In the dominant energy regime, i.e., for temperatures below the Onsager critical temperature, an external magnetic field promotes a spontaneous global spin alignment, leading to a ferromagnetic configuration. While the exchange interaction drives local spin alignment between neighboring sites, achieving a fully ordered state requires a larger energetic contribution (see Figure 4 and Figure 5). In contrast, the texture that minimizes the magnetostatic energy in the dipolar case does not favor global alignment; instead, it induces magnetic domains with radial anisotropy, producing two-dimensional patterns similar to sources or sinks. This competition effectively shifts the critical temperature $T_C$ toward higher values. In simple terms, the exchange term reinforces the Zeeman contribution, whereas the dipolar interaction competes with it, defining the dominant texture of the system. To visualize these effects, Figure 14 presents the temperature-dependent magnetization $M\left(T\right)$ obtained from the free energy of the Ising model for $L=3$ and $Q=2$, under two distinct energetic scenarios:

$$\begin{array}{cc}\hfill E& =E^{\left(\mathrm{J}\right)}-B\sum _{i}cos{\theta }_i,\hfill \\ \hfill E& =E^{\left(\mathrm{D}\right)}-B\sum _{i}cos{\theta }_i.\hfill \end{array}$$

(8)

where B denotes the external magnetic field, expressed in units of the exchange or dipolar coupling strength, and varied within the range $B_{min}=0.2$ to $B_{max}=2.0$.

The figure shows that the Zeeman term enhances spin alignment at low temperatures across the entire field range. For the exchange–Zeeman case (blue-green curves), the system saturates rapidly even at weak magnetic fields, reaching maximal magnetization at low temperature. In contrast, for the dipolar–Zeeman case (red-green curves), saturation is not achieved within the same field interval. Above a characteristic threshold $B^{*}$, the competition between dipolar anisotropy and the external field prevents global alignment, leading instead to partially ordered or domain-like spin configurations [37,38,39].

3.6. Low-Temperature Shoulder and Secondary Peak in $C(T,\alpha )$

The evolution of the specific heat for dipolar Ising and clock clusters reveals nontrivial finite-size structure beyond the main ferro–paramagnetic peak. In particular, the presence of low-temperature shoulders and secondary maxima reflects how competing interactions reshape the low-energy spectrum as the dipolar ratio $\alpha =D/J$ is varied. These features are summarized in Figure 15, which condenses the temperature locations of all local maxima of $C(T,\alpha )$ as functions of $\alpha$.

Figure 15(a) focuses on the Ising case ($Q=2$) for lattice sizes $L=3,4,5$. Each symbol marks a local maximum of $C(T,\alpha )$ identified from the full temperature dependence at fixed $\alpha$. The main ferro–paramagnetic peak is present for all $\alpha$ and shifts monotonically to higher temperatures, showing only a weak dependence on L, consistent with bulk-like excitations controlling this crossover. In contrast, the low-temperature peak appears only in restricted intervals of $\alpha$ and is strongly size dependent.

Importantly, the absence of points in certain $\alpha$ intervals does not indicate missing data but rather signals that no local maximum exists in that temperature range. In those intervals, $C(T,\alpha )$ instead develops a local minimum (or vanishes within numerical resolution), implying that the system crosses over smoothly between low-lying energy manifolds without accumulating spectral weight at finite temperature. As will be discussed in detail below, these gaps are closely connected to changes in the structure and degeneracy of the ground state induced by dipolar interactions.

Figure 15(b) shows the evolution of the specific-heat maxima as a function of the dipolar ratio $\alpha$ for a $3\times 3$ lattice and several values of Q. The dominant peak shifts monotonically to higher temperatures as $\alpha$ increases for all Q, reflecting the growing energetic weight of long-range dipolar interactions. In addition, a secondary low-temperature maximum appears over a broad interval around $\alpha \sim 1$. The occurrence of two points with the same color at a fixed $\alpha$ directly reflects this double-peak structure of $C(T,\alpha )$.

For $Q>4$, this behavior is consistent with a BKT-like crossover regime of the dipolar clock model, where competing low-energy configurations and vortex-like excitations give rise to distinct thermal response scales in finite systems. In the Ising case ($Q=2$), the secondary peak remains. Still, it appears shifted and less pronounced, indicating that it arises from finite-size rearrangements of the discrete energy spectrum rather than from a genuine collective crossover or thermodynamic transition.

For the clock models at fixed $L=3$, as shown in Figure 15(b), two distinct $\alpha$ intervals can be identified in which a secondary low-temperature peak emerges, particularly for $Q\ge 4$. These regions correspond to regimes where several low-energy configurations become nearly degenerate, enhancing thermal fluctuations at temperatures well below the main peak of $C(T,\alpha )$. Outside these intervals, the low-temperature structure disappears, again leading to gaps in the peak map that reflect smooth spectral reorganizations rather than sharp thermodynamic features.

Overall, Figure 15 should be interpreted as a map of thermal response scales encoded in $C(T,\alpha )$ rather than as a continuous phase diagram. The presence or absence of secondary peaks reflects how finite-size level crossings, boundary effects, and interaction competition organize the low-energy spectrum. In the following discussions, we complement this analysis with contour plots of $C(T,\alpha )$ and a detailed discussion of the ground-state structure, which together clarify the microscopic origin of these discontinuities.

Figure 16 synthesizes the two complementary mechanisms that govern the specific heat landscape of the model in the absence of an external magnetic field. In the low-temperature limit, the contour plot reveals a sequence of sharp minima along horizontal cuts at $T\approx 0$. These zeros of the specific heat $C(T,\alpha )$ occur at well-defined values of $\alpha$ and correspond to critical points where the ground state changes discontinuously due to level crossings in the energy spectrum. Because thermal excitations are suppressed at these isolated points, the specific heat vanishes as $T\to 0$, producing the pronounced dips highlighted in the lower inset.

In addition to these ground-state transitions, the figure also exhibits clear signatures of thermal crossovers at finite temperature. For intermediate and large values of $\alpha$—where the long-range dipolar interaction becomes comparable to or dominates the exchange term—the system develops a nontrivial temperature dependence. As illustrated by the vertical cuts (upper right panel), $C(T,\alpha )$ shows pronounced maxima arising from the thermal population of excited configurations with competing spin arrangements.

Since the dipolar interaction favors antiferromagnetic alignment, whereas the exchange interaction stabilizes a ferromagnetic arrangement, their competition gives rise to a progressive reorganization of the dominant spin correlations as temperature increases, before the system ultimately crosses over to the paramagnetic regime. This interplay accounts for the broad high-temperature peaks observed at large $\alpha$, including those persisting in the Ising limit. Overall, the $C(T,\alpha )$ map provides a unified view of the ground-state transitions and finite-temperature crossovers shaping the thermodynamic behavior of the system.

Figure 17 provides a direct and transparent confirmation of the critical values of $\alpha$ identified from the specific heat analysis. The dashed branches represent the energies of the lowest competing spin configurations, while the red curve tracks the lowest among them as a function of $\alpha$. The two crossings at ${\alpha }_1\simeq 1.19$ and ${\alpha }_2\simeq 1.91$ mark the points where the ground state changes discontinuously. These are precisely the same values at which the specific heat exhibits zeros in the limit $T\to 0$ in Figure 16, confirming that the low-temperature anomalies originate from level crossings in the energy spectrum.

An essential feature of these crossings is that the level degeneracy is not symmetric around the critical points. As a consequence, the corresponding Schottky anomaly does not produce a symmetric peak in the specific heat on both sides of the transition. Before the first crossing, the energy gap and its degeneracy give rise to a well-defined activated peak. However, beyond the critical value, the ordering of the levels is reversed, and the accessible low-lying excitations change, modifying both the position and amplitude of the Schottky contribution. This explains the asymmetric shape of the low-temperature peaks and the distinct behavior observed before and after each critical point $\alpha$.

Figure 18(a) shows how the low-temperature specific-heat structure evolves as the number of local spin states Q increases for a $3\times 3$ lattice. For $Q=2$ and $Q=4$, the curves exhibit a sequence of sharply defined peaks located at the same values of $\alpha$ identified in Figure 16 and Figure 17. These peaks originate from ground-state crossings, and their positions coincide with the critical values at which the lowest-energy configuration changes discontinuously. In contrast, for $Q=6$ and $Q=8$ the behavior changes qualitatively: the sharp low-temperature features disappear, and the specific heat displays only broad thermally activated maxima. This indicates that, for larger Q, the energy landscape no longer produces ground-state crossings within the relevant range of $\alpha$. The increased number of accessible configurations smooths out the level structure, removing the abrupt degeneracy changes that generate zeros of the specific heat as $T\to 0$. Thus, only the $Q=2$ and $Q=4$ cases retain the critical behavior associated with level crossings, whereas systems with $Q\ge 6$ exhibit purely thermal signatures.

Figure 18(b) highlights a striking dependence of the low-temperature specific heat on the parity of the lattice size L in the Ising case ($Q=2$). For the odd lattices $L=3$ and $L=5$, the curves display exactly two critical values of $\alpha$, consistent with the two ground-state crossings identified in the energy analysis. Immediately after each zero of the specific heat, a Schottky-type anomaly develops; however, these peaks are strongly asymmetric. This lack of symmetry reflects the fact that the degeneracies of the competing energy levels are not equivalent on both sides of the crossing. Once the ground state changes, the multiplicities of the first excited configurations differ, producing an imbalance in the thermal population that manifests as a distorted Schottky profile. The behavior for the even lattice $L=4$ is markedly different. In this case, three critical points appear within the same interval of $\alpha$, leading to three zeros in the low-temperature limit of the specific heat. Moreover, the associated Schottky peaks are significantly more symmetric. Such symmetry indicates that the levels involved in each crossing possess comparable degeneracies, so the activated contribution rises and falls in a nearly balanced manner around the critical point.

Altogether, panels (a) and (b) of Figure 18 show that the symmetry (or asymmetry) of the Schottky anomaly acts as a direct fingerprint of the degeneracy structure of the low-lying spectrum. Odd lattices, with unbalanced multiplicities across the crossings, yield asymmetric peaks, while even lattices tend to produce more symmetric anomalies due to the pairing of degeneracies. This establishes a clear connection between lattice parity, ground-state structure, and the thermodynamic signatures of the model.

Table 2. Summary of low-temperature specific-heat features for the Ising case ($Q=2$) as a function of lattice size L. The number of critical dipolar ratios $N_{\mathrm{crit}}$ corresponds to the zeros of $C(T\to 0,\alpha )$, while the Schottky-peak shape and degeneracy pattern reflect the structure of the low-lying spectrum discussed in Figure 18 and Figure 17.

L	Q	$N_{\mathrm{crit}}$	Schottky-peak shape	Degeneracy pattern
3	2	2	Strongly asymmetric	Unbalanced low-lying degeneracies
4	2	3	Nearly symmetric	Approximately paired degeneracies
5	2	2	Asymmetric	Unbalanced degeneracies (odd parity)

Table 2. Summary of low-temperature specific-heat features for the Ising case ($Q=2$) as a function of lattice size L. The number of critical dipolar ratios $N_{\mathrm{crit}}$ corresponds to the zeros of $C(T\to 0,\alpha )$, while the Schottky-peak shape and degeneracy pattern reflect the structure of the low-lying spectrum discussed in Figure 18 and Figure 17.

L	Q	$N_{\mathrm{crit}}$	Schottky-peak shape	Degeneracy pattern
3	2	2	Strongly asymmetric	Unbalanced low-lying degeneracies
4	2	3	Nearly symmetric	Approximately paired degeneracies
5	2	2	Asymmetric	Unbalanced degeneracies (odd parity)

4. Conclusions

In this work, we have carried out a detailed and fully controlled thermodynamic analysis of small dipolar Q-state clock lattices, with particular emphasis on the $3\times 3$ geometry, where exact enumeration allows the complete microscopic spectrum and its degeneracy structure to be resolved without approximation. This approach demonstrates that, even in the Ising limit ($Q=2$), the inclusion of long-range dipolar interactions is sufficient to generate nontrivial low-temperature caloric features. In particular, additional low-temperature peaks appear in the specific heat at well-defined values of the interaction ratio $\alpha =D/J$, signaling a reorganization of the low-energy manifold that is absent in the pure exchange model.

A central result of this study is the identification of sharp critical values of $\alpha$ that manifest as zeros of the specific heat in the limit $T\to 0$. These minima coincide exactly with discontinuous changes in the ground state driven by level crossings in the energy spectrum, establishing a direct and unambiguous correspondence between microscopic spectral rearrangements and macroscopic thermodynamic signatures. Contour maps of $C(T,\alpha )$, together with an explicit tracking of the lowest energy branches, confirm that each zero of $C(T\to 0,\alpha )$ originates from such a ground-state crossing.

By comparing different lattice sizes, we have shown that the shape of the Schottky-like anomaly following each critical point encodes detailed information about the degeneracy structure of the competing low-lying levels. For odd lattices ($L=3$ and $L=5$), the post-critical peaks are strongly asymmetric, reflecting unbalanced multiplicities across the crossings. In contrast, the even lattice ($L=4$) exhibits three critical values of $\alpha$ accompanied by nearly symmetric Schottky anomalies, indicating that the levels exchanging stability possess comparable degeneracies. These results highlight lattice parity as a key organizing principle for the pairing of degeneracies and provide a natural explanation for the symmetry or asymmetry of the observed caloric response.

For larger clock symmetries, the phenomenology changes qualitatively. While systems with $Q=2$ and $Q=4$ retain sharp critical values of $\alpha$ associated with ground-state crossings, for $Q\ge 6$ the energy landscape becomes sufficiently smooth that such crossings are suppressed within the explored parameter range. As a result, the low-temperature anomalies disappear and the specific heat displays only broad, purely thermal maxima. The absence of a genuine Berezinskii–Kosterlitz–Thouless transition in the finite clusters studied here is not a consequence of dipolar interactions alone, but rather reflects the combined effect of finite system size and the underlying discrete clock symmetry, which preclude the development of true topological order. Nevertheless, the evolution observed for $Q\ge 6$ bears a qualitative resemblance to the crossover behavior associated with vortex–antivortex physics in the two-dimensional XY model, in the sense that increasing angular freedom progressively smooths degeneracy-driven features.

Taken together, our findings demonstrate that the dipolar Q-clock model, even at the smallest accessible scales, exhibits a rich, highly structured interplay among long-range anisotropy, ground-state degeneracy, and finite-size geometry. The combined analysis of specific-heat extrema, residual entropies, and degeneracy patterns provides a coherent microscopic interpretation of the thermal anomalies induced by dipolar interactions. These results establish small dipolar clusters as controlled benchmark systems for understanding finite-size thermodynamics in van der Waals nanomagnets, patterned dipolar arrays, and mesoscopic magnetic textures. They also provide a firm foundation for future extensions to larger clusters and to dynamical or field-driven regimes, where the same competition between exchange and dipolar couplings is expected to generate similarly intricate thermal landscapes.