We expect the ISAT transition to be somewhat shifted away from \((\omega,\tau)=(1,3)\) due to considering a finite-size ensemble. If one considers the vertical line \(\omega=1\), one notices that the maximal eigenvalue peaks at a value of \(\tau\) somewhat greater than \(3\). We remind that the phase transition for \(\omega=1\) is a strong second-order phase transition where the specific heat diverges with an exponent \(0.68(5)\). For \(\omega<1\), there exists a line of even stronger peaks that join with the peak at \(\omega=1\), from which one can infer that there exists a strong phase transition on varying \(\tau\) for each \(\omega<1\). This is borne out by finite-size scaling analysis. Specifically, we have studied the model when \(\omega=0.5\): for this value of \(\omega\) we have simulated the model using a one-parameter flatPERM algorithm up to length \(n = 1024\), running \(7.9 \cdot 10^6\) iterations, and collecting \(2.7 \cdot 10^{10}\) samples at the maximum length. We find that the specific heat divergence is commensurate with a first-order transition with a linear divergence. To test this assumption of a first-order transition, we consider the distribution of the number of contacts for various values of \(\tau\) near the peak of the specific heat. Figure \ref{fig:double_peaks} shows a clear bimodal distribution, confirming the first-order character of the transition.

It is hence likely that there exists a line of first-order phase transitions at values of \(\tau\) near \(3\) for each value of \(\omega\lesssim 1\). We will return to the question of whether the line of first-order transitions extends all the way to \(\omega=1\) below.

Returning to our two-parameter data, let us first note that when \(\tau=0\) our model is the ISAW model and hence there is a weak \(\theta\)-like transition at \((\omega, \tau) = (\omega_c,0)\) with \(\omega_c \sim 1.94\) (see \cite{Caracciolo_2011} and references therein for recent estimates of \(\omega_c\)), which is reflected in the density plot by a broad peak near \(\omega = 2\) on varying \(\omega\) when \(\tau=0\). A line of such weak peaks extends to larger values of \(\tau\). When \(\tau=1\), the model becomes the previously studied INNSAT \cite{Bedini_2013_INNSAT}, which also demonstrated a weak \(\theta\)-like transition with exponent estimates encompassing ISAW values. We thus conjecture that the entire line lies in the \(\theta\)-universality class.

Now if the suggestion that the ISAT collapse corresponds to an infinite order multi-critical point \cite{Nahum_2013} is correct, it is then natural and simplest to conjecture that the line of first-order transitions meets with the line of \(\theta\)-like transitions at that point.

We obtain an indication of where the two lines might join by considering the ratio between the two eigenvalues of the covariance matrix \(H_{ij} = \partial^2_{ij} \log Z_n\) where \(i, j \in \{\omega, \tau\}\). This is based on the heuristic argument that the component \(\tau\tau\) (respectively \(\omega\omega\)) is the specific heat associated to multiply visited sites (nearest-neighbour contacts) and will dominate the spectrum when the transition is driven by this interaction. When the two eigenvalues coincide is then argued to indicate the presence of a higher order critical point. The density plot of the eigenvalue ratio is shown on Figure \ref{fig:eigenvalue_plot}. One clearly observes a unique point close to the ISAT collapse point, where the two eigenvalues have the same magnitude.