Phase transition properties via partition function zeros: The Blume-Capel ferromagnet revisited

TL;DR

This paper revisits the analysis of phase transitions in the Blume-Capel model using partition function zeros, showing that small systems and off-critical simulations can reliably estimate critical exponents, enhancing efficiency and sustainability.

Contribution

It extends previous work by demonstrating that accurate critical exponent estimation is possible with small system sizes and off-critical simulations, reducing computational effort.

Findings

Small system sizes suffice for reliable exponent estimates.

Simulations off the nominal transition point are effective.

Method improves efficiency and reduces computational footprint.

Abstract

Since the landmark work of Lee and Yang, locating the zeros of the partition function in the complex magnetic-field plane has become a powerful method for studying phase transitions. Fisher later extended this approach to complex temperatures, and subsequent generalizations introduced other control parameters, such as the crystal field. In previous works [Moueddene et al, J. Stat. Mech. (2024) 023206; Phys. Rev. E 110, 064144 (2024)] we applied this framework to the two- and three-dimensional Blume-Capel model -- a system with a rich phase structure where a second-order critical line meets a first-order line at a tricritical point. We showed that the scaling of Lee-Yang, Fisher, and crystal-field zeros yields accurate critical exponents even for modest lattice sizes. In the present study, we extend this analysis and demonstrate that simulations need not be performed exactly at the nominal transition point to obtain reliable exponent estimates. Strikingly, small system sizes are sufficient, which not only improves methodological efficiency but also advances the broader goal of reducing the carbon footprint of large-scale computational studies.

Tables

Table 1. Table 1 : Expected value y h IM = 1.875 y_{h}^{\rm IM}=1.875 . Best fits (in terms of χ s 2 \chi^{2}_{s} ) are highlighted in bold. For Δ = 0.991 \Delta=0.991 , the L = 56 L=56 data were excluded due to poor quality.

$\Delta$	$T$	$y_h^{\mathrm{IM}}$	${\chi }_s^2$
0.5	1.5640	$1.8898\left(1\right)$	207.05
	1.5657	$1.8764\left(1\right)$	$0.18$
	1.5658	1.8756(1)	0.25
	1.5659	1.8748(1)	1.77
0.991	1.400	1.8822(1)	55.57
	1.4005	1.8779(1)	8.33
	1.4006	1.8770(1)	3.89
	1.4007	1.8761(1)	1.16
1.5	1.1510	1.8685(1)	19.33
	1.1506	1.8725(1)	3.58
	1.1504	1.8745(1)	1.91
	1.1503	1.8755(1)	2.57
1.75	0.9580	1.8490(2)	132.62
	0.9558	1.8758(2)	2.19
	0.9559	1.8746(2)	1.00
	0.9560	1.8735(2)	0.42
	0.9561	1.8723(2)	0.45

Equations

$$\mathcal{H}=-J\sum _{(i,j)}{\sigma }_i{\sigma }_j-H\sum _i{\sigma }_i+\Delta \sum _i{\sigma }_i^2=E_J-HM+\Delta E_{\Delta },$$

$$Z_N(\beta J,\beta H,\beta \Delta )=\sum _{E_J=-dN}^{dN}\sum _{M=-N}^N\sum _{E_{\Delta }=0}^N\rho (E_J,M,E_{\Delta }),x^{E_J}{y}^M{z}^{E_{\Delta }},$$

$$\mathrm{Im}\left(H^{(j)}\right)\sim L^{-y_h},$$

$$\mathrm{Im}({\beta }^{(j)})\sim L^{-y_t},$$

$$\mathrm{Im}({\Delta }^{(j)})\sim L^{-y_g},$$

$$y_t^{\mathrm{IM}}=1,y_h^{\mathrm{IM}}=\frac{15}{8},$$

$$y_t^{\mathrm{TP}}=\frac{9}{5},y_h^{\mathrm{TP}}=\frac{77}{40},y_g^{\mathrm{TP}}=\frac{4}{5}.$$

$$⟨⟨E^n(\beta )⟩⟩=(-1{)}^n\frac{{\partial }_{\beta }^n\ln Z(\beta )}{N}.$$

$$\mathrm{Re}[{\beta }^{(1)}-\beta ]$$

$$|{\beta }^{(1)}-\beta {|}^2$$

$$\mathrm{Im}[H^{(1)}]\approx \pm \frac{1}{\beta }\sqrt{2n(2n+1)\left|\frac{⟨⟨M^{2n}(0)⟩⟩}{⟨⟨M^{2(n+1)}(0)⟩⟩}\right|}.$$

$$⟨⟨E^n⟩⟩\sim L^{-d+n{y}_t^{\mathrm{IM}}}.$$

$$E_{\Delta }=\Delta \sum _i{\sigma }_i^2.$$

$$G_L\left(r_i^{(j)}\right)=a_1{r}_i^{(j)}{L}^{a_2}+a_3,$$

$$G_L({\beta }_i^{(1)})=a_1{\beta }_i^{(1)}{L}^{\frac{d}{y_t}}+a_3,$$

$$G_L({\Delta }_i^{(1)})=a_1{\Delta }_i^{(1)}{L}^{\frac{d}{y_t}}+a_3.$$

Full text

††: Journal of Statistical Mechanics

1 Introduction

In this work, we build upon the extensive studies of the critical and tricritical properties of the two-dimensional (2D) Blume-Capel model [1, 2]. This model generalizes the Ising model by allowing spins to take values $\sigma_{i}=\pm 1,0$ instead of just $\pm 1$. Along with the usual nearest-neighbor ferromagnetic exchange interaction $J$ and coupling to an external magnetic field $H$, the Blume-Capel model introduces a crystal-field interaction $\Delta$. The latter acts as a chemical potential that controls the density of zero-spin states: larger values of $\Delta$ favour configurations where $\sigma_{i}=0$ is more prevalent. The introduction of the uncoupled zero-spin state, like invisible states in the Potts model [3] that only contribute to the entropy, enriches the phase diagram, giving rise to both second-order and first-order phase transitions that meet at a tricritical point, see figure 1. This added complexity makes the Blume-Capel model a particularly valuable framework for investigating multicritical phenomena in statistical physics [4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15].

Since it was first proposed, the Blume-Capel model has been investigated through mean-field theory, perturbative expansions, and numerical simulations on a variety of lattices, primarily in two and three dimensions; see, e.g., [6, 16, 17, 18]. The vast majority of studies have focused on the square lattice, employing a broad range of methods. These include real-space renormalisation [19], Monte Carlo simulations and Monte Carlo renormalisation-group approaches [20, 21, 22, 23, 24, 25, 26, 27, 28, 29, 30, 31], $\epsilon$-expansion techniques [32, 33, 34, 35], high- and low-temperature series expansions [36, 37, 38], and transfer-matrix calculations [39, 40, 25, 41, 42].

Here, we extend our analysis of partition function zeros–a technique we previously showed to be highly effective for probing critical properties of the 2D Blume-Capel model, even with very modest system sizes [43] (this constraint was also applied in three dimensions in [44, 45] and on complete graphs in [46]). Our strategy, detailed in [47], remains unchanged: restrict simulations to extremely small lattices and extract maximal information from the zeros of the partition function after analytic continuation into the complex plane of the magnetic field (Lee–Yang zeros), temperature (Fisher zeros), or the crystal-field parameter.

Throughout this study, we adopt the same notation as in the referenced article [43], to which we refer readers for detailed definitions. The Hamiltonian of the Blume-Capel model is given by [1, 2]

[TABLE]

and the corresponding partition function for a system of $N$ spins in $d$ dimensions on a simple hypercubic lattice can be expressed in terms of the density of states $\rho(E_{J},M,E_{\Delta})$ as

[TABLE]

where $x=e^{-\beta J}$, $y=e^{\beta H}$, and $z=e^{-\beta\Delta}$.

To illustrate the method, let us consider the Lee-Yang approach. We fix the temperature and crystal field at their critical values, $\beta_{\rm c}$ and $\Delta_{\rm c}$—for instance, along the second-order transition line or at the tricritical point—and examine the dependence of the partition function on the magnetic field. For a finite system, factoring out the smallest power of $y=e^{-\beta_{\rm c}H}$ transforms $Z_{N}(\beta_{\rm c}J,\,\beta_{\rm c}H,\,\beta_{\rm c}\Delta_{\rm c})$ into a polynomial in $y$. Extending $H$ to complex values, this polynomial has roots $H^{(j)}$, ordered by increasing imaginary part. At each zero $H^{(j)}$, the free-energy density, $-\beta_{\rm c}^{-1}\ln Z_{N}(\beta_{\rm c}J,\,H^{(j)},\,\beta_{\rm c}\Delta_{\rm c})$, is singular. In the thermodynamic limit, these singularities manifest at real parameter values: the lowest zeros condense onto a continuous curve that pinches the real axis precisely at the critical point, $H_{\rm c}=0$. In the Lee-Yang case, symmetry under $M\to-M$ further constrains all zeros to lie on the imaginary axis (in the Lee–Yang theorem, the $\mathbb{Z}_{2}$ symmetry ($M\to-M$) ensures that zeros come in reciprocal pairs, while ferromagnetic interactions further constrain all zeros of the partition function–in the fugacity $y=e^{-\beta H}$–to lie exactly on the unit circle, i.e., on the imaginary axis in the complex magnetic field plane. The validity of the theorem can nevertheless be compromised in some cases [48]).

At a critical (or tricritical) point, the finite-size scaling of the lowest zeros follows, to leading order, a power law:

[TABLE]

where $y_{h}$ is the magnetic renormalisation-group exponent for the universality class considered.

The extension to Fisher zeros and crystal-field zeros requires accounting for the fact that the real parts converge to non-zero values of the corresponding parameters:

[TABLE]

Here, $y_{t}$ and $y_{g}$ are associated with even scaling fields (the smaller of which corresponds to a correction to scaling), while $y_{h}$ corresponds to the magnetic sector, i.e., an odd scaling field under the $\mathbb{Z}_{2}$ transformation $\forall i,\sigma_{i}\to-\sigma_{i}$. The renormalisation-group eigenvalues along the second-order line of the Ising model (IM)’s universality class are denoted $y_{i}^{\rm IM}$, whereas those at the tricritical point (TP) are denoted $y_{i}^{\rm TP}$. The known values in two dimensions are:

[TABLE]

These values serve as references, and we will demonstrate that they can be accurately recovered even using simulations of very small system sizes.

In the following, we explore various regions of the 2D Blume-Capel phase diagram, shown in figure 1. For clarity, we will return to this phase diagram at several points throughout the article, each time explicitly highlighting the region under investigation.

2 Numerical observations along the transition lines of the Blume-Capel model

In this section, we visualize representative spin configurations at various points in the phase diagram.

Along the line $\Delta=0$ (figure 2), the system exhibits a clear evolution with temperature. At low temperatures, it is in a ferromagnetic, ordered phase dominated by either $+1$ spins (red online) or $-1$ spins (blue). As the temperature approaches the Ising transition point, clusters of red and blue spins become increasingly fragmented and fractal-like. Neither spin orientation dominates, and the overall magnetization approaches zero. Above the transition, in the paramagnetic regime, small red and blue clusters appear and average out, while zero spins (white) remain rare throughout this region.

As the crystal-field parameter $\Delta$ increases, white (zero-spin) sites become progressively more common, and white clusters start to appear in the configurations. Near the tricritical point (figure 3), the ordered phase can be predominantly red, blue, or white, with white clusters often dominating, particularly close to the tricritical temperature. Even in the surrounding paramagnetic region, white spins remain abundant, highlighting the growing influence of the crystal field.

In the first-order transition regime (figure 4), the change is much sharper. In the paramagnetic phase, zero spins overwhelmingly occupy the lattice, whereas in the ordered phase, red or blue spins form large, well-defined clusters that dominate the system. This dramatic change in spin composition and structure illustrates the distinct character of first-order transitions. However, since the simulations are performed on finite lattices, the transitions appear slightly rounded and are not as abrupt as in the thermodynamic limit.

3 Ising universality class

3.1 Renormalisation-group scaling dimensions along the critical line

Although the universality class of the 2D Blume-Capel model has been extensively studied–as illustrated by the phase diagram presented in the introduction–we begin here by analyzing the scaling behaviour of three types of zeros: Fisher zeros in the complex-temperature plane, zeros in the crystal-field plane, and Lee-Yang zeros in the complex magnetic-field plane. In this section, we focus solely on the finite-size scaling of the imaginary part of the first zero at four distinct points along the critical line.

The simulations cover system sizes $L=16$ to $64$. This allows us to (1) confirm that all four points belong to the same universality class–namely, the 2D Ising class, with $y_{t}^{\rm IM}=1$ and $y_{h}^{\rm IM}=15/8$–and (2) validate our zero-based analysis approach, which is developed further throughout this work.

In figure 5 (left panel), we show the locations in the phase diagram where the simulations were performed, with coordinates taken from the literature. The three panels below illustrate the finite-size scaling of the imaginary parts of the first Fisher zero, the first crystal-field zero, and the first Lee–Yang zero. The first two confirm $y_{t}^{\rm IM}\approx 1.000$, while the third confirms $y_{h}^{\rm IM}\approx 1.875$, in excellent agreement with the 2D Ising universality class. Detailed numerical data are reported in table 1. We note that the simulation at $\Delta=1$ yields slightly lower-quality fits, which we attribute to the less precise knowledge of the critical line at this point; accordingly, we tuned $\Delta$ to $0.991$ for this simulation.

3.2 Role of off-critical observations

In this section, we revisit a key insight by Deger et al [49, 50], namely that critical properties can be accurately extracted without restricting simulations to the immediate vicinity of the transition. This is noteworthy, since finite-size studies typically focus on points of maximal divergence–such as susceptibility peaks–close to the critical point.

To illustrate this, we return to simulations at $\Delta=0$ and examine the first Fisher zero across a broad temperature range, employing histogram reweighting to scan values of $\beta$. For each real $\beta$, we plot $\operatorname{Re}(\beta-\beta^{(1)})$ and identify the temperature where this quantity crosses zero. Around this crossing, both the real and imaginary parts of the first Fisher zero, $\operatorname{Re}(\beta^{(1)})$ and $\operatorname{Im}(\beta^{(1)})$, exhibit plateaus, remaining essentially constant. This behaviour is illustrated in figure 6, where the left panel shows the wide temperature interval explored through histogram reweighting in the phase diagram.

The emergence of such plateaus explains why finite-size scaling of the Fisher zero can be reliably performed over a wide window, $0.80\,\beta_{\rm c}<\beta<1.20\,\beta_{\rm c}$, using the cumulant method, even for very small systems ($L=5-8$). This robustness allows critical exponents to be extracted from data collected far from the nominal transition, thereby confirming the conclusion of Deger et al [49, 50]: critical exponents can be accurately determined even when the system is away from the phase transition. An application of this result is presented in figure 7, where we fit $\operatorname{Re}(\beta^{(1)})$ versus the system size $L$ to equation (4) at two different temperatures, $\beta_{\rm c}$ and $\beta_{r}^{(1)}$. Both fits extrapolate to the same critical value $\beta_{\rm c}$, with a relative difference smaller than $0.009\%$.

Finally, it is also interesting to look at the behaviour of the Lee-Yang zeros with respect to the temperature. Figure 8 shows that for $\beta<\beta_{\rm c}$ they are not close to the real axis suggesting that there is a Lee-Yang gap even in the thermodynamic limit, while for $\beta>\beta_{\rm c}$ it seems that the zeros are getting closer and closer to the real axis.

3.3 Cumulant method

In [49, 50, 51, 52], Deger and co-authors demonstrated that the first zero of the partition function can be numerically extracted from high-order cumulants of observables–either magnetization or energy-denoted $\langle\!\langle M^{n}\rangle\!\rangle$ or $\langle\!\langle E^{n}\rangle\!\rangle$. They showed that, in the limit of large cumulant order $n$, these expressions converge extremely rapidly.

Fluctuations of a quantity, e.g., $E$, with respect to a control parameter, e.g., $T$ or $\beta$, are defined as derivatives of $\ln Z$:

[TABLE]

This relation leads to closed-form expressions for the zeros. For the Fisher zeros, one obtains:

[TABLE]

Extension to other zeros is obvious. Formulas (9) and (10) hold in the asymptotic regime $n\gg 1$ and enable determination of both the real and imaginary parts of the first zero from four consecutive cumulants of $E$, measured at a fixed $\beta$. As emphasized above, this fixed value may lie within a finite window–often broad, especially in the Fisher case–around the critical point $\beta_{\rm c}$.

For the Lee-Yang zeros, odd cumulants $\langle\!\langle M^{2n+1}\rangle\!\rangle$ vanish due to symmetry. Equation (10) then simplifies to:

[TABLE]

Further discussion can be found in [47, 43].

High-order cumulants are numerically delicate: they tend to vanish rapidly or diverge strongly. It is therefore essential to assess their stability before using them to infer zeros. An illustration of energy cumulants up to $n=9$ is shown in figure 9. Two distinct behaviours appear depending on parity: even-order cumulants display a minimum or maximum near $T_{\rm c}$, whereas odd-order cumulants vanish at $T_{\rm c}$. In this chapter, cumulants are obtained from the central moments of the energy, while in [45] they were computed directly from magnetization averages. Small discrepancies can therefore arise; however, the central-moment method consistently yields more stable results.

Since cumulants can be generated up to large $n$, their internal consistency can be checked via finite-size scaling at criticality. The expected scaling reads:

[TABLE]

In two dimensions, this reduces to $\langle\!\langle E^{n}\rangle\!\rangle\sim L^{\,n-2}$. On a log-log scale, this implies a straight line with slope $n-2$, as confirmed in figure 10.

Separate fits for even and odd cumulants show acceptable convergence even for $n\geq 1$, with particularly precise agreement for even cumulants when $n\geq 2$.

The close match between measured slopes and the theoretical expectation validates the cumulant approach: high-order cumulants from relatively small lattices already reproduce the correct scaling behaviour. Their rapid convergence provides accurate estimates of the first Fisher zero, while exponential convergence with order has been confirmed in other models (Ising, Potts), attesting to the robustness of the cumulant method. This has been exploited in tables 5 and 6 in [43] for both Fisher and Lee-Yang zeros.

4 Tricritical point singularities

4.1 Off-critical scaling

In this section, we investigate the stability range of the first zero near the tricritical point. In this regime, it is more appropriate to analyse crystal-field zeros rather than Fisher zeros, since the transition line at $(T_{\rm t},\Delta_{\rm t})$ is oriented nearly perpendicular to the crystal-field axis in the phase diagram (see figure 11). Accordingly, the crystal field $\Delta$ serves as a natural parameter for probing the transition. In the right panel of the figure, we first locate the region where the real part $\operatorname{Re}(\Delta-\Delta^{(1)})$ vanishes, and then examine both $\operatorname{Re}(\Delta^{(1)})$ and $\operatorname{Im}(\Delta^{(1)})$ as functions of $\Delta$. Although the stability observed here is less pronounced than in the case of the Fisher zero at $\Delta_{\rm c}=0$, the first crystal-field zero nevertheless remains remarkably stable within the interval $0.95,\Delta_{\rm t}\leq\Delta\leq 1.05,\Delta_{\rm t}$, which represents a fairly broad window.

4.2 Cumulant scaling

As for the case $\Delta_{\rm c}=0$, we now focus on the tricritical point by extracting the thermal exponent from the scaling of the cumulants of the $\Delta$-contribution to the energy, defined as

[TABLE]

The results of the cumulant scaling are presented in figure 12. Compared to the critical case, the exponents–both odd and even–deviate more noticeably from the expected tricritical values, suggesting that either corrections to scaling must be included or larger system sizes are required. The goodness of fit, measured by the reduced chi-square, is approximately $\chi^{2}_{s}\approx 2$ in each case.

4.3 First crystal-field zero scaling

We now turn to the finite-size scaling analysis of the real and imaginary parts of the first crystal-field zero, determined either directly from the partition function or via the cumulant expansion. Figure 13 presents the scaling of the imaginary part $\operatorname{Im}(\Delta^{(1)})$, evaluated at both $\Delta=\Delta_{\rm t}$ and $\Delta=\Delta_{r}^{(1)}$, and compared against the partition function method for system sizes $L=16$–$64$. The agreement is excellent: at $\Delta=\Delta_{r}^{(1)}$ we obtain $y_{t}=1.8035(22)$ with $\chi^{2}_{s}=2.2$, in almost perfect correspondence with the partition function result $y_{t}=1.8034(21)$. At $\Delta=\Delta_{\rm t}$ the fit yields $y_{t}=1.8091(22)$ with $\chi^{2}_{s}=1.51$, which is likewise consistent within statistical uncertainties.

For the real part, the results are presented in figure 14. Unlike the imaginary component, fitting the scaling of $\operatorname{Re}(\Delta^{(1)})$ with $\Delta_{\rm t}$ fixed yields a poor estimate of the exponent $y_{t}$. However, when the tricritical exponent is fixed to $y_{t}^{\rm TP}=1.80$, the shift behaviour of the real part provides an excellent determination of $\Delta_{\rm t}$. Specifically, we obtain $\Delta_{\rm t}=1.96590(1)$ with $\chi^{2}_{s}=1.94$ at $\Delta=\Delta_{\rm t}$, and $\Delta_{\rm t}=1.96590(1)$ with $\chi^{2}_{s}=1.95$ at $\Delta=\Delta_{r}^{(1)}$, in remarkable agreement between the two approaches.

4.4 First Lee-Yang zero scaling

The finite-size scaling analysis of the first Lee–Yang zero, obtained via the cumulant method up to order $n=3$, has been extensively studied in [43]. In the present work, we extend this analysis to $n=6$ and demonstrate that the convergence continues to improve with higher-order cumulants; see figure 15.

5 Density of the first zeros

The density of partition function zeros provides valuable insight into both the nature and strength of phase transitions. This line of investigation was systematically pursued by Janke and Kenna [53], who demonstrated its effectiveness across a wide range of systems, even when only relatively small lattice sizes were available.

Subsequently, Alves et al [54] compared this classification scheme with an alternative approach introduced by Borrmann et al [55], which relies on the linear scaling of the limiting density of zeros. Their analysis showed that finite-size scaling of modest system sizes is already sufficient to distinguish the order of the transition in the $4$- and $5$-state Potts models. More recently, the density of Lee-Yang zeros in the three-dimensional Blume–Capel model was investigated in [45].

In this framework, the cumulative distribution of zeros can be written as

[TABLE]

which encapsulates key information about the phase transition. At the transition, the constant term $a_{3}$ must vanish, while the exponent $a_{2}$ determines the nature of the transition. For a first-order transition one expects $a_{2}\approx 1$, with the prefactor $a_{1}$ directly related to the latent heat $\Delta e$ (for Fisher zeros) or to the spontaneous magnetisation (for Lee–Yang zeros). By contrast, for a second-order transition one finds $a_{2}=2-\alpha=d/y_{t}$ for Fisher zeros, and $a_{2}=d/y_{h}$ for Lee–Yang zeros. Here, $r^{(j)}_{i}$ denotes the imaginary part of the $j$-th zero.

At the critical point $\Delta_{\rm c}=0$, the finite-size scaling of the Fisher zero density is analysed through

[TABLE]

whereas at the tricritical point the analysis focuses on the scaling of the crystal-field zeros:

[TABLE]

This latter case was first investigated in [47]. Figure 16 shows the Fisher zero density at the critical point (left panel) and the crystal-field zero density at the tricritical point (right panel), for system sizes $L=16-64$. In both cases, we perform two separate fits. In the first, $a_{3}$ is kept as a free parameter to test whether it vanishes–a necessary condition for the presence of a phase transition. For both the critical and tricritical analyses, the fitted values of $a_{3}$ are indeed consistent with zero. Having established this, we then set $a_{3}=0$ in a second round of fits, which allows us to extract the critical exponents with improved accuracy.

At the critical point, the fit yields $d/y_{t}=2.0053(35)$, in excellent agreement with the exact value $d/y_{t}^{\text{IM}}=2$. At the tricritical point, where the theoretical prediction is $d/y_{t}^{\text{tri}}\approx 1.111$, we obtain the numerical estimate $d/y_{t}=1.1090(13)$, again in strikingly good agreement.

At the tricritical point, we extend previous analyses by investigating crossover effects in its vicinity ($\Delta_{\rm t}=1.966$, $T_{\rm t}=0.608$) through the behaviour of the Lee–Yang zeros. As demonstrated in [43, 45], Lee–Yang zeros are highly sensitive to variations in the external parameters, making them a powerful probe of the three distinct regimes expected in the Blume–Capel phase diagram near tricriticality. To explore this sensitivity, we analyse two nearby temperatures–one slightly below the tricritical point, $T=0.6075$, and one slightly above, $T=0.6083$–using the same parameter values employed in the analysis of figure 8 in [43]. In all cases, the crystal-field parameter is fixed to its tricritical value $\Delta_{\rm t}$.

As in the previous analyses, the fitting procedure starts by checking whether $a_{3}$ vanishes, following the approach illustrated in figure 17 at criticality. Once this condition is confirmed, we fix $a_{3}=0$ to obtain a more precise estimate of the critical exponent. At the critical point, this procedure yields $d/y_{h}=1.0637(2)$, in excellent agreement with the Ising universality class prediction, $d/y_{h}^{\text{IM}}\approx 1.0666$. At the tricritical point, the estimate $d/y_{h}=1.0393(32)$ closely matches the theoretical value $d/y_{h}^{\text{tri}}\approx 1.03896$. Examining nearby temperatures reveals crossover effects: for $T=0.6075$, the fitted exponent $d/y_{h}=1.0041(25)$ clearly indicates a first-order transition, consistent with the expected value $d/y_{h}=1$. Conversely, at $T=0.6083$, the extracted value $d/y_{h}=1.0650(38)$ aligns remarkably well with the Ising universality class, highlighting the sensitivity of the Lee–Yang zeros to the external parameters and the transition between different scaling regimes near tricriticality.

Finally, we extend the analysis beyond the leading zero to include higher-order zeros, which also carry valuable information on the critical behaviour of the system [56]. While the cumulant method provides the leading zero, the higher-order zeros were extracted by tracking sign changes of the partition function in the complex crystal-field plane. Our study focuses on the tricritical point for system sizes $L=8$, $32$, $56$, and $64$, where the first three crystal-field zeros were determined. Their finite-size scaling behaviour is presented in the left panel of figure 18, alongside the corresponding density of crystal-field zeros in the right panel. Because the zeros are correlated at each lattice size, independent fits were performed for each zero.

The analysis of the second and third crystal-field zeros provides highly consistent estimates of the tricritical exponent, yielding $y_{t}^{\rm TP}=1.7986(7)$ from $\Delta_{i}^{(2)}$ and $y_{t}^{\rm TP}=1.7981(5)$ from $\Delta_{i}^{(3)}$. Likewise, the corresponding density-of-zeros analysis, shown in figure 18 (right panel), produces results in excellent agreement with the literature value $d/y_{t}^{\rm TP}\approx 1.11$, specifically $d/y_{t}=1.1120(4)$ for the second zero and $d/y_{t}=1.1123(3)$ for the third zero.

6 Conclusions

In this work, we have extended the study of partition function zeros as a tool for analyzing critical and tricritical behaviour in the two-dimensional Blume–Capel model. Building on our previous studies [43, 47], we examined Lee–Yang zeros in the complex magnetic-field plane, Fisher zeros in the complex-temperature plane, and crystal-field zeros, systematically exploring their finite-size scaling properties along both the second-order critical line and near the tricritical point. Our results show that accurate estimates of critical exponents can be obtained even from simulations performed away from the nominal transition points, thereby confirming and generalizing the findings of Deger et al [49, 50]. Notably, reliable results are achievable with surprisingly modest lattice sizes, significantly reducing the computational effort required.

Overall, our study highlights the versatility and efficiency of the partition function zero approach, especially when combined with the cumulant method. By showing that reliable critical information can be extracted from small systems and off-critical simulations, this work not only confirms theoretical predictions for the Blume–Capel model but also provides a framework for computationally efficient investigations of other complex systems.

Acknowledgments

The work of N G F was supported by the Engineering and Physical Sciences Research Council (grant EP/X026116/1 is acknowledged).