Generalized Bose-Hubbard Hamiltonians exhibiting a complete non-Hermitian degeneracy

TL;DR

This paper introduces a method to construct a broader class of non-Hermitian Bose-Hubbard Hamiltonians with complete exceptional point degeneracies, overcoming previous technical obstacles and enabling detailed spectral analysis.

Contribution

It presents a novel approach to circumvent technical challenges in modeling non-Hermitian Bose-Hubbard systems with EPN degeneracies, expanding the class of eligible models.

Findings

Broader class of Hamiltonians with EPN degeneracies identified

Method demonstrated through symbolic computations for larger N

Spectral analysis near and far from EPN regimes performed

Abstract

The conventional toy-model constructions of phase diagrams often use various versions of the standard Hermitian Bose-Hubbard Hamiltonians $H$. These studies were recently extended to cover several non-Hermitian PT-symmetric versions of the model. A key technical merit of the extensions has been found in the existence of Kato's phase-transition-related exceptional-point spectral degeneracies of the $N$th order (EPN). In this setting, an ill-conditioned nature of $H$ near its EPN singularity represented a decisive technical obstacle. In our paper we show how the obstacle can be circumvented, leading to a broader class of the eligible models with tridiagonal and symmetric complex Hamiltonians $H^{(N)}(z)$ exhibiting the EPN degeneracies. The implementation of the method (requiring the use of computer-assisted symbolic manipulations at the larger integers $N$) is illustrated by sampling the real or complex spectra of energies near as well as far from the EPN dynamical regime.

Tables

Table 1. Table 1: Numerical confirmation of the occurrence of EP6 in the non-BH Hamiltonian H a ( 6 ) ​ ( z ) subscript superscript 𝐻 6 𝑎 𝑧 H^{(6)}_{a}(z) .

precision $p$	min ${\varrho }_{mn}$	max ${\varrho }_{mn}$
10	7.2 ${10}^{-6}$	2.9 ${10}^{-5}$
20	6.4 ${10}^{-7}$	2.6 ${10}^{-6}$
30	2.7 ${10}^{-10}$	1.1 ${10}^{-9}$
40	1.4 ${10}^{-13}$	5.5 ${10}^{-13}$

Equations

$$H^{(2)}_{(BH)}(\gamma)=\left[\begin{array}[]{cc}-i{\it\gamma}&1\\ \vskip 6.0pt plus 2.0pt minus 2.0pt\cr 1&i{\it\gamma}\end{array}\right]\,,\ \ \ \ H^{(3)}_{(BH)}(\gamma)=\left[\begin{array}[]{ccc}-2\,i\gamma&\sqrt{2}&0\\ \vskip 6.0pt plus 2.0pt minus 2.0pt\cr\sqrt{2}&0&\sqrt{2}\\ \vskip 6.0pt plus 2.0pt minus 2.0pt\cr 0&\sqrt{2}&2\,i\gamma\end{array}\right]\,,\ldots$$

$$E_k^{(BH)}(\gamma )=(1-N+2k)\sqrt{1-{\gamma }^2},k=0,1,\dots ,N-1.$$

$$H^{(N)}(\pm 1)\sim J^{(N)}(E)=\left[\begin{array}[]{ccccc}E&1&0&\ldots&0\\ {}0&E&1&\ddots&\vdots\\ {}0&0&E&\ddots&0\\ {}\vdots&\ddots&\ddots&\ddots&1\\ {}0&\ldots&0&0&E\end{array}\right]\,$$

$$H^{(4)}(z,A,B)=\left[\begin{array}[]{cccc}-3\,i\sqrt{z}&\sqrt{B}&0&0\\ \vskip 6.0pt plus 2.0pt minus 2.0pt\cr\sqrt{B}&-i\sqrt{z}&\sqrt{A}&0\\ \vskip 6.0pt plus 2.0pt minus 2.0pt\cr 0&\sqrt{A}&i\sqrt{z}&\sqrt{B}\\ \vskip 6.0pt plus 2.0pt minus 2.0pt\cr 0&0&\sqrt{B}&3\,i\sqrt{z}\end{array}\right]\,$$

$$H^{(5)}(z,A,B)=\left[\begin{array}[]{ccccc}-4\,i\sqrt{z}&\sqrt{B}&0&0&0\\ \vskip 6.0pt plus 2.0pt minus 2.0pt\cr\sqrt{B}&-2\,i\sqrt{z}&\sqrt{A}&0&0\\ \vskip 6.0pt plus 2.0pt minus 2.0pt\cr 0&\sqrt{A}&0&\sqrt{A}&0\\ \vskip 6.0pt plus 2.0pt minus 2.0pt\cr 0&0&\sqrt{A}&2\,i\sqrt{z}&\sqrt{B}\\ \vskip 6.0pt plus 2.0pt minus 2.0pt\cr 0&0&0&\sqrt{B}&4\,i\sqrt{z}\end{array}\right]\,$$

$$A_{(BH)}^{(4)}=4,B_{(BH)}^{(4)}=3,A_{(BH)}^{(5)}=6,B_{(BH)}^{(5)}=4$$

$$A_{(non-BH)}^{(4)}=64,B_{(non-BH)}^{(4)}=-27,A_{(non-BH)}^{(5)}=-54,B_{(non-BH)}^{(5)}=64.$$

$$H^{(6)}({{z}},A,B,C)=\left[\begin{array}[]{cccccc}-5\,i{{\sqrt{z}}}&\sqrt{C}&0&0&0&0\\ \vskip 6.0pt plus 2.0pt minus 2.0pt\cr\sqrt{C}&-3\,i{{\sqrt{z}}}&\sqrt{B}&0&0&0\\ \vskip 6.0pt plus 2.0pt minus 2.0pt\cr 0&\sqrt{B}&-i{{\sqrt{z}}}&\sqrt{A}&0&0\\ \vskip 6.0pt plus 2.0pt minus 2.0pt\cr 0&0&\sqrt{A}&i{{\sqrt{z}}}&\sqrt{B}&0\\ \vskip 6.0pt plus 2.0pt minus 2.0pt\cr 0&0&0&\sqrt{B}&3\,i{{\sqrt{z}}}&\sqrt{C}\\ \vskip 6.0pt plus 2.0pt minus 2.0pt\cr 0&0&0&0&\sqrt{C}&5\,i{{\sqrt{z}}}\end{array}\right]$$

$$s^3+\left(-A+35z-2C-2B\right)s^2+\left(B^2+2AC+28Cz+259z^2+C^2-34Az+2BC-44Bz\right)s+$$

$$+30C{z}^2+C^{2}z-A{C}^2+25B^{2}z-225A{z}^2+150B{z}^2+225z^3+10BCz-30ACz=0$$

$$\underset{z\to z^{(EPN)}}{\lim }{E}_{\pm k}(z,A^{(N)},B^{(N)},\dots ,Z^{(N)})=0,k=1,2,\dots ,J$$

$$s^J+P_1(A,B,\dots ,Z)s^{J-1}+\dots +P_{J-1}(A,B,\dots ,Z)s+P_J(A,B,\dots ,Z)=0$$

$$P_m(A^{(N)},B^{(N)},\dots ,Z^{(N)})=0,m=1,2,\dots ,J.$$

$$(A_{(BH)}^{(6)},B_{(BH)}^{(6)},C_{(BH)}^{(6)})=(9,8,5).$$

$$416{\xi }^4+20909{\xi }^3+22505{\xi }^2+28734375\xi -48828125=0.$$

$${\xi }_a=-65.80360706245132477179785808904814860530$$

$${\xi }_b=1.693394621288372898472626362413820064872.$$

$$A_a^{(6)}=673.7717872,B_a^{(6)}=-253.5822865,C_a^{(6)}=-65.80360706.$$

$$A_b^{(6)}=107.5337579,B_b^{(6)}=-37.96027355,C_b^{(6)}=1.693394621.$$

$$H^{(7)}(z,A,B,C)=\left[\begin{array}[]{ccccccc}-6\,i{\sqrt{z}}&\sqrt{C}&0&0&0&0&0\\ \vskip 6.0pt plus 2.0pt minus 2.0pt\cr\sqrt{C}&-4\,i{\sqrt{z}}&\sqrt{B}&0&0&0&0\\ \vskip 6.0pt plus 2.0pt minus 2.0pt\cr 0&\sqrt{B}&-2\,i{\sqrt{z}}&\sqrt{A}&0&0&0\\ \vskip 6.0pt plus 2.0pt minus 2.0pt\cr 0&0&\sqrt{A}&0&\sqrt{A}&0&0\\ \vskip 6.0pt plus 2.0pt minus 2.0pt\cr 0&0&0&\sqrt{A}&2\,i{\sqrt{z}}&\sqrt{B}&0\\ \vskip 6.0pt plus 2.0pt minus 2.0pt\cr 0&0&0&0&\sqrt{B}&4\,i{\sqrt{z}}&\sqrt{C}\\ \vskip 6.0pt plus 2.0pt minus 2.0pt\cr 0&0&0&0&0&\sqrt{C}&6\,i{\sqrt{z}}\end{array}\right]\,.$$

$$s^3+\left(-2B-2A-2C+56\right)s^2+\left(2AB+B^2+2BC+C^2+784-56B-104A+40C+4AC\right)s+$$

$$+2304-2BAC+24BC-96AC+72AB-2A{C}^2-1152A+576B+36B^2+4C^2+192C=0.$$

$$A_{(BH)}^{(7)}=12,B_{(BH)}^{(7)}=10,C_{(BH)}^{(7)}=6,$$

$$B_{\alpha ,\beta }^{(7)}=76-6r_{\alpha ,\beta },C_{\alpha ,\beta }^{(7)}=6r_{\alpha ,\beta },$$

$$r_{\alpha }=-59+10\sqrt{34}\approx -0.69048,r_{\beta }=-59-10\sqrt{34}\approx -117.3095$$

$$r_{\gamma }=3/2+1/2\sqrt{21}\approx 3.7912878,r_{\delta }=3/2-1/2\sqrt{21}\approx -0.7912878$$

$$A_{\gamma ,\delta }^{(7)}=36r_{\gamma ,\delta },B_{\gamma ,\delta }^{(7)}=28-54r_{\gamma ,\delta },C_{\gamma ,\delta }^{(7)}=18r_{\gamma ,\delta }.$$

$$P_1=-A-2D+84-2C-2B=0,$$

$$P_2=1974+2AC+2AD+D^2+50D+4BD+B^2+2CD+C^2-83A-142B+2BC-70C=0,$$

$$P_3=1402C-2CBD+50C^2+12916-A{C}^2+74B^2+108BC+68AC+682D-2CAD-2B^{2}D-$$

Full text

Generalized Bose-Hubbard Hamiltonians exhibiting a complete non-Hermitian degeneracy

Miloslav Znojil

Nuclear Physics Institute of the CAS, Hlavní 130, 250 68 Řež, Czech Republic

e-mail: znojil@ujf.cas.cz

Abstract

The method of construction of the tridiagonal and symmetric complex-matrix Hamiltonians $H^{(N)}(z)$ exhibiting an exceptional-point (EP) degeneracy of the $N$th (i.e., maximal) order at a preselected parameter $z=z^{(EPN)}=1$ is proposed and tested. In general, the implementation of the method requires the use of computer-assisted symbolic manipulations, especially at the larger matrix dimensions $N$. The well known ${\cal PT}-$symmetric $N$-by-$N$-matrix Bose-Hubbard Hamiltonians as well as their recent $N\leq 5$ non-Bose-Hubbard alternatives as obtained as special solutions. Several other $N\geq 6$ non-Bose-Hubbard models $H^{(N)}(z)$ exhibiting the maximal EPN degeneracy are also constructed and analyzed in detail. In particular, their $z-$dependent, real or complex spectra of energies are displayed and discussed near as well as far from $z^{(EPN)}$.

Keywords

non-Hermitian quantum dynamics; multilevel degeneracies; exceptional points of high orders;

1 Introduction

In 2008, Graefe with coauthors [1] pointed out that one of the most user-friendly features of the family

[TABLE]

of the special one-parametric ${\cal PT}-$symmetrized Bose-Hubbard Hamiltonians (relevant, say, for the study of Bose-Einstein condensation) is that their energy spectra are all known in closed form,

[TABLE]

In spite of the manifest non-Hermiticity of the Hamiltonians, the energy levels themselves remain real and non-degenerate (i.e., potentially observable [2]) inside an open and dimension-independent interval of $\gamma\in(-1,1)$, therefore.

One of the related and phenomenologically as well as mathematically most challenging problems emerges in the two limits of $\gamma\to\pm 1$, i.e., on the boundary of the interval of stability and potential observability of the energies. At these values, all of the Hamiltonian matrices (1) cease to be diagonalizable, acquiring the canonical form of the $N$ by $N$ Jordan block,

[TABLE]

with $E=0$. This means that in the language of functional analysis the two extreme values $\gamma=\pm 1$ of the parameter may be interpreted as the Kato’s exceptional points (EP, [3]) of order $N$ (EPN). In these two limiting cases the spectrum becomes fully degenerate. In a broader context of general physics such a specific EP extreme is often called non-Hermitian degeneracy [4]. In the narrower area of unitary quantum systems the phenomenon is better known under the more intuitive, widely accepted name of the spontaneous breakdown of ${\cal PT}$ symmetry [5].

The authors of Ref. [1] emphasized that one of the most important phenomenological consequences of the existence of the EPN $\gamma=\pm 1$ singularities should be seen in the variability of the scenarios in which the degenerate spectrum “unfolds” under the influence of perturbations. In *loc. cit. * the phenomenologically most useful (viz., diagonal-matrix) choice of these perturbations has even been shown tractable analytically, by non-numerical means. A slightly more general version of the perturbative model-building strategy has been discussed in Ref. [6].

The maximality of the non-Hermitian EPN degeneracy (3) seems to follow from the highly specific choice of model (1) and, in particular, from its Lie-algebraic origin and symmetries [1]. In our recent paper [7] we decided to test such a conjecture. We managed to disprove it when, for the sake of simplicity, we reparametrized $\gamma\to\sqrt{z}\in[0,1)$ (so that just a unique EPN value of $z^{(EPN)}=1$ had to be taken into consideration), and when we only admitted the off-diagonal deformations of the model. Even though we restricted our attention just to the first two nontrivial, two-parametric deformations of the respective original Bose-Hubbard Hamiltonians $H^{(4,5)}_{(BH)}(\gamma)$ of Ref. [1], viz., to matrices

[TABLE]

and

[TABLE]

we were able to conclude, due to the comparatively elementary nature of such a generalization, that the maximal non-Hermitian EPN degeneracy of the spectrum can be achieved, at $z^{(EPN)}=1$, not only via the conventional Bose-Hubbard choice of parameters

[TABLE]

but also via its symmetries violating and strongly deformed alternatives

[TABLE]

Unfortunately, the applicability of the construction remained restricted to the latter two examples. The method we used did not seem to admit any immediate extension beyond $N=5$ (cf. section 2 and paragraph 2.1 below). Still, the potential physical relevance of the exceptional points of the higher orders [8] forced us to search for an amendment of the method. The search succeeded, and its results will be reported in what follows.

The overall idea of the innovated construction yielding the new, non-BH Hamiltonians will be explained and, choosing $N=6$, illustrated in section 2. Several characteristic features of its extension beyond $N=6$ will be then discussed in section 3. In subsequent section 4 one of the key technicalities (viz, the necessity of a reliable numerical proof of the maximality of the degeneracy) will finally be identified and resolved. The description of the parameter-dependence of the energy spectra far from the EPN singularity will also be discussed, in separate section 5, in some detail. After a thorough discussion of some terminological and experimental aspects of our results in section 6, the paper will be concluded by a concise summary in section 7.

2 Three-parametric non-BH deformation at $N=6$

2.1 The method of Ref. [7] and its

failure

The feasibility of the constructions of the toy-model Hamiltonians with property (3) using the method of Ref. [7] remained restricted to $N\leq 5$. Indeed, the next, $N=6$ model with the three-parametric candidate

[TABLE]

for the deformed, non-Bose-Hubbard Hamiltonian may be assigned the secular equation

[TABLE]

[TABLE]

where we abbreviated $E^{2}={s}$. Obviously, this equation only defines all of the bound state energies $E_{\pm k}=\pm\sqrt{s_{k}}$, $k=1,2,3$ in terms of the conventional but, in this particular application, prohibitively complicated Cardano formulae for the three roots $s=s_{k}$.

The localization of the EP6 non-Hermitian degeneracy becomes difficult. The main problem is that the Cardano formulae define the real roots as superpositions of complex numbers. This feature of the construction (i.e., the necessity of a guarantee of the exact mutual cancelation of the respective imaginary components of the root) converts the construction of the canonical representation (3) of the Hamiltonians in question into a purely numerical task. Unfortunately, such a numerical task is ill-conditioned. In other words, the explicit canonical-representation approach of Ref. [7] fails. Its applicability remains restricted to the smallest matrix dimensions $N\leq 5$. In what follows, an amended, “implicit”, effective alternative treatment of the problem will be developed and applied at a few sample matrix dimensions $N>5$, therefore.

2.2 Maximal degeneracy condition and the Gröbner-basis

technique

The double symmetry (i.e., the symmetry with respect to its two main diagonals) of our complex matrix (8) is reflected by the up-down symmetry $E_{\pm k}=\pm\sqrt{s_{k}}$ of the spectrum. Thus, at any even matrix dimension $N=2J$, our key requirement of the existence of the complete EPN degeneracy

[TABLE]

(say, at the conveniently chosen $z^{(EPN)}=1$) implies that the secular equation

[TABLE]

must acquire, in the EPN limit, the utterly elementary form $s^{J}=0$. In other words, the $J-$plet of the EPN-compatible parameters $A^{(N)},B^{(N)},\ldots,Z^{(N)}$ must satisfy the $J-$plet of polynomial equations

[TABLE]

The basic tool of the iterated-elimination solution of the similar sets of equations is provided by the construction of the so called Groebner basis. Such a construction is available via the commercially available symbolic-manipulation software – in what follows we shall use MAPLE [9].

2.3 EPN-admitting non-BH Hamiltonians at $N=6$

Once we return to the $J=3$ secular equation (9) the entirely routine application of the Gröbner-basis solvers leads to the five sets of solutions of the three coupled polynomial equations (12). The first, most elementary one just reproduces the well known Bose-Hubbard matrix of Ref. [1],

[TABLE]

In addition one reveals that there exist the other four sets of the EP6-compatible solutions which are all expressed in terms of the four roots $\xi$ of the following auxiliary polynomial equation

[TABLE]

For the sake of brevity we will drop here the discussion of the pair of the complex roots. Thus, we will only consider the remaining two real roots

[TABLE]

and

[TABLE]

In our calculations such a high-precision representation of these “seed” roots appeared necessary, mainly due to the perceivable loss of numerical precision caused by the mutual cancelations between the separate terms in the polynomials in question.

2.3.1 The $a-$subscripted non-BH deformation

In the first, $a-$subscripted case the conventional Gröbner-basis method gave us the desired EP6-compatible parameters with numerical values

[TABLE]

This enables us to complement the $N=4$ and $N=5$ “generalized Bose-Hubbard” models of Ref. [7] by their $a-$subscripted $N=6$ descendant. The explicit form of its one-parametric Hamiltonian $H^{(6)}_{(a)}(z)$ is obtained by the insertion of the EP6-compatible parameters (14) into the general four-parametric matrix (8). The global $z-$dependence of its spectrum is displayed in Fig. 1 and, in a magnified version near the EPN singularity $z^{(EP6)}=1$, in Fig. 2.

2.3.2 The $b-$subscripted non-BH deformation

In the second, $b-$subscripted case we arrive, using the same recipe, at the other, $b-$subscripted set of the EP6-compatible parameters

[TABLE]

Along the same lines as before we also obtain the new one-parametric Hamiltonian matrix $H^{(6)}_{(b)}(z)$ and the $z-$dependence of its spectrum (see Figs. 3 and 4 and also the comments in section 5 below).

Marginally let us add that for the fourth-order polynomial in (13) the exact expression for the root can be obtained, in closed form, via computer-assisted symbolic manipulations. Unfortunately, this type of result is hardly of any use in practical considerations. In contrast, the easy accessibility of the approximate values of the two real roots $\xi_{a,b}$ can still help us to study many relevant properties of the purely numerically constructed EPN-admitting Hamiltonians.

3 General EPN-admitting models

3.1 The deformed-Hamiltonian ansatz at $N=7$

The strategy of construction of the $J-$parametric deformed Bose-Hubbard Hamiltonians with odd $N=2J+1$ remains unchanged. For illustration let us consider the first nontrivial seven-dimensional Hamiltonian ansatz

[TABLE]

We reveal that one of its eigenvalues is constant, identically equal to zero, $E_{0}=0$. The rest of the spectrum, i.e., energies $E_{\pm k}=\pm\sqrt{s_{k}}$ with $k=1,2,3$ are expressed again in terms of the roots of the simplified secular equation

[TABLE]

[TABLE]

The mere inessential modification enables us now to apply the above-outlined method of the search for the EPNs to the model with $N=7$. Indeed, at any odd matrix dimension $N=2J+1$, the condition (10) of the complete EPN degeneracy remains unchanged. The same conventional choice of ${{z}}^{(EPN)}=1$ leads to the formally identical constraint (12). This specifies the $J-$plet of the EPN-compatible parameters $A^{(N)},B^{(N)},\ldots,Z^{(N)}$ at any $J$ and $N=2J+1$, in principle at least.

3.1.1 The $(\alpha,\beta)-$subscripted deformations

Once we recall the Groebner-basis constructive technique we may return to the present $N=7$ upgrade (16) of the $J=3$ secular equation. This enables us to construct the EP7-supporting Hamiltonians, indeed. Besides the known and expected Bose-Hubbard solution with

[TABLE]

we discover the existence of the other two subfamilies of the Hamiltonians. In the first subfamily we arrive at the fixed integer value of $A^{(7)}=-48$. The other two parameters become then defined, by formulae

[TABLE]

in terms of the two respective non-numerical, exactly known roots

[TABLE]

of quadratic equation $r^{2}+118\,r+81=0$.

3.1.2 The $(\gamma,\delta)-$subscripted deformations

In the second half of our exhaustive constructive analysis we would have to analyze the second subfamily of the EP7-responsible coupling constants. These appear to be defined in terms of the other pair of the auxiliary roots

[TABLE]

of the quadratic equation ${r}^{2}-3\,r-3=0$. The explicit closed definitions of the respective EP7-compatible generalized, deformed Bose-Hubbard Hamiltonians then read

[TABLE]

The more detailed analysis of the Hamiltonians defined by these parameters is left to the readers.

3.2 The EPN constructions at $N>7$

3.2.1 The case of even $N=2J=8$

The occurrence of the EP8 degeneracy can be, naturally, studied along the same lines as above. In the first step we have to evaluate the secular polynomial and require that at ${{z}}={{z}}^{(EP8)}=1$ the secular equation degenerates again to its trivial EP8 form $s^{J}=0$. This leads to the quadruplet of coupled polynomial relations

[TABLE]

[TABLE]

[TABLE]

[TABLE]

and

[TABLE]

[TABLE]

[TABLE]

These equations are, expectedly, satisfied by the Bose-Hubard parameters

[TABLE]

It is less elementary to find any other, non-Bose-Hubard solutions. It was necessary to use the algebraic-manipulation software [9]. This enabled us to accept the computer-assisted elimination strategy, to construct the Gröbner basis and to reduce the search for the EP8-supporting matrix elements of $H^{(8)}$ to the purely numerical search for the roots of a single polynomial $R(y)$ of degree $M(N)=M(8)=17$. The polynomial can still be displayed in the single printed page, $R(y)=$

$=153712881941946532798614648361265167-453762279414621179815552897029039797\,y+$

$+235326754101824439936800228806905073\,{y}^{2}-68875673245487669398850290405642067\,{y}^{3}+$

$+8129925258122948689157916436170874\,{y}^{4}-145759836636885012145070948315366\,{y}^{5}+$

$+2361976444746440513605248930610\,{y}^{6}+40525434802944282153115803370\,{y}^{7}+$

$+676326278232758784369966787\,{y}^{8}+62429137451114251409236415\,{y}^{9}+$

$+720991093724510065469933\,{y}^{10}+14670346929744822064505\,{y}^{11}+$

$+167556261648918275684\,{y}^{12}+917318495163561932\,{y}^{13}+$

$+3133529909492864\,{y}^{14}+4574211144896\,{y}^{15}-$

$-5932158016\,{y}^{16}+314432\,{y}^{17}$.

One should add that the search for the roots of this auxiliary “resultant” polynomial is a formidable numerical task. Using again the restriction to the mere real roots we obtained the following menu of seven eligible values

[TABLE]

[TABLE]

The insertion of any one of them would generate an independent EP8-supporting Hamiltonian $H^{(8)}$. These insertions as well as the derivation of their consequences remain routine.

3.2.2 The case of odd $N=2J=9$

At $N=9$ we have to deal with the quadruplet of polynomials

[TABLE]

[TABLE]

[TABLE]

[TABLE]

[TABLE]

[TABLE]

[TABLE]

The paradox of the simplification noticed after we moved from $N=2J$ to $N=2J+1$ at $J=3$ does not recur at $J=4$. The “resultant” polynomial $R(y)$ (again, of the 17th degree in $y$) does not factorize in any obvious manner. The whole process as sampled at $N=8$ must be repeated without any specific alterations. We omit the details here.

4 The EPN confluence of *all * of the

eigenvectors

4.1 The problem

Up to now we only studied the necessary conditions of the complete non-Hermitian degeneracy of the energy spectrum, $E_{n}\to E^{(EPN)}=0$ at all $n=0,1,\ldots,N-1$. Such a result must be complemented by the *elimination * of the *pathological * possibility of the replacement of Eq. (3) by its non-maximal-degeneracy alternatives like

[TABLE]

(with $N=4$ and $E=E^{\prime}=0$) or, in general, like

[TABLE]

(with $N=N_{1}+N_{2}$ and $N_{1}\geq 1$ and $N_{2}\geq 1$), or like

[TABLE]

(with $N=N_{1}+N_{2}+N_{3}$ and $N_{1}\geq 1$, $N_{2}\geq 1$ and $N_{3}\geq 1$), etc.

In Ref. [7] the disproof of the existence of pathologies at $z=z^{(EPN)}=1$ remained feasible, due to the absence of raounding errors at $N=4$ and $N=5$, via the explicit construction of the respective transition matrices $Q^{(N)}$. These matrices proved obtainable directly from the definition

[TABLE]

of the canonical-representation mapping (3). Such a recipe necessarily fails in the presence of any rounding error, i.e., as we saw, at any $N>5$. Even in the strictly EPN-compatible scenario, every numerically evaluated matrix $H^{(N)}(z)$ with $z\neq z^{(EPN)}$ then remains, in an arbitrarily small vicinity of $z^{(EPN)}$, diagonalizable. In the language of functional analysis one can say that the transition matrices $Q^{(N)}$ cease to exist at almost all parameters. Any numerical attempt of their construction must fail. The tests of the non-existence of the pathologies must be indirect.

4.2 The solution

Due to the ubiquitous presence of the numerical uncertainties we cannot construct the transition matrices $Q^{(N)}$ which mediate the isospectrality between our EPN-admitting matrix $H^{(N)}$ and its Jordan-block partner $J^{(N)}$. Formally this means that we will not be able to solve Eq. (22). In the EPN regime of interest, *any * standard numerical solver of such a problem would become wildly unstable. As an immediate consequence we lose the possibility of separating the EPN-related Hamiltonians from their isospectral alternatives with the Jordan blocks in subspaces, i.e., with a non-EPN, *incomplete * confluence of the eigenvectors.

In place of the exact Hamiltonians $H^{(N)}$ we only have access to their numerical representations $H=H^{(N)}+V$ containing a random, precision-dependent round-off perturbation $V={\cal O}(10^{-p})$. Up to the set of measure zero this makes our perturbed EPN-compatible Hamiltonians $H=H^{(N)}+V$, paradoxically, diagonalizable. Naturally, all of their normalized eigenvectors $|\psi_{n}\rangle$ are mutually almost parallel. This is also the property which enables us to test and confirm the occurrence of the EPN degeneracy even when $V\neq 0$.

We imagined that in many implementations of computer arithmetics the overall size of the round-off errors (i.e., of the exponent $p$ in the estimate of $V={\cal O}(10^{-p})$) can be varied. Whenever we amend the numerical precision $\sim{\cal O}(10^{-p})$ of the construction of our EPN-supporting Hamiltonian $H^{(N)}(z)$, we become able to distinguish between the correct, exhaustive, maximal degeneracy scenario (cf. Eq. (3)) and all of its incorrect, pathological alternatives characterized by one of the relations (19) or (20) or (21), etc.

A sample of such a test is presented in Table 1. In the test the numerical matrix of the form $H=H^{(6)}_{a}+V$ (cf. paragraph 2.3.1) with an unspecified, random round-off term $V={\cal O}(10^{-p})$ was assigned the numerically evaluated normalized eigenvectors $|\psi_{n}\rangle$. Their mutual confluence was then measured via an evaluation of the “non-overlaps” $\varrho_{mn}=1-|\langle\psi_{m}|\psi_{n}\rangle|$. In the given example they are clearly decreasing with the growth of $p$. The “no pathology” hypothesis (3) may be declared confirmed. The possibility of an *incompleteness * of the degeneracy of the $N-$plet of the eigenvectors of $H^{(N)}(z^{(EPN)})$ is persuasively excluded.

One should add that one of the important methodical merits of such a $p-$variation strategy is that it is rather robust. As long as it remains feasible at arbitrary matrix dimensions $N$, it enables us to confirm the absence of the pathological alternative scenarios in which the models with some hidden symmetry happen to be block-diagonalizable.

5 The spectra of energies $E=E(z)$ with $z\neq z^{(EPN)}$

5.1 $N=6$

In both of the one-parametric generalizations (7) of the conventional one-parametric $N=6$ Bose-Hubbard model one encounters, due to the absence of several conventional Lie-algebraic symmetries, a richer structure of the $z-$dependence of the spectrum. Both of the present deformed-model $N=6$ energy spectra can be perceived as the phenomenologically welcome immediate complements of their $N<6$ predecessors of Ref. [7].

5.1.1 The $a-$subscripted Hamiltonian

For the spectrum of the $a-$subscripted Hamiltonian one observes that the whole real line of $z$ is split into several subintervals of qualitatively different spectral form (cf. Fig. 1). In the leftmost subinterval, viz., for $z\in(-\infty,-z_{1})$ with $-z_{1}\approx-400$, the values of $\gamma=\sqrt{z}$ are purely imaginary. Hence, it is not too surprising that all of the six bound state energies remain strictly real.

At the emergent EP2 boundary $-z_{1}$ the innermost quadruplet of the energies merges and complexifies, pairwise, at $s\approx 2000$. Inside the subsequent interval $(-z_{1},1)$, therefore, the innermost quadruplet of the energies remains complex. We also find that the remaining auxiliary root $s_{3}=E^{2}_{\pm 3}$ is positive.

The latter value becomes negative (and keeps decreasing) in the next subinterval $(1,z_{2})$ with $z_{2}\approx 12.4$ (cf. the magnified picture in Fig. 2). In this interval the related energies $E_{\pm 3}$ become purely imaginary.

At the EP2-boundary $z_{2}$ one encounters again, at a positive value of $s\approx 104$, an unfolding of the two real roots. In the adjacent interval $(z_{2},z_{3})$ with $z_{3}\approx 15.5$ the roots $s_{1}$ and $s_{2}$ stay real and positive. They also keep dominating the third, negative real root $s_{3}$. The related four energies are real.

At the end of the above interval the middle real root $s_{2}$ changes sign. Inside the subsequent interval $(z_{3},z_{4})$ with $z_{4}\approx 835$ the further two energies $E_{\pm 2}=\sqrt{s_{2}}$ acquire the purely imaginary values. The maximal real root $s_{1}$ becomes also negative beyond $z_{4}$. In the ultimate, rightmost interval of $z\in(z_{4},\infty)$ this converts the further two energies $E_{\pm 2}=\sqrt{s_{2}}$ into purely imaginary quantities as well.

5.1.2 The $b-$subscripted Hamiltonian

Let us now turn attention to the other, $b-$subscripted Hamiltonian and to its spectrum (cf. Fig. 3). The observations and conclusions remain similar. First of all, the leftmost EP2 boundary gets merely shifted to $-z_{1}\approx=-66.527$ yielding the energy merger and subsequent complexification at $s\approx 346.2$.

In the second and third intervals $(z_{1},z_{2})$ and $(z_{2},1)$ with boundary $z_{2}\approx 0.01041$ we notice the change of the shape of the curve $s(z)$ near $z^{(EP6)}=1$. This is displayed in Fig. 4. The picture also implies the obvious classification of the energies. Four of them keep complex for $z\in(-z_{1},z_{3})$ with $z_{3}\approx 131.6044$ while the remaining two remain real. In the subinterval $z\in(z_{2},1)$ they get again converted into the purely imaginary quantities. Finally, in the rightmost interval of $z\in(z_{3},\infty)$ the model supports the two real and four purely imaginary energy eigenvalues.

We may conclude that besides the preselected EP4 separator of the real line of $z$, both of our “a” and “b” $N=6$ Hamiltonians also appeared to keep the parallels with the $N=4$ results of Ref. [7]. Thus, all of these generalized Bose-Hubbard-type models are found to share the emergence of a few remote EP2 separators as well as of the qualitatively different types of spectra supported by the individual, model-dependent subintervals of $z$.

5.2 $N=7$

In the vicinity of the most interesting EP parameter $z^{(EP7)}=1$ the local shape of the single real $\alpha-$subscripted curve $s_{\alpha}(z)$ is shown in Fig. 5. In comparison with its two $N=6$ predecessors (cf. Figs. 2 and 4 above), this shape is still different. The horizontal line which played just the auxiliary, eye-guiding role at $N=6$ must be now reinterpreted as representing also one of the energies.

In a more detailed analysis of the spectrum the points of the change of the sign of $s_{\alpha}(z)$ can be also localized exactly. In Fig. 5 we have $s_{\alpha}(z)=0$ not only at the EPN prescribed value $z=1$ but also at the small negative exact nodal point

[TABLE]

In the global picture of Fig. 6 we also find the third real zero of $s_{\alpha}(z)$ at $z_{-}\approx-44.39497$.

Naturally, after the change $\alpha\to\beta$ of sign in Eq. (17) the zeros of the real function $s_{\beta}(z)$ get shifted,

[TABLE]

A surprise still emerges at the second auxiliary root $r_{\beta}$. Qualitatively, the global shape of $s_{\beta}(z)$ parallels Fig. 6 – only a rescaling of the axes makes the main difference. In contrast, locally, a magnified detail as given in Fig. 7 shows a different, utterly unexpected new pattern near $z^{(EP7)}=1$ when compared with Fig. 5.

6 Discussion

6.1 Complex symmetric Hamiltonians

In the Heisenberg’s more than ninety years old formulation of quantum mechanics [10] an important role has been played by the finite-dimensional and real symmetric matrices $H=H^{T}$. Nowadays, quantum physics still relies heavily upon the inspiration provided by such an elementary mathematics. In the Kato’s influential monograph [3], for example, several abstract features of various advanced quantum Hamiltonians may still be found illustrated by the most elementary two by two matrices, real or complex, and symmetric, or not. Typically, multiple relevant physical questions (e.g., of the stability of quantum systems under small perturbations) as well as some related mathematical concepts may be found explained there in terms of these utterly elementary toy models.

The recent return of interest to less conventional elementary models (cf., e.g., the use of *complex symmetric * matrices in Refs. [11, 12, 13]) found its motivation in an innovative physical interpretation of EPs. Specifically, we may mention the relevance of EPs in the study of classical as well as quantum dynamical systems [14, 15, 16] and/or the key role of the unfolding of the EPs of the $N-$th order in the description of the dynamics of the Bose-Einstein condensation [1].

What is observed in the new research motivated by similar ideas is not only an enormous increase in the experimental activities (cf., e.g., their concise list in [17]) but also the growth of quality of theoretical investigations [18]. Typically, the interest of mathematicians in the properties of finite-dimensional complex symmetric matrices [19] is now being extended to the infinite-dimensional cases [20] and to the various non-Hermitian operators exhibiting the innovative Krein-space self-adjointness *alias * ${\cal PT}-$symmetry [2, 5, 18]. In parallel, the accessibility of the efficient, computer-assisted symbolic manipulation techniques opens also new horizons in our understanding of certain less conventional finite-dimensional models (cf., e.g., the lists of references in [1, 21, 22]).

6.2 Bose-Hubbard and non-Bose-Hubbard models

The latter studies inspired our theoretical investigation [7] where we considered the tridiagonal, ${\cal PT}-$symmetric and complex-symmetric $N$ by $N$ matrix Hamiltonians $H^{(N)}(A,B,\ldots,Z)$ up to $N=5$. The availability of the utterly elementary formulae for the multiparametric energies $E_{n}^{(N)}(A,B,\ldots,Z)$ enabled us to succeed in obtaining an exhaustive $N\leq 5$ classification of the subset of models $H^{(N)}(A,B,\ldots,Z)$ guaranteeing, at an ad hoc set of parameters $A^{(EPN)},B^{(EPN)},\ldots,Z^{(EPN)}$, the existence and guarantee of the confluence of *all * of the lower-order EPs into a single, maximal EPN singularity.

In our present paper we succeeded in extending these results to all $N$. The level of the complexity of our innovated, unusual Hamiltonians and of their spectra appears controlled, roughly speaking, by the integer part $J$ of $N/2$. A remark should be added concerning the presentation of the results at the larger dimensions $N$ when all formulae become rather long, not fitting a printed-page format. We, nevertheless, imagined that a detailed discussion of the models with $N\leq 9$ already provides a sufficiently comprehensible picture of the situation and, in particular, of the universality of our construction of the generalized Bose-Hubbard-type Hamiltonians, exhibiting still the presence of the $N-$th order exceptional-point degeneracy.

Among the current applications of the similar models in physics one can notice an intensification of interest in the so called open quantum systems in which one admits a non-unitarity of the evolution caused by an uncontrolled interaction with an environment [23]. It is not too surprising that one of the most natural theoretical formulations of such a situation is offered by a combination of the EP-related phenomenology with the real- or complex-symmetric-matrix mathematics. *Pars pro toto * we should mention here a well-balanced combination of the mathematical and phenomenological insight provided by the Bose-Hubbard multi-bosonic models studied in both of their real- and complex-symmetric versions – cf., e.g., [24] and [25], respectively.

6.3 Non-equivalent concepts of non-Hermitian physics

In conventional quantum mechanics the predictions concerning an observable $\Lambda$ are probabilistic, expressed in terms of matrix elements $\langle\psi(t)|\Lambda|\psi(t)\rangle$. The time-dependence of the states $|\psi(t)\rangle$ is controlled by Schrödinger equation

[TABLE]

Traditionally, the evolution is assumed unitary so that the Hamiltonian itself must be chosen selfadjoint, $H=H^{\dagger}$ [26]. A less conventional, non-unitary-evolution paradigm in which $H\neq H^{\dagger}$ is being developed in several alternative directions at present [18]. Opening several entirely new areas of research. Thus, whenever one speaks about non-Hermitian Hamiltonians (i.e., about operators $H$ such that $H\neq H^{\dagger}$ in some preselected Hilbert space ${\cal L}$), it is necessary to check their spectrum $\sigma(H)$. This enables us to distinguish between the so called quasi-Hermiticity (best known and used in nuclear physics [22], with $\sigma(H)\in\mathbb{R}$) and the less unusual genuine non-Hermiticity taking place whenever $\sigma(H)$ is *not * real [23].

In our present paper our attention was concentrated on the open-system scenario in which the Hilbert space ${\cal L}$ of ket vectors $|\psi(t)\rangle$ is considered physical. For this reason the non-unitarity of the evolution (caused by the non-Hermiticity of $H$) is accepted as natural, finding its explanation in a non-negligible interaction of the quantum system in question with a certain rather vaguely specified environment (cf., e.g., the original ideas of Feshbach [27] and Löwdin [28] as well as their multiple most recent theoretical developments and implementations as sampled, say, in reviews [23, 29]).

6.4 Exceptional points: Theory vs. experiment

The Kato’s [3] concept of exceptional points found multiple applications in the study of Hamiltonians which are non-self-adjoint in a preselected Hilbert space ${\cal L}$, physical or not. Typically, the pioneering simulation of the confluence of two quantum resonances in [11] offered one of the first *experimental * localizations of the EP2 anomaly in a classical microwave setup. What followed was a long series of other experimental studies of EP2s covering, say, their manifestations in electronic circuits [30], in exciton-polariton billiards [31] or in optomechanical systems [32]. Not too surprisingly, a “natural” transition to the analogous experimental EPN scenarios with $N\geq 3$ opened a number of new and challenging questions and obstacles [33]. In fact, one of the first physical systems making the EP3-related patterns experimentally accessible was only proposed, in the pioneering study [12], in 2012.

Also the progress in the underlying theory of the localization of EPNs is not too quick. Its characteristic nontrivial $N=3$ implementation as proposed in [34] consisted, for example, of a mere triplet of coupled wave-guides, with the gain and loss areas arranged in a ${\cal PT}-$symmetric manner. An experimentally feasible arrangement of the semiconductor wave-guides has been determined. The most elementary tridiagonal and complex symmetric choice of the matrix $H^{(3)}$ assigned to the system appeared realistic and sufficient for the purpose.

The idea of the study of the role of EPs in Bose-Hubbard-type models was followed also in [7]. The present further development and completion of such a project was still motivated by the generic phenomenological appeal of the models with higher-order EPs as well as by the feasibility of an experimental realization and verification of the EP-related theoretical hypotheses and measurable predictions. A short outline of the history was provided here in Introduction. In our present paper we only re-emphasized the relevance of the higher-order EPs and, in particular, of the guarantee of their occurrence in the fairly broad class of the generalized, Bose-Hubbard-inspired $N$ by $N$ matrix Hamiltonians $H^{(N)}$.

7 Summary

Our present extension of the non-BH constructions beyond $N=5$ was achieved only after a thorough innovation of the method. The main obstacle appeared represented by the presence of the round-off errors in the numerically constructed Hamiltonian matrices $H_{(non-BH)}^{(N)}(z)$ at the larger dimensions. Originally we were sceptical in this respect. A decisive progress has only been achieved when we turned attention to the computer-assisted constructions using the variable-precision arithmetics. Without an explicit construction of transition matrices $Q^{(N)}$, this enabled us to check the existence and maximality of the EPN anomaly.

In its present form, our arbitrary$-N$ construction of the non-BH Hamiltonians still originates from the energy-determining secular equation. The energy-degeneracy constraint (10) then yields the polynomial constraints equivalent to a nonlinear algebraic set of $J$ coupled polynomial equations for $J$ unknown quantities $A^{(EPN)},B^{(EPN)}$, $\ldots$, $Z^{(EPN)}$, real or complex [cf. Eq. (12)]. Such a set is, in general, solvable by the well known Gröbner-basis elimination technique. Thus, at any Hamiltonian-matrix dimension $N$ the Gröbner-basis elimination reduces the construction of the parameters $A^{(EPN)},\ldots$ to the localization of the roots $y_{k}$ of a single polynomial $R^{(N)}(y)$.

The viability of such an approach and algorithm has been confirmed here to work very comfortably up to the fairly large matrix dimensions $N$. A decisive technicality has been found to lie in the localization of the “seed” roots $y_{k}$. The reason is that the degree $M(N)$ of polynomial $R^{(N)}(y)$ grows quickly with $N$. For example, we had $M(N)=17$ in our $J=4$ samples of dimensions $N=8$ and $N=9$. Thus, the main criterion of the practical feasibility of the construction is given by the value of $M(N)$.

Incidentally, the quick increase of $M(N)$ makes it practically impossible to display any sample results in print at $N\geq 10$. Even the very length of the individual (integer) numerical coefficients in $R^{(N)}(y)$ would exceed the single-line capacity in such a case. Fortunately, whenever needed, it makes sense to keep the $N\geq 10$ results stored in the computer while asking just for the display of output containing the measurable (e.g., spectral) predictions.

In the physical, phenomenological model-building context the most characteristic appeal of our present EPN-supporting toy-model Hamiltonians is twofold. Firstly, these models offer a truly rich menu of the parameter-dependence of the energy spectra even far from the EPN merger. Secondly, the price to be paid for the flexibility remains acceptable. Thus, the key message delivered by our study is that one can deal with the non-BH Hamiltonian matrices at the unexpectedly large dimensions $N$. The way to the success has been found in keeping the influence of the round-off errors under control. This was shown to confirm that in the EPN limit of $z\to 1$ the non-BH quantum systems really *do * reach the *maximal * non-Hermitian degeneracy at which *all * of the eigenvectors of the Hamiltonian *do * simultaneously merge.

Acknowledgement

Work supported by GAČR Grant Nr. 16-22945S.