Entropy and Mutability for the q-State Clock Model in Small Systems

In this paper, we revisit the q-state clock model for small systems. We present results for the thermodynamics of the q-state clock model for values from q=2 to q=20 for small square lattices of L×L, with L ranging from L=3 to L=64 with free-boundary conditions. Energy, specific heat, entropy, and magnetization were measured. We found that the Berezinskii–Kosterlitz–Thouless (BKT)-like transition appears for q>5, regardless of lattice size, while this transition at q=5 is lost for L<10; for q≤4, the BKT transition is never present. We present the phase diagram in terms of q that shows the transition from the ferromagnetic (FM) to the paramagnetic (PM) phases at the critical temperature T1 for small systems, and the transition changes such that it is from the FM to the BKT phase for larger systems, while a second phase transition between the BKT and the PM phases occurs at T2. We also show that the magnetic phases are well characterized by the two-dimensional (2D) distribution of the magnetization values. We made use of this opportunity to carry out an information theory analysis of the time series obtained from Monte Carlo simulations. In particular, we calculated the phenomenological mutability and diversity functions. Diversity characterizes the phase transitions, but the phases are less detectable as q increases. Free boundary conditions were used to better mimic the reality of small systems (far from any thermodynamic limit). The role of size is discussed.


Introduction
The q-state clock model is the discrete version of the famous 2D XY model, which is probably the most extensively studied example showing the Berezinskii-Kosterlitz-Thouless (BKT) transition [1][2][3][4] in the presence of a frustrated quenched disordered phase [5][6][7][8]. The q-state clock model is often used as a reference model due to its peculiar critical behavior at the transition point and its universal features [9][10][11][12][13]. Instead of the exclusion of the explicit continuous symmetry essential for the BKT transition, it can also emerge from a system without explicit continuous symmetry [13]. The Hamiltonian of the q-state clock model can be written in many forms, and one of the simplest forms is the following expression, where no magnetic anisotropies or external magnetic fields are included: where J > 0 is the ferromagnetic coupling connecting pairs of nearest neighbors i and j; the discrete angle between the spin orientations is given by θ i,j (η, q) = 2πη/q, for η = {0, 1, .., q − 1}. While the exact XY model is recovered only in the limit of infinite q, it has been found that the BKT characteristics appear in the clock models when q ≥ 5 [14][15][16][17][18][19]. The nature of the phase transitions in the general clock model has been widely studied with different theoretical and numerical approaches. However, these studies have given mixed results for the characterization of transitions at the lower bound of q (for instance, see the summary of the related debates in [20]). In the present work, we want to consider the clock model as a generalization of the Ising model by establishing the similarities and differences that arise due to the increase of the degrees of freedom due to the local states at each site. All of the characterizations are aimed at the behavior of small systems: specifically, we focus on square lattices ranging from 3 × 3 up to 64 × 64. Free boundary conditions are preferred since they better represent the importance of surface states in small systems. The clock model systems under scrutiny range from q = 2 (equivalent to the usual Ising model) to q = 20.
There are various scattered results for the thermodynamics and phase transition of the clock model. Therefore, we report below a consistent compendium of its thermodynamic properties, such as internal energy, U; specific heat, C; entropy S; and magnetization. We also report the transition temperatures T 1 for transitioning from the ferromagnetic (FM) to the P phase for very small q values, continuing to the transition from the FM to the BKT phase for larger values of q, as well as T 2 for the transition from the BKT phase to the disordered paramagnetic (P) phase. These series of results are elaborated into a phase diagram of T C vs. q, where T C is determined from the C(T) curves.
The information content of a sequence was measured by the mutability, called as such for the first time by Vogel et al. [21] in a relationship with the characterization of the 2D phase transitions in the Edwards-Anderson spin system as an alternative to the Binder cumulant analysis [22]. Later, an appropriate information recognizer was proposed (named word length zipper (wlzip)) which optimizes recognition of digital information associated with properties of the system [23]. The method was later successfully applied to the reentrant phase diagram in the case of the 3D Edwards-Anderson model [24]. Successive applications of the information content methods addressed stock markets [25], pension funds [26], blood pressure [27], seismology [28], nematic transitions [29], and wind energy production [30].
In the present paper, we follow the approach of magnetic transitions [21,23] complemented with the definitions of mutability and diversity [25], to be elaborated below. The aim is to establish an alternative method to characterizing and distinguishing the phases present in the q clock model.
In the next section, we present the clock model and the main methods used to characterize it. A presentation of the results and discussion are given in Section 3, and the last section is devoted to conclusions.

General Definitions
Let us begin by considering the q-state clock model on a two-dimensional (2D) square lattice of dimensions L × L = N, where the local magnetic moment or "spin" S i at site i can point in any of the q directions in a given plane. S i is then a 2D vector, i.e., S i = (cos( 2π q k), sin( 2π q k)), where k = 0, 1, ..., q − 1, with equal probability for all q values. S i are dimensionless vectors of magnitude one.
The isotropic Hamiltonian for such a system can be written as: where J > 0 is the ferromagnetic exchange interaction to nearest neighbors i and j; the sum extends to all such pairs through the lattice, which is indicated by the symbol < i, j > under the summation symbol. B is an external field applied along one direction in the plane. In this Hamiltonian, J is one unit of energy and B is also measured in energy units. This form is completely equivalent to the one of Equation (1), when B = 0.

Exact Theoretical Approach for a Small System
Let us begin by considering the theoretical approach for the q-state clock lattice with L = 3, which was introduced in the previous section. The partition function can be expressed as: where the coefficients c n are the numbers of all of the possible spin configurations compatible with energy E n according to the Hamiltonian of Equation (2); in the case of a magnetic field B, the energies and the coefficients c n depend upon B; λ is the number of different energy levels. We express energy and temperature in the same units, so the Boltzmann constant is set to unity, k B = 1.
The coefficients c n = c(E n , q) for a small N = L × L = 3 × 3 lattice can be straightforwardly calculated by computing the energies for all q 9 spin configurations, which are shown in Table 1 for q = 4 and q = 6 and in Table 2 Table 1), some symmetry rules apply: energy distribution is symmetric around E n = 0, and the majority of the density of states occur for E n = 0 and c(−E n , q) = c(E n , q) holds. There are q 9 states to be spread among all available energies in the energy range E n = [ −12, 12]. Due to symmetry, the total number of different energies are λ = 11, 23, 47, 289, and 699 for q = 2, 4, 6, 8, and 10, respectively. Table 1. Coefficients c(E n , q) for a 3 × 3 lattice with free boundary conditions and B = 0 for q = 4 and q = 6. The first column enumerates energy levels, n, and the second and third columns give the corresponding energies E and degeneracies c(E, 4), respectively. The fourth and five columns show half of the energies and degeneracies for q = 6, respectively. The rest of the energies and degeneracies can be found using the following symmetry: c(−E n , 6) = c(E n , 6). For odd q values, the main symmetry of the Hamiltonian of Equation (2) for B = 0 is lost, because the spin inversion S i → −S i is not possible. This is simply because, for odd q values, if there is the symmetry S i , then −S i does not exist. Therefore, the energy distribution and the corresponding degeneracy is not symmetric around E n = 0. Thus, the highest possible energy is not the negative value of the ground state energy. However, as we will see, this lowering of symmetry does not affect the thermodynamic observables.
Once the partition function is known, the thermodynamic observables can be calculated directly. Thus, for the cases of internal energy U, specific heat C, and entropy S, they can be obtained by employing: The next table shows the results for q = 5, where the energies E n and degeneracies for n = 1 to n = 85 are explicitly given. In this case, for a given energy E n , their negative counterpart does not exist. Table 2. Coefficients c(E n , 5) for a 3 × 3 lattice with free boundary conditions and B = 0. In this case, there are no symmetries and, therefore, we present all E n values with their respective degenerations, from n = 1 to n = 85. Here, as explained in the manuscript, the symmetry S i → −S i is not present. We now use numeric simulations to calculate larger sizes for the same system.

Numerical Simulations
In addition to the theoretical calculations, most of the work presented here deals with numerical calculations based on Monte Carlo (MC) simulations. A square lattice L × L was chosen; free boundary conditions were imposed; a site was randomly visited, and the energy cost, ∆, of rotating the corresponding spin around q possible states was calculated: if the energy is lowered, the change of orientation is accepted; otherwise, only when exp(−∆/T) ≤ r is the spin rotation accepted, where r is a freshly generated random number in the range [0, 1] with equal probability. This is the usual Metropolis algorithm. A Monte Carlo step (MCS) is reached after N = L × L spin-rotation attempts. One of the main goals here is to report the sensitive temperatures, T 1 and T 2 , that define the transition from FM to BKT and from BKT to PM (disordered) phase for different systems.
For each lattice size and q value, a sequence of temperatures was defined in the range [0.02, 3] at steps of 0.02 for each temperature. Also, 5τ MCSs were performed: the first τ MCS was used to equilibrate at a fixed temperature T, while the next 4τ MCS was used to measure the observables every 20 MCSs, reaching a total of 2 × 10 5 = 200, 000 measurements. Unless otherwise specified, τ = 10 6 is used in the rest of the paper; this τ value gives stable results and leads to coincidence with the analytic expressions obtained as described above. The energy U(T) was computed according to the Hamiltonian given by Equation (2). With the internal energy known as a function of T, Equations (4)-(6) can be used to generate the thermodynamics. In parallel, the magnetization M(T) for each temperature can also be instantaneously measured. Alternatively, we can also use direct relationships based on the thermal treatment of the variables, as presented next.

Thermal Averages
The lattice average of the spin configuration, equivalent to the magnetization per site M, is given by the following expression: where S j is the value of the spin at site j at a given time t, and N = L × L is the total number of spins.
In this particular case, M is a vector of two components, M = (M x , M y ). Normally, the magnitude or absolute value of this vector is calculated, i.e., |M| = M 2 x + M 2 y . Then, the thermal average of the absolute value |M| is < |M| >, and it is given by where N c = 2 × 10 5 is the number of configurations used to perform thermal averages, as explained in the preceding section. Energy is the main value used in the Monte Carlo method to reach thermal equilibrium. Therefore, after τ MCSs, the internal energy U can be obtained by averaging the N c = 2 values for E k , where k runs over the accepted configurations after the Metropolis algorithm, namely: where every spin configuration is separated from the next one by 20 MCs. The energy per site is then U/N, which is the thermal average of the lattice average of the system energy. The specific heat is then calculated as proportional to the fluctuations of the energy as follows: The absolute entropy S can be calculated by calculating ∆S(T f , T i ) = S(T f ) − S(T i ) by numerical integration of the specific heat divided by the temperature, as follows: We know the entropy at zero temperature, because we know the energy degeneration at T = 0. For the case B = 0, only one spin configuration has minimum energy (every spin aligned to B); therefore, S(0) = 0. On the other hand, when the magnetic field is zero, there are q ferromagnetic spin configurations with equal energies. Therefore, S(0) = ln q and, hence, S(T) = ln q + ∆S(T, 0).
As we have seen, the thermal average of a physical quantity is a summation over N c quantities (normalized to N c ), like energies, magnetic moments, etc. However, the order of the sequence of these N c quantities is totally irrelevant to the result of this evaluation. Next, we introduce the mutability, which is a quantity (defined in terms of information theory) that can be calculated from any temporal sequence of quantities, like energy, magnetic moment, etc. Mutability is a quantity that depends on the sequence order, and, therefore, it contains information that differs from a mere average.

Information Theory, Mutability, and Diversity
Word length zipper (wlzip for short) is a compressor designed to recognize meaningful information in a data chain. As its name indicates, it recognizes "words" of precise length and precise location within the data chain. Then, it compresses less than other file compressors which precisely optimize this function. In the case of wlzip, the optimization is of the recognition of patterns beginning at a precise location and for a given number of digits. In this way, physical properties and/or other inherent properties of the system can be recognized.
The recognition of information within a file can render at least two parameters: mutability and diversity. They were defined a few years ago, including a working example, given in Table 2 of Reference [23]. Here, we very briefly review these definitions, beginning with the basic rules under which wlzip operates. Let us assume our data are contained in a vector file (one register per row) named "data.txt" whose number of rows is λ(data) and whose weight in bytes is w(data). We now create a map of the previous file in the following way: we go over each row of data.txt and, whenever this value is new, we add it as a new row to a new file which will be named map.txt. So, the process begins with the first element of data.txt which is also the first element in map.txt. Then, the second element of data.txt is considered: if this element is new, we add it as a second row in map.txt; if this element is the same as the previous one, we just add a digit 1 to the right of this element in map.txt (here, 1 means one position from last time this element appeared). The process continues like that so that each time a new record is detected, a new row in map.txt is created. This register is copied in map.txt at the beginning of a new row, and the first number that is written to its right is its absolute position in the data.txt file. Each time a register is repeated, its position with respect to its last appearance in the data.txt file is written to the right of the last annotation on this row in the map file. At the end of the process, map.txt contains as many rows as there were different values present in data.txt; the more a value is repeated, the wider the corresponding row is in map.txt. In a sense, this is like a histogram organized according to the appearance order. The number of rows in map.txt is λ * (data), while the weight of this file is w * (data).
The mutability µ(data) of the entire file data.txt is simply given by Similarly, the diversity of the file data.txt is also a ratio, namely, The previous definitions can also be dynamic at time t and given with respect to the ν records in data.txt counted from that one at time t. Then, we can define the dynamic parameters for data.txt: and As shown in the discussion accompanying Figure 1 of Reference [25], referred to the Ising model, mutability can be more appropriate for discussing the variability of a parameter with changing conditions (like an increase in temperature), while diversity is more appropriate for discussing critical phenomena, such as phase transitions. This is one of the issues discussed below in light of the results obtained from the clock model.

Results and Discussions
Let us begin by considering the theoretical approach to the q-state clock model according to Equation (2) for B = 0, considering an L = 3 lattice, as presented in the previous section. We successfully calculated the partition functions from q = 2 (Ising) to q = 10. From here, all the thermodynamic quantities of the system can be analytically obtained as a function of the temperature T. Figure 1 shows the results obtained for internal energy U(T), specific heat C(T), and entropy S(T). For the specific heat, C vs. T (center of Figure 1), we see that the two peaks appear for q ≥ 6. At q = 5, we see only one peak, but it is notoriously more skewed to the right. We will see that this is due to the small size of the lattice. Therefore, for small systems, the BKT phase appears only for q ≥ 6.
We can see here (Figure 1, right) the basic features of entropy for the low-and high-temperature limits. When no external field is applied, the energy ground state of the ferromagnetic q-state clock model has a degeneracy equal to q, independent of the total number of spins N; therefore, S(0) = ln q. On the other hand, at very high temperatures, the exchange interaction is overridden, and every spin has q degrees of freedom with equal probabilities: hence, the system degeneracy is equal to q N , and, therefore, Figure 2 shows the same observables in the presence of magnetic field B, namely, U(T, B), C(T, B), and S(T, B), for the case where q = 7, as an example. The shift of the transition temperatures due to the variations in the magnetic field is clearly visible. This is also an exact analytical result obtained by calculating straightforwardly the partition function of Equation (3). The magnetic field was applied along the (1, 0) direction, along which a spin orientation is always possible, regardless of the q value. For B > 0, the Zeeman term is added to the energy (see Equation (2)); therefore, the FM ground state energy is lowered by the external field. As an example for B = 2, the Zeeman energy is −18, and this is added to the ground state energy (−12) at zero field; therefore, the new ground state energy is −30 (red curve in the left of Figure 2).
We also observe that the shapes of specific heat curves, C vs. T, maintain the two peaks, but they shift to higher temperatures as the external field goes from B = 0 to B = 2. This behavior is easily explained by the external field favoring ordered phases (FM and BKT) over disordered ones, and, therefore, the transition temperatures increase with the strength of the external field.
The magnetic field breaks the degeneracy of the ground state so that S(0) = 0 for B = 0 instead of S(0) = ln q. However, at high temperatures, the Zeeman term in Equation (2) is overridden, and we are back to the same situation as in previous analysis for B = 0, namely, every spin has q degrees of freedom with equal probabilities: hence, the system degeneracy is equal to q N , and, therefore, the entropy S(T >> J) = N ln q, as observed in the right of Figure 2.

Monte Carlo Simulations
Next, we present and discuss the output from Monte Carlo simulations made for the q-state clock model in lattices of up to 64 × 64 with free boundary conditions. We began by simulating a 3 × 3 lattice for different q values and comparing these numerical results to the analytic ones presented in the previous section. We did not find a single difference, which is expected of course, but which also serves as a check for the computer programs used extensively in the simulations reported next. Thus, thermodynamic observables for lattices 10 × 10, 16 × 16, 32 × 32, and 64 × 64 are presented in Figures 3-6, respectively.    We observe an overall self-agreement in the shape of energy U, specific heat C, and entropy S as functions of temperature for the different lattices. For a given size, the peak of the specific heat occurs at a lower temperature as q increases, and then it splits in two peaks. As q continues to increase, the high-temperature peak remains at the same temperature, whereas the low-temperature peak tends toward lower temperatures (eventually to a temperature of zero) as q increases tending to infinity. This low-temperature peak, T 1 , is the transition from the FM phase to the BKT-like phase, which is characterized by vortex spin configurations and FM spin-spin correlated configurations, like waves. Both the vortex and spin waves are low-energy excitation that occurs as q increases. Therefore, for q → ∞, we expect that T 1 → 0. We would like to stress the size-dependence of the peak heights at fixed q. Indeed, for q = 2, 3, and 4, the unique peak corresponds to a second-order phase transition, hence, the specific heat C/N is expected to diverge with N at T = T 1 . However, for q ≥ 5, the two peaks should correspond to two BKT transitions, which are infinite-order; that is, there is no divergence in C/N at T 1 and T 2 . Both features are observed in Figures 3-6, and, furthermore, the constant height of the peak at T = T 2 for q = 5 could represent further evidence of a BKT-like transition. To be more specific on the characterization of the clock model, we discuss next a phase diagram including the three possible magnetic phases: FM (long-range order), BKT (short-range correlations), and PM (total disorder).

Phase Diagram
In this section, we show the phase diagram for the clock model as extracted from the specific heat. We collected the analytic results for the 3 × 3 lattice and the numerical results for the 10 × 10 and the 64 × 64 lattices. This is shown in Figure 7 for q ranging from 2 to 20 using squares, triangles, and circles, respectively. The texture and color underneath illustrate the instantaneous magnetic phases for the larger lattice and q = 9 just as examples for the kind of magnetic ordering present in each phase. Phase diagram for the q-state clock model for a 3 × 3 lattice (squares), a 10 × 10 lattice (triangles), and a 64 × 64 lattice (circles) with free boundary conditions. For temperatures under T 1 , the ferromagnetic (FM) ordered phase dominates. For q < 5, the transition from this FM ordered phase is to the disordered PM phase. For large enough lattices and q ≥ 5, the transition from the FM phase is to the partially ordered Berezinskii-Kosterlitz-Thouless (BKT) phase. Then, a transition at a higher temperature T 2 appears which separates the BKT phase from the paramagnetic (PM) phase; it can be noticed that T 2 is essentially constant with respect to q. The background colors blue, yellow, and red (from bottom to top) mark the FM, BKT, and PM phases, respectively. Spin orientations for one possible snapshot corresponding to the 64 × 64 lattice and q = 9 are given in the gray color over the background colors. Figure 7 deserve special attention. First, the lower critical temperature T 1 follows a monotonous decrease with q approaching zero asymptotically. Second, the higher critical temperature T 2 remains at a constant value for q ≥ 6. Third, for q = 5, there is just one critical temperature following the tendency of T 1 , but for L = 10 (and also for L = 64), the two transitions are clearly visible for q = 5. Fourth, an order parameter beyond the usual magnetization is necessary to distinguish the BKT phase from the PM disordered phase.

Several features in
To further specify the different phases, Figure 8 shows the specific heat for a 32 × 32 lattice as a function of temperature. In this figure, we show a two-dimensional (2D) order parameter at certain characteristic temperatures which clearly discriminate the three (FM, BKT, and PM) different phases. The order parameter we use is the 2D distribution of the variable M = (M x , M y ), as defined in Equation (7), that is, the spin-lattice average at time t after the thermalization process of τ MCSs.
This corresponds to the vector spin average for a given spin configuration at time t. This vector is then calculated after the thermalization, and every 20 MCSs of a total of N c times. By plotting all N c vectors, we generate a 2D distribution that clearly characterizes the different phases. In the FM phase, only certain directions of the spin are allowed. In Figure 8, constructed for q = 9, we see that at low temperature, the average magnetization vector points in nine directions which correspond to the nine-fold symmetry of the FM states with equal probability. The BKT-like phase is characterized by spin waves and vortex structures; therefore, the lattice average of the spin points to any of the 2π directions while conserving, to a great extent, its magnitude in every lattice average. Therefore, a ring structure is formed. In the disordered (PM) phase, every spin in the lattice points randomly to any direction; therefore, the lattice average magnetization distribution exhibits a 2D-Gaussian peak with decreasing magnitude as T increases; hence, the circle begins to be filled with a higher probability (red color) near the center. The color code goes from black (zero value) to purple, blue, green, yellow, and red in increasing order of probability for this 2D order parameter. This parameter was already introduced as a complex order parameter by Baek et al. [31].
The next figure (Figure 9) depicts snapshots of some spin configurations showing the spin arrangements that occur at different temperatures.
The next figure (Figure 10) shows the magnetization modulus as the thermal average of the spin-lattice average, as defined in Equation (8). This makes the thermal average magnetization (the 2D order parameter) rotate while keeping its overall magnitude. At T = 0.27, the FM phase is seemingly still present, as can be seen in the ring structure with nine maxima, as observed at T = 0.21. In T = 0.5, a pure BKT phase is observed as a uniform ring structure, but their radii are beginning to shrink, maintaining a smaller magnitude. At T = 1 is the onset of the disordered phase in which the magnitude of the magnetization decreases, filling the interior of the circle as a 2D Gaussian distribution with a zero average. At T = 1.51, the overall magnitude of the M parameter further shrinks. The color code of the 2D order parameter distribution goes from black to purple, blue, green, yellow, and red as values of the 2D distribution increase. The radii of the ring-like 2D distributions represent the thermal magnetization modulus as a function of temperature, and this is given in Figure 10.  Let us now consider the information content as an independent test to characterize these phase transitions. We shall concentrate on the simulations for 32 × 32 lattices, measuring µ and δ as defined in the Models and Methods Section. Figure 11 presents the information content results for the same energy series generated by the previously described MC algorithm; the case of a 32 × 32 lattice for q = 2 was chosen for this report. The open-squares curve (red) gives the results for mutability, and it presents a local maximum near T = 2.2, in agreement with the specific heat results (see Figure 5). The solid-circles curve (blue) is the result for diversity and shows a sharper absolute maximum at T = 2.2. The vertical dashed line illustrates the temperature 2.269, at which the Onsager solution predicts the transition for an extremely large system [32]. This is just a reference since our systems are finite and have free boundary conditions, so they do not correspond well with this theoretical solution. The previous results were obtained from sequences presenting energy data. However, there are more sensitive parameters associated with states of the system whose sequences present alternatives which can recognize this phase transition in a better way. These are the cases of the magnetic order parameters (magnetization, neighbor correlations, site-memory correlation) [25] which we initially tested and obtained encouraging results. To consider this additional study is beyond the goals of the present paper, so we will stick to the energy data results analyzed for diversity to directly compare with the results reported above.
In Figure 12, we present the results for δ in the cases for q = 3, q = 4, q = 6, and q = 8. The ordinates were multiplied by the appropriate constants to fit a common arbitrary scale since the meaningful information is in the temperatures that are marked by each maximum or the inflections of the curves. For q = 3, just one critical temperature is visible and in agreement with the maximum of the specific heat for this lattice size, as shown in Figure 5. A similar situation is observed for q = 4, with a maximum at a lower temperature of around 1.1. For q = 6, a maximum of around 1.1 and a "knee" just over 0.5 are visible, in agreement with the two maxima for this value of q reported in the figure for specific heat. In the case of q = 8, the maximum near 1.1 is clearly present, although it is a bit broader than that for q = 6. The "knee" is a barely visible tiny change in slope at T = 0.5. The results shown by the previous two figures, as well as similar ones for other intermediate values of q, show that the information method is able to recognize the transitions present in the clock model for low values of q when applied to the energy data vectors produced by the MC simulations. However, the phases are less detectable by this method as q increases.
The previous shortcomings of the information theory method could be due to the system itself and the way in which the transition is characterized. With respect to the former, it should be noticed that as q increases, the number of possible configurations increases, the energy intervals decrease, and the density of states is highly degenerate and tending to a continuum as q → ∞; to distinguish among energy values is now increasingly difficult. With respect to the latter, energy is not the best order parameter to characterize these transitions, and, eventually, different order parameters more oriented to the recognition of magnetic states of the system (magnetization, site memory, neighbor correlations) can render better results. At the moment, this is an open question, and it should be addressed in future work.

Conclusions
Using analytically derived expressions and Monte Carlo simulations, we explored the q-state clock model for square lattices with free boundary conditions which better mimic the properties of small systems for which this approach is intended. We calculated their thermodynamic properties and characterized the three magnetic phases present for this model. The corresponding magnetic phase diagrams were calculated for lattice sizes of up to 64 × 64. It turns out that there exists an FM below a critical temperature T 1 . For q ≥ 5 and lattice sizes over 10 × 10, a BKT-like phase appears for temperatures between T 1 and a second critical temperature T 2 , separating this BKT phase from the disordered paramagnetic phase. The BKT phase reflects partial magnetic ordering characterized by spin vortex configurations and zones with FM spin-spin correlation, reflecting curling, or wave-like ordering. The three phases can be well characterized using the spin-lattice average distribution. This 2D distribution shows characteristic patterns which clearly identify the three (FM, BKT, and PM) phases.
The entropy of the system always increases with temperature, showing subtle slope changes at the transition temperatures T 1 and T 2 . The low-temperature limit for the entropy is simply given by ln(q) in absence of a magnetic field, while it vanishes for any magnetic field that breaks the degeneracy of the ground state, yielding a singlet as a ground state. On the other hand, the high-temperature value for entropy tends asymptotically to N ln(q), thus reflecting that all degrees of freedom are equally probable.
The information theory method produces results that are in agreement with those of the specific heat of the system. It distinguishes the phase diagram presented in Figure 7 for low values of q, thus confirming previous results obtained by conventional treatments. However, this recognition is progressively lost as q increases. This is probably due to the fact that the information recognition points to the energy values which are shared by increasingly more states as q increases (contribution of the internal degrees of freedom to the density of states). Eventually, different simulations involving a more elaborated order parameter are needed to address the slight difference between the BKT phase and the disordered phase, both with vanishing magnetization. This work is in progress and should produce results in the near future.