ABC of ladder operators for rationally extended quantum harmonic oscillator systems
José F. Cariñena Mikhail S. Plyushchay



 























30003000
18pt0pt
0pt0pt
0pt0pt
32pt35by 460pt







 addtoresetequationsection.equation










=-2.0truecm 3000 3000 18pt 0pt 0pt 60pt 0pt 
0pt 32pt 36 by 470pt
addtoresetequationsection.equation

    𝒞𝒲𝔔ℋ𝒫𝒦𝒟ℒRQSN1𝔸𝔹ℤℝ𝕎ℙ𝕋ℂℍ̋0 grx⃗y⃗P⃗p⃗∇⃗c⃗A⃗P⃗b⃗ #1 #1
            [ 𝔰𝔬σπsncndnnsncndcdsdsechℋ𝒜𝒫ℒ𝒮𝒬𝒟𝒱𝒲#1red#1#1blue#1#1green#1
             
            ABC of ladder operators for
            rationally extended quantum harmonic oscillator
            systems 
             
                 José F. Cariñena^a and Mikhail S. Plyushchay^b
            
            [8pt]
            ^aDepartamento de Física Teórica, 
            Universidad de Zaragoza, 50009 Zaragoza, Spain
            ^bDepartamento de Física,
            Universidad de Santiago de Chile, Casilla 307, Santiago 2,
            Chile  
            
            [4pt]
             E-mails:  jfc@unizar.es, 
            mikhail.plyushchay@usach.cl
                
            =======================================================================================================================================================================================================================================================================================================
            
            
            The problem of construction of ladder operators for
             rationally extended quantum harmonic oscillator (REQHO)
             systems of a general form is investigated in the light of
            existence  of different  schemes of the
            Darboux-Crum-Krein-Adler transformations
            by which such systems
            can be generated from the quantum harmonic oscillator. 
            Any 
            REQHO system is characterized by  
            the number of separated states in its spectrum,
            the number of `valence bands' in which the 
            separated states are organized, and by the total
            number of the missing energy levels
            and their position.
            All these peculiarities  
            of a REQHO system
            are shown to be detected and reflected by 
            a  trinity (𝒜^±,
            ℬ^±, 𝒞^±)  
            of the basic  (primary) lowering and raising 
            ladder operators related 
            between themselves by certain algebraic identities with 
            coefficients   polynomially-dependent 
            on the Hamiltonian.
            We show that all the  
            secondary, higher-order ladder operators
            are obtainable  by a composition  of  
             the basic ladder operators of the trinity
            which form  the set of the 
            spectrum-generating operators.
            Each trinity, in turn,  can be constructed from
            the intertwining operators of the two 
            complementary minimal schemes of the 
            Darboux-Crum-Krein-Adler transformations.
            .5cm
            
            § INTRODUCTION
            
            
            There are two basic  exactly solvable 
            quantum mechanical systems 
            which reveal themselves directly or indirectly
            in association with other systems
            and play a fundamental role in many
            physical theories and applications.
            One of them is a free particle characterized
            by a continuous spectrum.   The other 
            one
            is a harmonic oscillator with its 
            infinite equidistant discrete
            spectrum of the bound states.     
            
            Free particle is essential, in particular,
            for  understanding  the properties of the soliton solutions 
            in  integrable systems.  Quantum reflectionless
            potentials represent a snapshot of soliton solutions
            to the classical KdV equation.
            Initially, reflectionless potentials  
            were obtained  with the help of 
            the method of the inverse scattering transform in
            solving the problem of theoretical
            construction of a solid dielectric 
            medium that is perfectly transparent to 
            electromagnetic radiation <cit.>.
            This  important class of  the systems can also be 
            generated from the quantum  free particle  
            by means   of Darboux  transformations (DTs) and 
            their generalization in the 
            form   of  Darboux-Crum transformations 
            (DCTs) 
            <cit.>.
            
            Notice here that reflectionless 
            systems 
            appear for instance in 
            the Gross-Neveu model <cit.>
            in the context of the hadron physics <cit.> and the 
            physics of conducting polymers <cit.>.
            Any quantum reflectionless system with n bound states in its spectrum
            is characterized by 
            the presence of a nontrivial Lax-Novikov integral that is a differential operator 
            of order 2n+1 having  a structure of a Darboux-dressed
            momentum operator of the free particle system.
            It is this integral of motion  that distinguishes the  states in  
            the doubly degenerate continuous part of the spectrum of a reflectionless system
            and detects all the non-degenerate bound states 
            as well as a state at the very 
            edge of the continuous part of the spectrum by annihilating them
            <cit.>.
            It is the same  operator that plays a fundamental role 
            in the theory of nonlinear integrable systems 
            <cit.>.
             The peculiarity  of reflectionless systems 
            also reveals itself in the nature of
            the quantum mechanical supersymmetry associated with them.
            Instead of a usual linear or non-linear 𝒩=2 
            supersymmetric structure which appears in 
            an extended system composed from a  pair
            of  quantum mechanical systems related by 
            a DT or DCT, the
            extended  system composed from a pair of reflectionless
            systems is described by the 
            exotic nonlinear 𝒩=4 supersymmetry
            generated by two pairs of supercharges 
             alongside with the two bosonic integrals 
            of motion constructed from  the Lax-Novikov integrals
            of the subsystems <cit.>. The appearance of the 
            exotic nonlinear supersymmetric
            structure associated with reflectionless systems 
            is traced  to the fact that any reflectionless Schrödinger system 
            can be related to a free particle system  not 
            only by one  but by two different DCTs due to 
            the presence of the momentum operator 
            in the quantum free particle system in the form of
            its integral of motion [A general  picture is more
            complicated, however.  In the case of a coincidence of some or all  
            discrete energy levels  of the two subsystems,  the supersymmetric 
            reflectionless partners can be related directly by a DCT  of a lower order,
            without a necessity to construct a chain of Darboux transformations via 
            a  free particle system. This happens  due to 
            the opening of a kind of a  `direct tunnelling channel' 
            that can appropriately be understood from the standpoint of 
            the picture of  soliton scattering, see <cit.>.]. 
                
              
            Some time ago, there  has been discovered a new very broad
            class of exactly solvable systems 
            which represent certain rational extensions 
            of the quantum harmonic oscillator 
            <cit.>.
            The eigenstates of rationally extended quantum harmonic oscillator 
            (REQHO) systems are given in terms
            of exceptional Hermite polynomials,
            and can be obtained from  the quantum harmonic oscillator (QHO) system 
            by an appropriate 
            DCT, or its further generalization in the form 
            of the Darboux-Crum-Krein-Adler transformation
             <cit.> [Exceptional Jacobi 
            and Laguerre polynomials 
            <cit.>
            can be associated 
            with Darboux-Crum transformed free particle 
            on finite interval (particle in infinite potential well) 
            and isotonic oscillator, respectively; the isotonic oscillator, in turn,
             can  be related 
            to the QHO by a singular Darboux transformation.].
            In what follows we shall refer to any 
            generalized 
            Darboux transformation 
            of the QHO with
            intertwining operators to be higher order differential operators
            as the Darboux-Crum-Krein-Adler transformation (DCKAT).
            Instead of a continuous spectrum of the free particle, 
            the QHO is characterized 
            by an infinite discrete spectrum of  bound states.
            In spite of such a radical difference, 
            the QHO  is also a very peculiar system
            because  its discrete spectrum
            is equidistant. As a consequence, instead of the 
            Hermitian momentum operator integral
             that encodes reflectionless nature of the free particle, 
            the QHO possesses a pair of
            Hermitian conjugate ladder operators
            which are spectrum-generating operators.  
            Similarly to the relation between reflectionless 
            and free particle systems,
            a given REQHO system can be obtained from  
            the QHO by different DCKATs. 
            Then one can expect  by analogy 
            with the pairs of reflectionless systems 
            related by DCTs and the exotic supersymmetric structure 
            associated with them that the 
            REQHO systems should be characterized 
            by some special properties.
            Particularly, it seems to be natural to  expect 
            the appearance of peculiarities 
            related to the ladder operators
            for such a family of quantum
            systems. 
            
            
            It is worth to mention  here that 
            finite-gap systems, a  limit case 
            of which corresponds to reflectionless systems,
            and the QHO
            are generated by the periodic 
            Darboux chains  <cit.>. The
            last construction also produces 
            Painlevé equations <cit.>, 
            that are  intimately related with isomonodromic 
            deformations of linear systems and integrability properties
            of nonlinear systems in partial derivatives
            <cit.>.
            The REQHO systems are  
            isomonodromic 
            deformations of the QHO <cit.>. 
            
            Some investigations on ladder operators 
            in REQHO systems have already been realized
            in <cit.>. 
            This has been done, however, for some particular examples
            of the REQHO systems or for some particular families 
            of such systems, while the problem of construction 
             of  ladder  operators and investigation of their properties
            for REQHO systems of a general
            form  remains open. 
            Note also that  in the indicated 
            works only some special aspects related to the ladder 
            operators of the REQHO systems were studied.
             In particular, in a recent paper 
             <cit.>  we have considered 
            the problem of  construction of the ladder operators 
            for the simplest case of the REQHO
             by exploiting the simplest DT 
            which relates the system to the QHO. 
            We also investigated there the discrete chains  
            related with the obtained ladder operators
            given by a pair of Hermitian conjugate  
            third order differential operators
               to be fermionic generators of the 
            polynomially deformed
            bosonized 𝔬𝔰𝔭(1|2) superalgebra.
            
            In the present article we investigate the
            problem of construction  of ladder
             operators 
            for  REQHO systems 
            of a general form in the light of
            existence of different DCKATs by which any such a system
            can be related to and generated from the QHO. 
            In this point there shows up 
            a similarity of the REQHO systems with 
            reflectionless systems related with a free particle.
             We  
            show that for any REQHO, there exists a trinity
            of the basic  (primary) pairs of the   lowering and raising 
            ladder operators. This trinity is proved to form the set of the 
            spectrum-generating operators which  
            detects and reflects all the peculiar properties of a
            given REQHO  system.  
            
            The paper is organized as follows.
            In the next Section we briefly review
            the properties of the DTs  
            and their generalization in the form of the 
            DCKATs. 
            In Section <ref> we discuss 
            general schemes of the DCKATs for generation 
            of the REQHO systems from the QHO
            by using different
            sets of physical and non-physical eigenstates of  the 
            latter.
            Any REQHO system can be characterized by the 
            total number of separated states in the low part of its 
            spectrum, by the number of  the 
            gapless `valence bands'
            in which the eigenvalues of separated states  are organized,
            and by the total number of missing energy levels 
            and their position in the low part of the spectrum.  
            In Section <ref> we consider
            in detail the problem of construction of 
            different  ladder operators
            for the simplest REQHO system with one separated state
            and two missing energy levels in a unique gap 
            that separates it  from the equidistant infinite part of the spectrum.
            We investigate there  the general properties
             and relations between the basic
             (primary) and  secondary, higher-order  ladder operators.  
            In Section <ref> we consider two more
            particular examples of the REQHO systems,
            one of which corresponds to a generalization of the system
            from Section <ref>. Another example
            corresponds to a REQHO system with one 
            valence band composed from the two
             energy levels separated from the equidistant 
             infinite part of the spectrum by  the gap of two missing 
             energy levels. The results and observations obtained
             in Sections  <ref> and <ref>
             are generalized then in Section <ref>
             for the case of the REQHO systems of a general form.
             Section <ref> is devoted to a summary of the
             obtained results,
             where we also indicate some problems 
             that could be interesting for 
            further investigation.
            
            
            
            
            
            § DARBOUX-CRUM-KREIN-ADLER TRANSFORMATIONS
            
            
            Let ψ_*(x) be a nodeless physical or non-physical
            real eigenfunction of a Hamiltonian operator H=-d^2/dx^2+V(x) 
            with eigenvalue E_*, Hψ_*=E_*ψ_*.
            We assume that potential V(x) is a non-singular real function on all the real line
            .
            Define the first order differential operators
            
            A≡ψ_*d/dx1/ψ_*=d/dx-𝒲 ,
                
            A^†=-1/ψ_*d/dxψ_* ,
            
            where 
            𝒲=ψ'_*/ψ_*,  
            ψ'_*=dψ_*/dx.
            They   factorize the shifted Hamiltonian,
            H-E_*=A^† A,
            whose potential  is given by
             V=𝒲^2+𝒲'+E_* in terms of the superpotential
              𝒲 and factorization energy E_*.
            The product with permuted operators  A and 
            A^† defines a supersymmetric partner
            Hamiltonian, AA^†≡H̆-E_*,
            H̆=-d^2/dx^2+V_*, 
            for which V̆=𝒲^2-𝒲'+E_*.
            The relation of the potential V̆(x) of the superpartner H̆ 
            to the potential V(x) of  the system H 
            can be rewritten in a more convenient form
             for a further generalization,
            V̆=V-2(lnψ_*)” .
            From factorization relations 
            it follows that  the operators A and 
            A^† intertwine the partner Hamiltonians,
            
            AH=H̆A ,     
            A^†H̆=H A^† . 
            
            As a consequence, A and 
            A^†  mutually map
            the eigenstates of the superpartners.
            Namely, if ψ(x;E) is a physical or non-physical
            eigenstate of H of eigenvalue E≠ E_*,
            Hψ(x;E)=Eψ(x;E),
            then Ψ(x;E)=Aψ(x;E) is an eigenstate
            of H̆ of the same eigenvalue and of the same physical or 
            non-physical nature.
            Vice versa, if Ψ(x;E) is an eigenstate
            of H̆ of eigenvalue E≠ E_*,
            then A^†Ψ(x;E) is an eigenstate
            of H of the same eigenvalue
            and of the same nature.
            
            If ψ(x) is a solution of the second order 
            differential equation Hψ(x)=Eψ(x) 
            for arbitrary value of E,  a second, linearly independent 
            solution of this equation is 
            ψ(x)≡ψ(x) ∫^x dξ/(ψ(ξ))^2 .
            
            If ψ(x) is a normalizable on  function, then  
            ψ(x) is not normalizable, and vice versa. 
            The properties of A and A^† as the operators 
            that mutually map the corresponding eigenstates 
            ψ(x) of H and Ψ(x) of 
            H̆  of eigenvalue E≠ E_* are valid also for  
            the associated  eigenstates ψ(x)
            and Ψ(x). 
            The case E=E_* in this context is different.
            The eigenstate ψ_* of H of eigenvalue E=E_* 
            constitutes  the kernel of the  operator A,
            Aψ_*=0. The same is valid for the 
            state Ψ_*(x)≡ 1/ψ_*(x), which constitutes the 
            kernel of A^†, A^†Ψ_*=0,
            and is the eigenstate of H̆ of the same eigenvalue E=E_*.
            In this special case of E=E_*, 
            the operator  A  transforms the state ψ_*
            into the state Ψ_*=1/ψ_*∈ (A^†):
            
                Aψ_*(x)=(ψ_*d/dx1/ψ_*)ψ_*(x)=
                ψ_*(x)d/dx∫^x dξ/(ψ_*(ξ))^2=
                1/ψ_*(x) .
            
            Analogously, the state 
            Ψ_*=(1/ψ_*)
            is transformed by A^† into ψ_* being a kernel of A.  
            More details can be found in <cit.>
            where the same notation is used.
            
            Let us stress that for our constructions below it is important 
            that 
            if ψ_* is a normalizable state,
            then the  state ψ_* is not normalizable, 
            and vice versa, and that the same
            property  is valid  for the pair of
            the states Ψ_*=1/ψ_* and Ψ_*. 
            Before  the discussion of the concrete quantum systems
            we, however, completely  neglect  the questions of 
            normalizability of the corresponding wave functions.
            Furthermore, here and in what follows we do not 
            preoccupy about normalization of the states, 
            and specify wave functions   
            modulo a constant multiplication factor.
            
            
            The DT construction can be generalized for the case of
            the DCKAT. The latter is  
            generated on the basis of several seed eigenstates
            ψ_i_1, ψ_i_2,…,
            ψ_i_n of H
            of  different eigenvalues 
            E_i_k, Hψ_i_k=E_i_kψ_i_k,
            k=1,…, n, with E_i_k≠ E_i_k' for i_k≠ i_k'.
            To get a nonsingular partner system H_n,
            these states  should be  such  that 
            their Wronskian _n(x)≡(ψ_i_1(x), ψ_i_2(x),…,
            ψ_i_n(x))=||ℱ(x)||,
            ℱ_ij=d^i-1/dx^i-1ψ_j,
            i=1,…, n, j=i_1,…, i_n,  is a nodeless function.
             The partner system H_n is then given  by the potential
             
             V_n(x)=V(x)-2(ln_n(x))”.
             
             The eigenstates ψ(x;E) of H are mapped into 
             the eigenstates Ψ(x;E) of H_n of the same eigenvalue E 
             via the relation 
             Ψ(x;E)=(ψ_i_1(x),…,ψ_i_n(x),ψ(x;E))/_n(x) .
             
            In the case of  n=1, Eq. (<ref>) reduces to 
            Eq. (<ref>) while (<ref>) 
            reduces to  the relation 
            corresponding to the case of the DT  presented in the form
             Ψ(x;E)=A_1ψ(x;E),
            with ψ_i_1=ψ_*, A_1=A and H_1=H̆.
             Furthermore, 
            for n>1
            Eq. (<ref>) can 
            be presented in a form that generalizes 
            the indicated DT's formula. 
            For this  we iteratively define
             two sequences of related 
            differential operators A_m, m=1,…,
            and 
            _m, m=0,1,…,
             as follows: _0≡ 1,
            
            A_m=(_m-1ψ_i_m)d/dx1/(_m-1ψ_i_m) ,
                
            m=1,2,…,
            _n=A_n A_n-1… A_1 .
             
            Let us denote  H_0≡ H
            and define  H_m-1=A_m^† A_m+E_i_m.
            We have  
             A_m A_m^†=H_m-E_i_m,
             and  obtain a generalization of 
            the intertwining relations (<ref>),
            A_m H_m-1=H_m A_m,
            A_m^† H_m=H_m-1 A_m^†.
            Then  we can present 
              (<ref>) 
            in the equivalent form
            Ψ(x;E)=_nψ(x;E) . 
            
            
            The n-th order differential 
            operators _n and _n^† 
            intertwine the partner systems H_0 and H_n,
            _n H_0=H_n_n ,    _n^† H_n=H_0_n^† .
            
            The products of  operators   _n and _n^† 
            turn out to be polynomials in the corresponding Hamiltonian operators
            with  roots equal  to the energies of factorization:
            _n^†_n=∏_k=1^n (H_0-E_i_k) ,    _n _n^† =∏_k=1^n (H_n-E_i_k) .
            
            
            
            For non-singular partners
            H_0 and H_n, some or all of the 
            intermediate second-order differential operators 
            H_m with m=1,…,n-1
            can be singular.
            If we change the order of the 
            seed states, the Wronskian is left invariant 
            modulo a possible multiplication by (-1),
            that does not change Eqs. (<ref>) and
            (<ref>). Differential operator _n 
             is not changed either, and  from (<ref>) 
             one can conclude   that  
            the kernel of _n is spanned by 
             the complete set of the seed states, 
             (_n)=span {ψ_i_1,…,ψ_i_n}. 
            However, under  permutation in the order
            of the seed eigenstates 
            the first order differential operators
            entering into the factorized form (<ref>) 
            of _n  are  changed. 
            As a result, the nature of some or of all
            of the intermediate
            systems H_m can be changed from 
            a non-singular (singular) to a singular (non-singular).
            
            
            
            
            
            § GENERATION OF  THE REQHO SYSTEMS FROM THE QHO
            
            
            Let us turn now to the  specific example of the QHO
            system 
            given by the  Hamiltonian operator 
            H_osc=-d^2/dx^2+x^2.
            Its bound  eigenstates are described 
            by (not normalized here)  wave functions 
             ψ_n(x)=H_n(x)e^-x^2/2
            which correspond to eigenvalues
            E_n=2n+1, n=0,1,…, where H_n(x) are 
            Hermite polynomials. 
            The change of variable x→ i x
             generates a change
            of the sign of the Hamiltonian,
            H_osc→ -H_osc,  and so,  transforms physical eigenstates ψ_n(x) 
             into non-physical eigenstates  
            ψ^-_n(x)=ℋ_n(x)e^x^2/2
             of H_osc of eigenvalues E_n^-=-(2n+1), n=0,1,…, where 
             ℋ_n(x)=H_n(i x).
            Unlike H_n(x),  the polynomials ℋ_n(x) with even index have 
            no real zeros
             while the unique real zero of these polynomials with odd index is 
             at x=0.  More details can be found in <cit.>.
             
            The choice  of the ground state
            ψ_0=e^-x^2/2 
             as a seed eigenfunction ψ_* for the DT in (<ref>)
            generates the first order differential operators
            
            A_1=d/dx+x≡ a^- ,    
            A^†_1=-d/dx+x≡ a^+ .
            
            They factorize the shifted Hamiltonian, 
            a^+a^-=H_osc-1, and satisfy the commutation relation
            
            [a^-,a^+]=2 .
            
            As a consequence  we have 
            
            [H_osc,a^±]=± 2 a^± .
            
            The two relations in (<ref>)  can be rewritten equivalently in the form
            
            a^-H_osc=(H_osc+2)a^- ,    
            a^+(H_osc+2)=H_osc a^+  ,
            
            and mean that a^- and a^+ are  the ladder lowering and 
            raising operators of the QHO.
            Relations in (<ref>) can also be interpreted
            as that a^- and a^+ intertwine 
            the QHO system H_osc=H_0 with the SUSY-partner system 
            H_1=H_osc+2
            which is just the shifted
            QHO.
            Since a^-ψ_0=0, 
            from the point of view of the DT
            one can consider the shifted system
            H_1 as the QHO with the removed ground 
            state. But since a^-ψ_n=ψ_n-1,
            the partner system H_1
            is  the same QHO 
            with the Hamiltonian shifted for +2. 
            In this DT picture, the ground state ψ_0 of H_osc
            is created by application of A_1^†=a^+
            to the function  
            ψ^-_0=(1/ψ_0)
            that  is  a non-physical eigenstate of H_1 of 
            eigenvalue E=1.
            
            The DCKAT generated on 
            the basis of the set of the seed eigenstates 
            ψ_0,…,ψ_n-1,ψ_n,
            n=1,…,
            produces the intertwining operators
            _n+1^-=(a^-)^n+1 and 
            ^+_n+1=_n+1^†=(a^+)^n+1,
            which intertwine the QHO with the 
            SUSY-partner system
            H_n+1=H_osc+2(n+1),
            _n+1^-H_osc=H_n+1_n+1^-,    _n+1^+ H_n+1=H_osc_n+1^+ .
            
            The case n=0 here reproduces  
            the relations of the DT generated by
            the choice  ψ_*=ψ_0.
            Note that in the case of n=1 
            the permutation of the seed states
            in the DCKAT construction,
            (ψ_0,ψ_1)→ (ψ_1,ψ_0),
            gives rise to a singular 
            operator 
            
            A_1=ψ_1d/dx1/ψ_1=
            d/dx+x-1/x≡ a^-_iso .
            
             This operator 
            acting on the second 
            chosen seed eigenstate
            ψ_0 gives a function -1/xe^-x^2/2=-1/ψ_1^-, 
            which according to 
            (<ref>) generates the second 
            factorization operator
            A_2=d/dx+x+1/x that also is singular.
            The product of these two singular operators
            gives the same second-order non-singular intertwining 
            operator [For the discussion 
            of related phenomena in a context of the quantum 
            second-order supersymmetry 
            anomaly and coupling-constant metamorphosis see ref. <cit.>]
             as  the scheme
            with the (ψ_0,ψ_1) pair:
            A_2A_1=(a^-)^2=_2^-. 
            The intermediate Hamiltonian 
            H_1 in this case is a singular at x=0 operator
            corresponding to the (shifted)  quantum isotonic oscillator,
            A_1^† A_1=H_osc-3,
            A_1 A_1^†=H_1-3,
            which is given by the  potential 
            V_1(x)=x^2+2/x^2+2.
            We also have here the relations 
            A_2^† A_2=H_1-1,
            A_2A^†_2=H_osc+3.
            
            
            The choice of  a nodeless non-physical state ψ_*=ψ^-_0 
            corresponding to factorization energy E_*=-1
            in the DT construction gives 
            
            A_1=d/dx-x=-a^+ ,    
            A_1^†=-a^- .
            
            We obtain the same ladder operators of the QHO,
            but now they will intertwine the QHO with the
            SUSY-partner system H_1=H_osc-2, 
            
            a^+H_osc=(H_osc-2)a^+ ,    
            a^-(H_osc-2)=H_osc a^-  .
            
            Here the action of A_1 on the non-physical eigenstate
            ψ^-_0=(1/ψ_0)
             of H_0=H_osc of 
            eigenvalue E=-1 produces a physical ground state ψ_0
            for the partner system  H_1. 
            
            In the case of the DCKAT generated
            on the basis of the seed eigenstates
            ψ_0^-,…,ψ_n-1^-,ψ_n^-,
            n=1,…,
            we obtain a partner system with the added 
            n+1 bound states in the lower part of the spectrum
            of the QHO which is described by the shifted
            QHO Hamiltonian H_n+1=H_osc-2(n+1). 
            
            The described 
            DCKATs based on the choice of the seed eigenstates
            ψ_0,…,ψ_n or  ψ^-_0,…,ψ^-_n
            reflect the property of the special shape invariance of the 
            QHO system.
            
            0.1cm
            
            Before we proceed further, let us summarize  
            briefly
            the main features of the non-singular 
            DT  and DCKAT schemes
            based on other choices of the sets of 
            physical, ψ_n, and non-physical, ψ^-_n,
            eigenstates of the QHO as the seed eigenstates
            <cit.>.
            This will generalize  the preceding discussion and will
            allow us  to generate the REQHO systems. 
            
            0.1cm
            
            Consider the DCKAT generated on the basis of  the 
            physical eigenstates ψ_i_1,
            ψ_i_2, …, ψ_i_n, and non-physical 
            eigenstates ψ_j_1^-,
            ψ_j_2^-, …, ψ_j_l^-,
             of the QHO.
            The peculiarity of the physical and non-physical eigenstates in both
             families is that they have  
            a  form of polynomials  multiplied by exponential 
            functions e^-x^2/2 and e^x^2/2, respectively.
            As a result we obtain a quantum system described
            by a potential to be a rational function.
            In order 
            a new quantum system generated by the DCKAT be non-singular, 
            the complete set of the seed  states has to be composed from the 
            blocks of the states
            (ψ_0^-,…,ψ_j_1-1^-,ψ_j_1^-), (ψ_j_2^-,…,
            ψ_j_2+l_2^-),…,(ψ_j_r^-,…,
            ψ_j_r+l_r^-),
            (ψ_0,…,ψ_i_1-1,ψ_i_1),
            (ψ_i_2,…,
            ψ_i_2+2m^+_2+1),…,(ψ_i_s,…,
            ψ_i_s+2m^+_s+1),
            where j_2=j_1+2m^-_1+3, j_k+1= j_k+2m^-_k+3,
            i_1< i_2, i_k+2m^+_k+1<i_k+1,
            and m^+_k, m^-_k and l_k can take values 0,1,….
            The total number of such blocks can be arbitrary and 
            blocks that include the states ψ_0 and ψ^-_0 
            can be absent.
            Up to a possible  constant shift, 
            the  generated system will have the gapped spectrum 
            of the QHO with the deleted levels  appearing
            at the positions of energies 
             corresponding to  physical states ψ_i
            in these blocks, and with new, added energy levels 
            appearing at the positions of
            energies of non-physical eigenstates ψ^-_j.
            In other words, the inclusion of physical states ψ_i into the 
            generating set of the seed states eliminates the energy levels,
            while the inclusion of non-physical states ψ^-_j introduces 
            corresponding additional energy levels into the spectrum.
            Each gap in the spectrum of the resulting system 
            contains an even number of the missing energy levels. 
            As in the simplest examples we considered above
            with the partner system to be the same but the shifted 
            QHO, the same REQHO system
            can be produced by DCKATs based on different
            choices of the sets of the seed eigenstates. 
            Different choices of the seed states 
            generate different  intertwining operators
            which relate REQHO system with the 
            QHO. As a consequence, as we shall see, 
             there exist different 
            ladder operators for the 
            same REQHO, which possess  different 
            properties. 
            Below we first consider some simple concrete examples
            of the REQHO systems. This will allow us
            to investigate in detail the families
            of the DCKATs associated with a given REQHO system, to 
            identify different ladder operators, and 
            to study their properties as well as  
            to establish the relations between them.
            Then the results will be developed 
            for the  case of the REQHO systems of a general form. 
            
            
            
            
            § SIMPLEST REQHO SYSTEM AND ITS LADDER OPERATORS
            
            
            A  simplest REQHO system 
            can be produced by
            taking a nodeless non-physical eigenstate of 
            energy E=-5 of H_osc, 
            ψ^-_2=(1+2x^2)e^x^2/2 .
             
            The first order differential   operators 
            
            A^-≡ψ_2^-d/dx1/ψ^-_2=d/dx-x-4x/2x^2+1 ,
                 A^+=(A^-)^† 
            factorize the shifted QHO Hamiltonian 
            
            A^+A^-=-d^2/dx^2+x^2+5=H_osc+5≡ H .
            
            Their permuted product  generates a simplest REQHO system, 
             
            A^-A^+=-d^2/dx^2+x^2+3 +82x^2-1/(2x^2+1)^2≡H̆ .
            
            We have the  intertwining relations
            
            A^-H=H̆A^- ,    
            A^+H̆=H A^+ ,
            
            from which it follows that the systems H and H̆ are almost 
            isospectral, and  operators (<ref>) 
            provide a  map  between the eigenstates 
            of the QHO and the REQHO systems.
            The excited eigenstates 
            of H̆ are the bound states 
            Ψ_n(x)=A^-ψ_n(x) ,    
            E_n=6+2(n-1) ,    
            n=1,2,…,
            
            where ψ_n(x) are the QHO eigenstates (<ref>).
            In correspondence with the general relation (<ref>),
            the ground state and its energy are
            Ψ_0=A^-ψ^-_2=1/ψ^-_2 ,    
            E_0=0 .
            
            In correspondence with
            a general picture described above, this state constitutes the 
            kernel of the operator A^+.
            The ladder operators for H̆ can be constructed 
            by the Darboux-dressing 
            of the ladder operators a^± of the QHO,
            𝒜^±=A^-a^± A^+ .
            
            We have 
            
            [H̆,𝒜^±]=± 2𝒜^± ,
            
            and
            𝒜^+𝒜^-=H̆(H̆-2)(H̆-6) ,    𝒜^-𝒜^+=H̆(H̆+2)(H̆-4).
            
            Due to the last factor in (<ref>) and in correspondence 
            with relations (<ref>), 
            the ground-state of
            the REQHO of zero energy, Ψ_0=1/ψ_2^-,
             is annihilated 
            by both ladder operators 𝒜^- and 𝒜^+.
            Another peculiarity is 
            that the kernel of the lowering operator 
            𝒜^- also contains 
            the first excited physical state
            Ψ_1=A^-ψ_0 of energy E=6,
            and the non-physical state
            A^-ψ^-_1 which is the eigenstate of H̆ 
            of the eigenvalue E=2.
            The three-dimensional  kernel
            of the  lowering ladder operator 𝒜^-
            is therefore
             (𝒜^-)=span {Ψ_0,A^-ψ^-_1,
            Ψ_1} .
            
            Besides the 
            ground state Ψ_0,
            the kernel of 𝒜^+   includes 
            the  two
            states A^-ψ^-_3 and A^-ψ^-_0,
            which are non-physical  eigenstates of H̆ 
            of the eigenvalues -2 and 4,
             (𝒜^+)=span {Ψ_0,A^-ψ^-_3,
            A^-ψ_0^-} .
            
            
            We denote (α_1) this scheme based on the DT with
            generating function ψ^-_2, 
            
            (α_1)={ψ^-_2}.
            
            Up to a global shift,  the same REQHO system can also 
            be produced
            by means of any of the DCKAT  schemes
            
            (α_2)={ψ^-_0,ψ^-_3} , …, 
            (α_n+1)={ψ^-_0,ψ^-_1,…,
            ψ^-_n-1,ψ^-_n+2} .
            
            The presence of the first n states 
            ψ^-_0,ψ^-_1,…,
            ψ^-_n-1
            in the scheme
            (α_n+1) gives rise to 
            the addition  of the corresponding energy levels 
            into the spectrum of the QHO while the inclusion 
            of the state ψ^-_n+2 results finally in the generation  
            of the shifted gapped REQHO system H_n+1=H̆-2n.
            The intertwining operators between 
            H defined by (<ref>) and H_n+1  in this case are 
            _n+1=A^-(a^+)^n ,    _n+1^†=(a^-)^nA^+ ,
              
             where  A^- and A^+ are 
             given by Eq. (<ref>).
             Combining the intertwining operators of the schemes
             (α_n+1) and (α_1),
             we can construct the 
             higher-order  ladder operators
             
              A^-_n+1^† =A^-(a^-)^nA^+≡𝒜_n^- ,    _n+1A^+=A^-(a^+)^nA^+=(𝒜_n^-)^†≡𝒜_n^+ ,
             [H̆, 𝒜_n^±]=± 2n 𝒜_n^±,
             where 𝒜^±_1=𝒜^±.
             This in particular
              means  that the third-order differential  operators
             (<ref>), which have the nature of the 
             Darboux-dressed QHO operators a^±, 
             can also be considered as the ladder operators
             generated via a composition of the intertwining 
             operators corresponding to the schemes 
             (α_1) and (α_2),
               𝒜^-=A^-(a^-A^+),
             𝒜^+=(A^-a^+)A^+.
             The use of the scheme (α_n) with n>2 
             instead of (α_2) in such a composition 
             provides  us therefore with the  
             (α_1)-Darboux-dressed form 
             (<ref>)
             of the higher-order
             ladder operators (a^±)^n of the QHO.
            The following relations 
            between 𝒜^- and 𝒜^-_n  are valid, 
              
             (𝒜^-)^n=∏_j=1^n-1(H̆+2j)·𝒜^-_n=
              𝒜^-_n·∏_j=1^n-1(H̆-2j) ,
             
            and analogous relations are obtained from them 
            for 𝒜^+ and 𝒜^+_n 
            by  the Hermitian conjugation.
            In 
            a more general case the composition of the intertwining operators of the 
            schemes 
            (α_n) and (α_m) with n>m generates 
            the higher-order ladder operators 
            𝒜^±_n-m,
            _m+1_n+1^†=∏_j=1^m(H̆-4-2j)·𝒜_n-m^- ,
                _n+1_m+1^†= 𝒜_n-m^+·∏_j=1^m(H̆-4-2j) .
            
            
            
            The REQHO system (<ref> ) 
            can also be generated 
            via the DCKAT based on  the 
            physical eigenstates ψ_1 and ψ_2.
            We denote this scheme, 
            which eliminates two neighbour energy levels
            E=8 and E=10 
            in  the spectrum of the shifted 
            QHO (<ref>),
             as (β_2):
            
            (β_2)={ψ_1,ψ_2} .
            
            We denote  ^±_2 the second-order
             intertwining operators _2 and _2^†
            constructed 
            on the basis of these two states according to (<ref>), 
            _2^-≡ A^-_isoa^-_iso ,    _2^+≡ (_2^-)^†=a^+_isoA^+_iso .
            
            The operator a^-_iso is defined in Eq. (<ref>),
            and a^+_iso=(a^-_iso)^†.
            The result of the action of the 
            operator a^-_iso on the
            QHO's eigenstate ψ_2 can be presented 
            in terms of its physical and non-physical  eigenstates in the 
            form
            a^-_isoψ_2=-ψ_0ψ^-_2/ψ^-_1≡ϕ.
            The 
            first-order differential operators 
            A^±_iso  are  generated by the function ϕ(x),  
            according to  (<ref>),
            
            A^-_iso=ϕ(x)d/dx1/ϕ(x)=d/dx
            +x+1/x-4x/1+2x^2 ,    
            A^+_iso=(A^-_iso)^† .
            
            We have 
            
            a^+_isoa^-_iso=H-8 ,    
            a^-_isoa^+_iso=H_iso ,
            
            where 
            
            H_iso=-d^2/dx^2+x^2+2/x^2-1
            
            is the shifted isotonic oscillator 
            to be singular at x=0,
            and H corresponds to the shifted
            Hamiltonian of the QHO defined in (<ref>).
            We also have the relations
            
            A^+_isoA^-_iso=H_iso-2 ,    
            A^-_isoA^+_iso=H̆-4 ,
            
            where H̆ is the Hamiltonian 
            of the REQHO system defined in Eq. (<ref>).
            From (<ref>) and (<ref>)
            we find
            _2^- H=(H̆+6)_2^- ,    _2^+ (H̆+6)=H_2^+ ,
            
            and 
            _2^+ _2^-=(H-8)(H-10) ,    _2^-_2^+=(H̆-2)(H̆-4) .
            
            By the construction, (_2^-)=span {ψ_1,ψ_2}. 
            One can also see that 
            (_2^+)=span {A^-ψ^-_0,A^-ψ^-_1}=
            span {^-_2ψ_1,
            ^-_2ψ_2}=
            span {ψ_1/_2,ψ_2/_2},
            where _2(x)=(ψ_1,ψ_2)(x)=-ψ^-_2(x)  e^-3/2x^2.
            0.1cm
            
            The DCKATs corresponding to the 
            (α_1)- and (β_2)-schemes  are
            in some sense complementary. The (α_1)-scheme introduces   
            effectively a new energy level into the 
            spectrum of the QHO  below its ground-state
            energy at the  distance equal to the tripled
            distance between equidistant energy levels.
            The (β_2)-scheme 
            makes a similar job but by deleting the 
            first two excited energy levels in the spectrum 
            of the QHO. Since the Wronskian 
            in the DCKAT in the 
            latter scheme includes the additional 
            exponential factor e^-3/2x^2 in comparison
            with the structure of the non-physical eigenstate
            ψ^-_2,  this 
            produces the additional constant shift +6 
            in the potential generated by means of 
            relation (<ref>)
            that is reflected in  Eq. (<ref>),
            cf. Eq. (<ref>).
            From Eq. (<ref>) it follows that if ψ(x;E)
            is an eigenstate of H of energy E, 
            then ^-_2ψ(x;E) is an 
            eigenstate of H̆ of energy (E-6).
            As a consequence,   
            all the spectrum of the system generated  by 
            the (β_2)-scheme will be shifted 
            for  -6 in comparison 
            with the spectrum of the 
            REHQO system (<ref>) produced via
            the (α_1)-scheme.
            In correspondence with this picture, 
            the ground-state 
            Ψ_0 of H̆ can alternatively
            be constructed from the QHO
            ground-state ψ_0, Ψ_0=^-_2ψ_0,
            cf. (<ref>).
            The excited states Ψ_n+1 of
            energy E_n+1=6+6n with n=0,… 
            can be presented in the alternative to  (<ref>) form
             Ψ_n+1=^-_2ψ_n+3.
            All this gives a possibility 
            for the construction of another pair of  ladder
            operators for the REQHO system (<ref>)
            with the properties rather different to
             those of the ladder operators we obtained
            by using only the (α)-schemes. For this we take now
            the composition of the intertwining operators
            of the (α_1)- and (β_2)- schemes
            to construct the operators
            𝒞^-=_2^-A^+ ,    𝒞^+=A^-^+_2 .
            
            Instead of (<ref>) and (<ref>),
            they satisfy the relations 
            
            [H̆,𝒞^±]=± 6 𝒞^± ,
            
            and 
            𝒞^+𝒞^-=H̆(H̆-8)(H̆-10) ,    𝒞^-𝒞^+=(H̆+6)(H̆-2)(H̆-4) .
            
            Like 𝒜^±, these are  third-order
             differential
            operators of the nature of  ladder operators.
            However, acting on eigenstates of the REQHO system 
            H̆, they change the energies not in 2 but in 6.
            In this aspect they are somewhat  similar to the 
            higher-order ladder operators 𝒜^±_3,
            which are
            fifth-order differential operators discussed
            above. The essential
            difference of 𝒞^- 
            from 𝒜^-_n and in particular 𝒜^-
             is that 
            in correspondence with the
            first relation from (<ref>),
            the kernel of 𝒞^- 
            is composed only from physical eigenstates
            of H̆,  
             (𝒞^-)=span {Ψ_0,Ψ_2,Ψ_3} .
            
            The energies 0, 8 and  10 of these states are the roots
            of the third-degree polynomial in the first identity in (<ref>). 
            Also, unlike 𝒜^+,
             the kernel of the raising operator 𝒞^+
            is composed only from the
            non-physical eigenstates of H̆, 
             (𝒞^+)=span {A^-ψ^-_5,A^-ψ^-_1,A^-ψ^-_0}=
            span {^-_2ψ^-_3,^-_2ψ_1,^-_2ψ_2}.
            
            In correspondence with the second 
            relation in (<ref>),
            the eigenvalues of the states in (<ref>) are
            -6, 2 and  4.
            From the point of view of the structure
            of the kernels and commutation relations
            (<ref>), 
             the ladder operators
            𝒞^± are 
            similar to the 
            third-order  differential 
            operators (a^±)^3 in the QHO system.
            However, unlike 𝒞^-,  the operator  
            (a^-)^3  annihilates the three lowest physical 
            eigenstates of 
            the QHO of the three subsequent  values of energy. 
            The first exited state Ψ_1 of the REQHO 
            system of energy E=6 does not belong  to 
            the kernel of 𝒞^- and is annihilated by 
            (𝒞^-)^2: 
            𝒞^-Ψ_1=Ψ_0, (𝒞^-)^2Ψ_1=0.
            
            The following relations can be 
            established by comparing  
            the kernels of the operators on both sides of
            the  equalities,
            
            A^+^-_2=-(a^-)^3 ,    ^+_2A^-=-(a^+)^3 .
            
            These  and their analogous  relations  for
            other REQHO systems will play important role
            in what follows. From them one can find 
            in particular the operator identities 
            a^-^+_2=-(a^+)^2A^+,
            (a^-)^2^+_2=-a^+A^+(H̆-2),
            as well as the Hermitian conjugate ones.
            
            
             
             One can introduce  additionally the operators
             a^± inside the factorized 
             structure of the  operators 𝒞^±.
             In this way one can construct the operators
             𝒞^-_n+1=_2^-(a^-)^nA^+ ,    𝒞^+_n+1=A^-(a^+)^n^+_2 ,  
              n=0,… , 
              
              with the implied identification  𝒞^±_1=𝒞^±
              for n=0. 
              They satisfy the relation [H̆,𝒞^±_n+1]=
              ± (6+2n)𝒞^±_n+1.
              These operators can be treated either as 
              the QHO operators (a^±)^n dressed
              by the intertwining generators 
              of the (α_1) and (β_2) schemes,
              or as the operators produced by
              intertwining operators (<ref>) from the
              (α_n+1) scheme and those
              from the same (β_2) scheme. 
              The kernel of 𝒞^-_n+1 is 
             composed only by
            the physical eigenstates of H̆, while
            the kernel of   𝒞^+_n+1 is 
            spanned only by its 
            non-physical eigenstates.
            For instance, (𝒞^-_2)=span {Ψ_0,Ψ_1,Ψ_3,Ψ_4}.
            With the help of identities (<ref>) we also find that
            (𝒞^±)^n=(-1)^n+1𝒞^±_3(n-1)+1, 
            n=1,… .
            
            
            Analogously to 𝒜^±, we also
            introduce   the operators
            ℬ^±=^-_2 a^±^+_2 .
            
            Unlike the third-order ladder operators 
            𝒜^± and 𝒞^±, 
            the ℬ^±
            are fifth-order differential operators. 
            They satisfy the relations 
            
            [H̆,ℬ^±]=± 2ℬ^±
            and 
            ℬ^+ℬ^-=
            H̆(H̆-2)(H̆-6)(H̆-4)^2 ,    ℬ^-ℬ^+=
            H̆(H̆+2)(H̆-4)(H̆-2)^2 .
            
            The kernel of ℬ^- involves
            two  physical 
            and three non-physical eigenstates 
            of H̆,
             (ℬ^-)=span {Ψ_0,Ψ_1,A^-ψ^-_0,A^-ψ^-_1,A^-ψ^-_0} .
            
            The eigenvalue E=4 of the 
            non-physical eigenstates 
            A^-ψ^-_0 and A^-ψ^-_0 
            in the kernel of 
            ℬ^- 
            corresponds to the
            double root of the last factor in the first relation in  
            (<ref>). 
             The kernel of the increasing ladder operator
             ℬ^+
             includes only one physical eigenstate,
             (ℬ^+)=span {Ψ_0,A^-ψ^-_0,A^-ψ^-_1,A^-ψ^-_3,A^-ψ^-_1} .
            
            The eigenvalue E=2 of the 
            non-physical eigenstates 
            A^-ψ^-_1 and A^-ψ^-_1 
            in the kernel of the increasing  ladder operator 
            ℬ^+
            corresponds to the
            double root of the last factor in the 
            second relation in  
            (<ref>). 
            By analogy with (<ref>), 
            one can consider the 
            higher-order ladder operators
            ℬ^±_n=^-_2(a^±)^n^-_2 ,
            [H̆,ℬ^±_n]=± 2n
            ℬ^±_n,
            with the identification ℬ^±_1=ℬ^±.
            The  lowering operator ℬ^-  can be   related 
            to the ladder operators 𝒜^-
            and 𝒞^-
            via the identities 
            ℬ^-=𝒜^-(H̆-4) ,    ℬ^-_2=(𝒜^-)^2 ,    ℬ^-_3 =-𝒞^-(H̆-2)(H̆-4) ,
            ℬ^-_4=-𝒞^-_2(H̆-2)(H̆-4),
            ℬ^-_5=-𝒞^-(𝒜^-)^2,
            ℬ^-_6=(𝒞^-)^2(H̆-2)(H̆-4),
            etc.  The  increasing operator 𝒞^+  is related 
            to  𝒜^+
            and 𝒞^+ via   the conjugate identities.
            Similarly to (<ref>), for degrees n>1 
            of ℬ^- we have 
            
            (ℬ^-)^n=∏_j=1^n-1(H̆-2+2j)
            (H̆-4+2j)·ℬ^-_n ,
            
            and an analogous  relation for (ℬ^+)^n.
            
            
             A generalization of the (β_2)-scheme 
              corresponds to the family of the DCKAT schemes
               (β_3)={ψ_0,ψ_2,ψ_3} , …,
            (β_n+2)={ψ_0,…,ψ_n-1,ψ_n+1,ψ_n+2} .
            
            In the case of the scheme 
            (β_n+2), the intertwining operators  
            constructed according to the prescription (<ref>) 
            are _n+2^-≡^-_2(a^-)^n
            and ^+_n+2=(a^+)^n_2^+.
            The   schemes (β_n+2) 
            do not give anything new but 
            allow us to re-interprete 
            the already discussed higher-order
            ladder operators 𝒞^±_n+1
            as those produced
            via the composition of the (β_n+2) and (α_1) schemes,
            𝒞^-_n+1=^-_n+2A^+, 
            𝒞^+_n+1=A^-^+_n+2.
            Analogously, ℬ^-_n+1=^-_n+2^+_2,
            ℬ^+_n+1=^-_2^+_n+2.
            
            Besides the two infinite families (α_n),
            n=1,…,  and (β_n), n=2,…,
            which involve  as the seed eigenfunctions
             either only  non-physical 
            or only physical eigenstates of the QHO,
            there are two additional, `intermediate' schemes 
            which simultaneously include  eigenstates  of   
            both types. These are 
            the schemes
            
            (γ_2)={ψ_0,ψ^-_1} ,    
            (γ_3)={ψ_0, ψ_1,ψ^-_0} .
            
            The intertwining operators 
            in the case of the scheme (γ_2)
            are the second-order differential operators
            _2=A^-a^- and _2^†=a^+A^+, while 
            in the scheme  (γ_3),
            the intertwining operators are 
            the third-order  differential operators
            _3=A^-(a^-)^2 and _3^†=(a^+)^2A^+.
            These operators have a structure similar 
            to that of  the intertwining operators in 
            the family of the schemes (α_n).
            With their help we do not 
            obtain anything essentially new for the construction 
            of the ladder operators for the REQHO system
            H̆ in comparison 
            with the already discussed structures.
            Indeed, employing the superposition 
            of the intertwining operators from the (β_2)-scheme
            and either (γ_2)- or (γ_3)-  schemes, 
            one can construct the ladder operators 
            𝒞^-_-n≡^-_2(a^+)^n A^+ ,    𝒞^+_-n≡ A^- (a^-)^n  ^+_2 ,    
            n=1,2 .
            
            Here  n=1 and n=2 correspond, respectively,
            to the (γ_2)- and  (γ_3)- schemes.
            Operators 𝒞^±_-1 are 
            fourth-order differential operators,  
            while 𝒞^±_-2 are  fifth-order differential 
            operators.
            They, however,  are not independent 
            but  can be expressed in terms of 
            the already constructed intertwining operators.
            Namely, we have, in particular,  
            𝒞^±_-1=-𝒜^±_2 ,    𝒞^-_-2=-(H̆-2)𝒜^- ,    𝒞^+_-2=-(H̆-4)𝒜^+ .
            
            These relations can be established by  comparing  the  
            kernels of  𝒞^-_-1 and 
            𝒞^-_-2,
            (𝒞^-_-1)=span {Ψ_0,A^-ψ^-_0,Ψ_1,Ψ_2}, 
            (𝒞^-_-2)=span {Ψ_0,A^-ψ^-_1, A^-ψ^-_0,
            A^-ψ^-_0, 
            Ψ_1},
            with the kernels of 
            𝒜^-_2
            and of the operator 
            (H̆-2)𝒜^-=𝒜^-(H̆-4), respectively,
            and by comparison of the signs before the leading 
            derivative terms in the corresponding pairs of
            operators.
            Due to the identities  (<ref>), 
            a generalization of the operators (<ref>)
            for n>2  does not give us 
            anything new  since 
            𝒞^±_-3=-H̆(H̆-2)
            (H̆-4).
            This last relation as well as relations 
            (<ref>) can also be obtained by
            employing the identities (<ref>).
            We also  have the operator identities 
            𝒜^+𝒞^+_n=(H̆-2)𝒞^+_n+1 ,    𝒜^-𝒞^-_n=-𝒜^-_n+3 ,    
            n=1,…,
            𝒜^+𝒞^-=-(H̆-6)𝒜^-_2 ,    𝒜^-𝒞^+=-(H̆+2)𝒜^+_2 ,
            
            as well as the Hermitian conjugate 
            relations.
            
            Let us look  in more detail at 
            the already 
            mentioned similarity between the operators
            𝒞^± and 𝒜_3^±.
            Using the first identity from (<ref>) and 
            factorization relation (<ref>), we obtain 
             𝒜^-_3=A^-(a^-)^3A^+
            =-H̆𝒞^- and 
            𝒜^+_3=-(H̆-6)𝒞^+.
            These relations are similar to those in (<ref>).
            Employing Eq. (<ref>),
            we find  that 
            
            (𝒜^-)^3=-H̆
            (H̆+2)(H̆+4)𝒞^-=
            -𝒞^-(H̆-6)
            (H̆-4)(H̆-2) .
            
            This relation from the point of view 
            of the kernels of the involved 
            operators  corresponds to the following picture.
            The kernel of the operator 𝒜^- 
            is formed by  the two physical eigenstates
            Ψ_0 and Ψ_1
            and by one non-physical eigenstate A^-ψ^-_1.
            We have also the relations <cit.>𝒜^-(A^-ψ^-_1)=Ψ_0,
            𝒜^-(A^-ψ^-_0)=A^-ψ^-_1,
            𝒜^-Ψ_2=Ψ_1,
            𝒜^-(A^-ψ^-_0)=A^-ψ^-_1,
            𝒜^-(A^-ψ_0)=A^-ψ^-_0,
            𝒜^-Ψ_3=Ψ_2. As a result we obtain
             (𝒜^-)^3=span {Ψ_0,
            A^-ψ^-_1,
             Ψ_1, 
             A^-ψ^-_1,
            A^-ψ^-_0,
             Ψ_2, 
             A^-ψ^-_0,
            A^-ψ_0,Ψ_3} .
            
            The pairs  of states (A^-ψ^-_1, 
            A^-ψ^-_1), 
            (A^-ψ^-_0, ψ^-_0)
            and (Ψ_1=A^-ψ_0, A^-ψ_0)
            constitute, respectively,   
            the kernels of the factors (H-2),
            (H-4) and (H-6)
            in (<ref>).
            The remaining three physical eigenstates  
            Ψ_0, Ψ_2 and Ψ_3 in (<ref>)
            correspond to the kernel of the 
            operator 𝒞^-.
            
            According to  Eq. (<ref>) 
            and its conjugate version,  the ladder operators
            𝒞^± can be generated 
            by 𝒜^±.
            Then with taking into account Eq. (<ref>)
            and all the described relations,  
            we conclude that in the case of the simplest 
            REQHO system given by the Hamiltonian
            H̆ defined in (<ref>)
            all the set of the  ladder operators
            can be obtained, in principle, 
            from the  compositions of  
            the 
            ladder
            operators ℬ^+
            and  ℬ^-.
            
            In conclusion of this section 
            let us note, however,  that in comparison with the 
            QHO picture, the peculiarity of the system
            (<ref>) in particular is that
            its  ground-state Ψ_0 
            cannot be achieved 
            from physical states
            by action of
            the lowering  operators 𝒜^- and ℬ^-
            which are  differential operators of orders 3 and 5.  
            Like the first-order differential operator
            a^- in the QHO,  the ladder operators
             𝒜^- and ℬ^- 
            decrease the energy values of H̆ 
             in 2, but they produce 
            the ground-state by acting on the 
            non-physical eigenstate A^-ψ^-_1
            of the eigenvalue E=2.
            The ground-state Ψ_0 of zero energy can be achieved,
            however,
            by application of the lowering operator
            𝒞^-, which is a third-order differential operator,
              to the physical eigenstate
            Ψ_1  with eigenvalue E=6.
            
            
            We also notice here  that the ladder operators 
            (<ref>) for the REQHO system 
            of the simplest form  (<ref>) were constructed 
            (without  employing   the Darboux-dressing procedure)  in  
             <cit.> where this model was introduced and 
             investigated  for the first time.
            Later these ladder operators were constructed,
            particularly,  in <cit.>,
            <cit.> and recently in <cit.>
            via the  Darboux-dressing prescription 
            based on the non-physical seed state
            we used here.
            The ladder operators (<ref>) we constructed
            on the basis of the (α_1)- and (β_2)- schemes 
            by employing the analogy with  reflectionless 
            quantum systems <cit.>
            where the corresponding Lax-Novikov integrals of motion can be
            generated either by Darboux-dressing of the
            free particle momentum operator or by `gluing'
            two different intertwining 
            operators that act in the opposite directions. 
            The same last mentioned method 
            also allows ones to generate the Lax-Novikov integrals 
            for periodic finite-gap systems where
            the Darboux-dressing mechanism  can not be applied, see  
            <cit.>. 
            In reflectionless and finite-gap systems, 
            however, the 
            two glued  intertwiners always 
            are differential operators of the 
            `opposite', even and odd, differential orders,
            but both intertwine the two corresponding 
            partner systems
            without additional relative displacement.
            It is because of the relative  displacement 
            in intertwining relations (<ref>) 
            and  (<ref>) that here we obtain the
            ladder operators 𝒞^± 
            while in reflectionless  and finite-gap systems
            analogous procedure generates the integrals 
            of motion.
            Within the same framework we used here and 
            based on  employing the
            two schemes with physical and non-physical seed states
            of the quantum harmonic oscillator,  
            the ladder operators  (<ref>) were introduced
            earlier in <cit.> 
             (but without exploiting the 
            indicated analogy with generation of the Lax-Novikov integrals)
            and derived later in <cit.> 
            for  more general families related to 
            multi indexed exceptionnal orthogonal polynomials.
            By another method such ladder operators were
            introduced  for the system (<ref>) 
            even earlier in <cit.>. 
            The ladder operators (<ref>)  constructed here
            by the Darboux-dressing procedure based on physical seed states 
            seems were not discussed earlier in the literature.
            
            With subsequent  analysis we shall see that the 
            trinity 
            of the basic ladder operators (𝒜^±,
            ℬ^±, 𝒞^±) 
            admits a natural generalization for the case of
             REQHO systems of a general form, and 
            that each pair of the conjugate lowering and raising 
            ladder operators  
            detects and reflects some specific properties 
            of a corresponding quantum system.
            
            
            
            § TWO FURTHER EXAMPLES OF THE REQHO SYSTEMS
            
            
            A  REQHO system generated by 
            the  DT
            based on the non-physical 
            state ψ^-_2n with n>1 is 
            similar to the considered 
            REQHO system generated by
            the DT based on ψ^-_2.
            In this case the gap in the spectrum 
            of the REQHO corresponds,
            up to a global shift,  to the missing 
            2n energy levels  with  E=3,…,4n+1
            in the spectrum of the QHO. 
            We also have here the two infinite families of
            the DCKAT schemes
            of the structures which generalize those of the 
            case n=1.  For instance, in the case of n=2, 
            we have ψ^-_4=(4x^4+12x^2+3)e^x^2/2, and
            the two infinite families of the schemes are
            
            (α_1)={ψ^-_4} ,    
            (α_2)={ψ^-_0,ψ^-_5} ,    
            (α_n+1)={ψ^-_0,…,ψ^-_n-1,ψ^-_n+4} ,
            
            and
            
            (β_4)={ψ_1,ψ_2,ψ_3,ψ_4} ,    
            (β_n+4)={ψ_0,…,ψ_n-1,
             ψ_n+1,ψ_n+2,ψ_n+3,ψ_n+4} .
             
             In addition, we have the `intermediate' schemes whose sets 
             of seed states include both physical and non-physical 
             eigenstates of the QHO.
             These are
             
             (γ_2)={ψ_0,ψ^-_3} ,
               
              (γ_3)={ψ_0,ψ_1,ψ^-_2} ,
               
             (γ_4)={ψ_0,ψ_1,ψ_2,ψ^-_1} ,
               
             (γ_5)={ψ_0,ψ_1,ψ_2,ψ_3,ψ^-_0} .
              
            The scheme (α_1)  generates
            the intertwining operators A^-=ψ^-_4d/dx1/ψ^-_4,
            and A^+=(A^-)^†.
            They allow us to construct ladder operators that are third-order differential
            operators,  the Darboux-dressed 
            ladder operators of the QHO, 𝒜^±= A^-a^± A^+.
            They satisfy the relations of the form 
            (<ref>), 
            [H̆,𝒜^±]=± 2𝒜^±,   with
            H̆≡ A^-A^+=
            -d^2/dx^2+x^2+7+324x^6+4x^4+3x^2-6/(4x^4+12x^2+3)^2 ,
            
            and similarly to the 
            already considered case,
            here both ladder operators 𝒜^± annihilate the
            ground state of H̆ which is Ψ_0=1/ψ_4^-=A^-ψ^-_4.
            Besides the ground state Ψ_0 
            of  energy E=0,  the kernel of 𝒜^- contains, 
            the first excited state Ψ_1=A^-ψ_0 of  
            energy E=10
            and one non-physical eigenstate A^-ψ^-_3 of 
            energy E=2.
            The kernel of 𝒜^+ contains besides the ground state 
            Ψ_0 also two non-physical eigenstates 
            A^-ψ^-_0 and A^-ψ^-_5
            of H̆ of eigenvalues E= 8 and E=-2. 
            The intertwining operators corresponding to the 
            DCKAT scheme (β_4)={ψ_1,ψ_2,ψ_3,ψ_4} are the 
            fourth order differential operators
            constructed in accordance with 
            Eq. (<ref>),
            which by analogy with the already 
            considered case 
            we denote here as
            ^-_4 and ^+_4=(^-_4)^†.
            Then we define another pair of   ladder operators
            via a composition of the intertwining  operators
            of this (β_4)-scheme and of the 
            (α_1)-scheme,
            𝒞^-=^-_4 A^+, 𝒞^+=
            A^-^+_4.
            Unlike the previously discussed case of the REQHO system
            (<ref>),
            these are fifth-order differential operators,
            which satisfy the relations [H̆,𝒞±]=
            ± 10𝒞^±.
            The ladder operator 𝒞^- annihilates
            five physical eigenstates of H̆,
            which are the ground state 
            Ψ_0 and the states Ψ_j+1=A^-ψ_j, j=1,2,3,4,
            with the  energy values E_0=0 and
            E_j+1=10+2j.
            Like in the REQHO system we considered before, 
             the first excited state Ψ_1=A^-ψ_0 of energy
            E_1=10 here does not belong to the kernel of 
            the decreasing ladder operator 𝒞^-,
            and we have 𝒞^-Ψ_1=ψ_0, 
            (𝒞^-)^2Ψ_1=0.
            The kernel of 𝒞^+
            is composed only from non-physical eigenstates.
            Yet another pair of the ladder operators corresponds to
            ℬ^±=^-_4a^±^+_4,
             which are differential operators of order 9.
             Like 𝒜^±, they satisfy the 
             relations 
             [H̆,ℬ^±]=± 2ℬ^±.
             The kernel of the  lowering ladder operator 
             ℬ^- is spanned by two physical eigenstates 
             Ψ_0 and Ψ_1
             of energies 0 and 10, and seven non-physical eigenstates
             of H̆ of eigenvalues 8 (twice), 6 (twice), 4 (twice)
             and  2,
              (ℬ^-)=span {Ψ_0,Ψ_1,A^-ψ^-_0,
             A^-ψ^-_0,Ψ_1,A^-ψ^-_1,
             A^-ψ^-_1,
             Ψ_1,A^-ψ^-_2,
             A^-ψ^-_2,A^-ψ^-_3}.  
             The kernel of ℬ^+ is spanned by the 
             ground-state Ψ_0 and by eight non-physical eigenstates
             of H̆.
            
            Other,  
            secondary ladder  operators can be constructed 
            by introducing the QHO ladder operators (a^±)^n
            inside the factorized structures 
            of the basic ladder operators
            𝒜^±, ℬ^±
            and 𝒞^±,
            or by considering compositions
            of the intertwining operators corresponding to
            (<ref>), (<ref>) and (<ref>) schemes 
            analogously to
             how it was done for the simplest REQHO system.
             The  secondary, higher-order ladder operators can also be
             generated by taking the products of  the 
             basic ladder operators.
             Analogously to (<ref>), we also have here 
             the relations 
             A^+^-_4=-(a^-)^5, 
            ^+_4A^-=-(a^+)^5.
            As an analog of relation (<ref>)
            we have 
            
            (𝒜^-)^5=
            -𝒞^-(H̆-10)
            (H̆-8)(H̆-6)(H̆-4) ,
            
            and the relation 
            ℬ^-=(H̆-2)(H̆-4)(H̆-6)𝒜^-
            
            is analogous here to the first relation in (<ref>). 
            As in the case of the simplest REQHO system
            considered in the previous section,
            for the  REQHO system described by the Hamiltonian(<ref>)  
            all the ladder operators can be generated and extracted 
            from the powers of the basic ladder operators ℬ^±.
             0.1cm
             
            Let us consider yet  another  example of the REQHO system
            in which  two states are separated by a gap 
            from the infinite equidistant part of the spectrum.
            A simplest system of such a nature
             can be
             generated by employing  the  minimal (α)-scheme 
            
            (α_2)={ψ^-_2,ψ^-_3} .
            
            Let us shift the Hamiltonian of the QHO for +7 
            and denote  H=H_osc+7,  for which the potential is 
            V(x)=x^2+7 and the spectrum is E_n=8+2n, n=0,1,….
            The Wronskian here is
             _2(x)=e^x^2(3+4x^4).
            The second order DCKAT based on (<ref>) produces 
            the partner system which is  the REQHO system 
            H̆=-d^2/dx^2+V̆(x),
            where in accordance with (<ref>),
            V̆(x)=3+x^2
            +32x^24x^4-9/(3+4x^4)^2 .
            
            Its gapped spectrum is 
             
            E_0=0 ,   
            E_1=2 ,   
            E_2+n=8+2n ,   
            n=0,1,… .
            
            The intertwining second order differential operators
            constructed via (<ref>) on the basis of the seed 
            QHO eigenstates (<ref>) 
             we denote as _2^- and
            _2^+=(_2^-)^†.
            The kernel of ^-_2 is spanned by the states
            ψ^-_2 and ψ^-_3, 
            while the kernel of ^+_2 is spanned by 
            the lowest physical eigenstates of H̆
            of the energies 0 and 2,
            which can be obtained from 
            the QHO non-physical eigenstates ψ^-_2
            and ψ^-_3,
            Ψ_0=^-_2ψ^-_2,
            Ψ_1=^-_2ψ^-_3.
            The operators ^+_2 and ^-_2 
            satisfy the relations
            ^+_2^-_2=H(H-2) ,    _2^-^+_2=H̆(H̆-2) ,
            
            and
            ^-_2H=H̆^-_2,
            ^+_2H̆=H_2^+.
            We construct the ladder operators 𝒜^±  for H̆
            by the Darboux-dressing of the QHO operators a^±,
            𝒜^± =
            ^-_2 a^±_2^+ .
            
            These fifth-order differential operators satisfy 
            the relations 
            [H̆,𝒜^±]=± 2𝒜^±,
            and 
            𝒜^+𝒜^-=H̆(H̆-2)^2(H̆-4)(H̆-8) ,    𝒜^-𝒜^+=(H̆+2)(H̆)^2(H̆-2)(H̆-6) .
            
            The kernel of the lowering operaror is
            (𝒜^-)=span {Ψ_0,Ψ_1,Ψ_2, Ψ_1, _2^-ψ^-_1}.
            Here the first three states are the lowest three 
            physical eigenstates of H̆ of 
            energies E=0, 2 and 8,
            respectively, and the two last states are non-physical eigenstates 
            of energies E=2 and 4. 
            The indicated energies correspond 
            to the roots of the polynomial in the first equality in (<ref>).
            The kernel of the increasing operator 
            𝒜^+ is spanned by the two lowest eigenstates
            Ψ_0 and Ψ_1 and by the three non-physical eigenstates:
            (𝒜^+)=
            span {Ψ_0,Ψ_1,^-_2ψ^-_4,Ψ_0,^-_2ψ^-_0}.
            The energies of these five eigenstates 
            correspond to zeros of the polynomial in 
            the second equality in  (<ref>).
            
            The same system, up to a global shift, 
            can also be generated 
            via  the complementary minimal (β_2)-scheme
            
            (β_2)={ψ_2,ψ_3} .
            
            The Wronskian of the seed states in this case 
            is 
            _2(x)=(ψ_2,ψ_3)=
            e^-x^2(3+4x^4).
            The potential calculated according to (<ref>)
            is the potential (<ref>) shifted for +8.
            Let us denote the corresponding intertwining second order
            differential operators constructed on the basis of these seed states
            as ^-_2 and ^+_2=(^-_2)^†.
            They satisfy the relations 
             ^-_2H=(H̆+8)^-_2 ,    ^+_2(H̆+8)=H^+_2 ,
             
             and
             ^+_2^-_2=(H-12)(H-14) ,    ^-_2^+_2=(H̆-4)(H̆-6) .
            
            Note that here 
            _2^-(x) =-_2^-(ix).
            The kernel of ^-_2 is spanned by the seed states
            (<ref>), whereas the kernel of ^+_2 
            is spanned by non-physical 
            eigenstates of H̆ of energies
            E=4 and E=6:
            (^+_2)=span {^-_2ψ_2,^-_2ψ_3}=
            span {^-_2ψ^-_0,^-_2ψ^-_1}.
            The ladder operators 
            𝒞^-=^-_2^+_2 ,    𝒞^+=^-_2^+_2
            
             are differential 
            operators of the order four.
            They  obey the relations 
            
            [H̆,𝒞^±]=± 8 𝒞^± ,
            
            and 
            𝒞^+𝒞^-=H̆(H̆-2)
            (H̆-12)(H̆-14) ,    𝒞^-𝒞^+=(H̆+8)(H̆+6)
            (H̆-4)(H̆-6) .
            
            The lowering operator 𝒞^- annihilates
            the four physical states Ψ_0,Ψ_1,Ψ_4,Ψ_5
            of the energies
            E=0,2,12,14.
            The kernel of 𝒞^+ is spanned only by non-physical 
            eigenstates 
            of energies E=-8,-6,4,6,
            (𝒞^+)=span {^-_2ψ^-_3,^-_2ψ^-_2,^-_2ψ_2,
            ^-_2ψ_3}=
            span {^-_2ψ^-_7,^-_2ψ^-_6,^-_2ψ^-_1,
            ^-_2ψ^-_0}.
            Changing the order of the second order operators 
            in the factorized form of the ladder operators in 
            (<ref>),  we 
            obtain the operator identities
            ^+_2^-_2=(a^+)^4 ,    ^+_2^-_2=(a^-)^4 ,
            
            cf. (<ref>).
            Again, these relations  can be verified   by comparing the 
            kernels of the corresponding operators.
            
            Yet another pair of the ladder operators corresponds to 
            differential operators of order 5,
            ℬ^± =
            ^-_2 a^±_2^+ .
            
            They satisfy relations [H̆,ℬ^±]=
            ± 2ℬ^±, and 
            ℬ^+ℬ^-=H̆(H̆-4)(H̆-6)^2(H̆-8) ,    ℬ^-ℬ^+=(H̆+2)(H̆-2)(H̆-4)^2(H̆-6) .
            
            In correspondence with the first relation in (<ref>), 
            the kernel of the lowering operator is spanned 
            by the ground-state Ψ_0 and the eigenstate Ψ_2 at the bottom 
            of the equidistant infinite part of the spectrum as well as by 
            the three non-physical eigenstates
            of H̆,
            (ℬ^-)=span {Ψ_0,Ψ_2, _2^-ψ^-_0,
            _2^-ψ^-_1,_2^-ψ^-_0}.
            The kernel of the increasing operator ℬ^+
            is spanned by the first excited physical eigenstate
            Ψ_1 of energy 2, and by four non-physical eigenstates 
            of H̆ of eigenvalues  -2, 4 (twice) and 6,
            (ℬ^+)=span {Ψ_1,_2^-ψ^-_3,
            _2^-ψ^-_1,_2^-ψ^-_1,_2^-ψ^-_0}.
            Here ℬ^- and  𝒜^-
            are related by the operator identity
            H̆ℬ^-=(H̆-4)𝒜^- ,
            
            and the relation between the increasing ladder operators 
             ℬ^+ and  𝒜^+
             is given by Hermitian conjugation of 
             (<ref>). 
             Note that in comparison with the first relation 
             in (<ref>) and relation (<ref>)
             in the examples of the REQHO systems 
             with one separated energy level here 
             the relation (<ref>) contains 
             a Hamiltonian-dependent factor before the operator 
             ℬ^-.  Coherently with this,
             in the REQHO system under consideration 
             the lowering operator ℬ^- in comparison 
             with 𝒜^- annihilates only 
             the lowest state Ψ_0 
              in the separated part of the spectrum.
             The raising operator ℬ^+ annihilates 
             another separated state Ψ_1 in comparison 
             with both separated states Ψ_0 and Ψ_1 
             annihilated by 𝒜^+.
             This difference can be understood  if we 
             note that the polynomial in the first identity in (<ref>)
             does not have the root 2.
             Rewriting  relation (<ref>) in the equivalent form
             ℬ^-(H̆-2)=𝒜^-(H̆-6),
             we see then that the annihilation of the state
             Ψ_1  from the kernel of 𝒜^-
             is provided by the factor (H̆-2) on the left hand side 
             of the identity. In the same way one can understand 
             the difference in the kernels of the raising operators
              ℬ^+ and 𝒜^+ by looking at 
              the roots of the polynomial in the second identity 
              in (<ref>) and by taking into account 
              the identity relation 
              ℬ^+H̆=𝒜^+(H̆-4)
              to be conjugate to (<ref>).
             
               
            
            The same (up to a global shift) system H̆ given by 
            the potential (<ref>) can be produced by using the
            higher-order  (α)- and (β)-schemes,
            
            (α_n+2)={ψ^-_0,…,ψ^-_n-1,ψ^-_n+2,ψ^-_n+3} ,    
            (β_n+2)={ψ_0,…,ψ_n-1,ψ_n+2,ψ_n+3} ,
            n=1,…, 
            and by the  two intermediate 
            (γ)-schemes,
            
            (γ_3)={ψ_0,ψ^-_1,ψ^-_2} ,    
            (γ_4)={ψ_0,ψ_1,ψ^-_0,ψ^-_1} .
            
            The secondary, higher-order ladder operators can be generated here
            in the same way as  for the REQHO systems
            with one added gapped bound state.
            Here the relation
            
            (𝒜^-)^4=
            𝒞^-(H̆-8)^2
            (H̆-6)^2(H̆-4)^2(H̆-2)(H̆-10) .
            
            is analogous to the  relations (<ref>) and 
            (<ref>),
             and shows that the ladder operators
            𝒞^± can be  generated by
            the  ladder operators 𝒜^±.
            
            
            
            
            
            § REQHO SYSTEMS OF A GENERAL FORM AND THEIR LADDER OPERATORS
            
            
            We generalize now our analysis  of the particular examples
            for the case of  REQHO systems 
            of a  general form, for which 
            we construct    the ladder operators and investigate their
            properties.
            
            Each REQHO system can be generated 
            by employing 
            the appropriate DCKAT based on 
            any of the schemes from the two infinite families.
            The (α)-type schemes
            include only 
            non-physical eigenstates of the QHO   
            chosen as
            the seed states.
            The (β)-type schemes involve 
            only the corresponding physical eigenstates
            of the QHO.
            Besides, there exists also a  finite number
            of intermediate (γ)-type schemes
            which simultaneously 
            use the eigenstates of both types.
            Let us denote _n_+^- and _n_+^+
             the mutually conjugate intertwining operators 
            constructed on the basis of the (α)-type scheme with
            a minimal number  n_+ of  seed non-physical eigenstates. 
            They are
             differential operators of order n_+.
            Analogously, let us denote the 
            intertwining operators constructed on the basis
            of the (β)-scheme with a minimal number 2n_-
            of seed physical eigenstates as 
            _2n_-^- and _2n_-^+.
            By the construction we have 
            ^-_n_+ψ^-_j_s=0 and ^-_2n_-ψ_i_s=0,
            where ψ^-_j_s are n_+ non-physical eigenstates 
            of the QHO
            which are the seed states in the  (α_n_+)-scheme,
             ψ^-_j_s∈ (α_n_+), while
             ψ_i_s are 2n_-  physical eigenstates which are used 
             in the (β_2n_-)-scheme, 
             ψ_i_s∈ (β_2n_-). 
             Such two minimal (α_n_+)-
            and (β_2n_-)-schemes are complementary 
            similarly to the  schemes (α_1) and (β_2)
            in the case of the simplest REQHO system
            we considered in detail above.
             The corresponding 
            Wronskians in these two schemes  have the form 
            _n_+=exp(1/2x^2n_+)ϕ(x)
            and 
            _2n_-=cexp(-x^2 n_-)ϕ(x),
            where ϕ(x) is some nodeless 
            polynomial 
            function, 
            and  c is some constant.
            We fix the additive constant shifts 
            in  the Hamiltonians H of the QHO  
            and H̆ of the REQHO 
            in such a way that ^-_n_+H=H̆^-_n_+
            and  that the ground-state Ψ_0 of H̆ 
            has 
            zero  energy, E_0=0.  Then 
            the lowest state Ψ_n_+=^-_n_+ψ_0 in the equidistant  
            infinite part of the spectrum of H̆ will be characterized
            by the energy value 
            
            E_n_+=2(n_++2n_-)≡ 2Δ
            that also will be the energy of the ground-state ψ_0 of 
            the shifted QHO Hamiltonian H,
            Hψ_0=E_n_+ψ_0.
            Then we have H=H_osc-1+2Δ.
            For the other basic lowering 
            intertwining operator ^-_2n_-
            we have the relation
             ^-_2n_-H=(H̆+2Δ)^-_2n_- .
             
            This  means in particular that if
             ψ(x;E) is an eigenstate of H of energy E and 
             if ^-_2n_-ψ(x;E)≠ 0,
             then this latter state  will be
             eigenstate of H̆
             of the eigenvalue (E-2Δ),
             H̆ (^-_2n_-ψ(x;E))=
             (E-2Δ)^-_2n_-ψ(x;E).
             
             
            In terms of the operators _n_+^±
            and _2n_-^±
            we construct the three pairs  of the 
            basic ladder operators 
            𝒜^±=_n_+^- a^±_n+^+ ,
            ℬ^±=_2n_-^- a^±_2n-^+ ,
            
             and 
            𝒞^-=_2n_-^-^+_n_+ ,    𝒞^+=^-_n_+_2n_-^+ .
            
            The operators 𝒜^± and ℬ^± are
            differential operators of  
            orders 2n_++1 and 4n_-+1, respectively, while
            the ladder operators 𝒞^± 
            are  differential operators of  order 
            n_++2n_-.
            These basic ladder operators  satisfy the relations 
            
            [H̆,𝒜^±]=
            ± 2 𝒜^± ,     
            [H̆,ℬ^±]=
            ± 2 ℬ^± ,     
            [H̆,𝒞^±]=
            ±  2Δ 𝒞^± .
            
            We also have the operator identities 
            𝒜^+𝒜^-=𝒫_𝒜(H̆) ,    𝒜^-𝒜^+=𝒫_𝒜(H̆+2) ,
            ℬ^+ℬ^-=𝒫_ℬ(H̆) ,    ℬ^-ℬ^+=𝒫_ℬ(H̆+2) ,
            
            and
            𝒞^+𝒞^-=𝒫_𝒞(H̆) ,    𝒞^-𝒞^+=𝒫_𝒞(H̆+2Δ) ,
            
            where 
            𝒫_𝒜(H̆)=
            (H̆-2Δ)
            𝒫_(H̆-2)𝒫_(H̆) ,  𝒫_ℬ(H̆)=
            H̆𝒫_(H̆+2Δ-2)𝒫_(H̆+2Δ) ,
            𝒫_𝒞(H̆)=𝒫_(H̆)
            𝒫_(H̆) .
            
            The polynomial 𝒫_ of order n_+ is  defined here by
            𝒫_(H̆)≡^-_n_+^+_n_+=H̆∏_i=1^n_+-1(H̆-E_i),
            where E_i, i=1,… ,n_+-1 , are nonzero
            eigenvalues of the 
            corresponding excited separated (gapped) physical eigenstates 
            Ψ_i of H̆. 
            Together with zero energy E_0=0 of the ground state Ψ_0,
            the energy values of the n_+ gapped physical 
             eigenstates of H̆ are the shifted by the  constant 
            2Δ=2(2n_-+n_+) energies of the corresponding 
            non-physical eigenstates
            ψ^-_j_s  of H which appear as the seed
            states in the minimal  (α_n_+)-scheme.
            The permuted product of the intertwining 
            operators gives here
              ^+_n_+^-_n_+=𝒫_(H).
            The polynomial 
            𝒫_ of order 2n_- is  defined 
            via the relation
            𝒫_(H)≡^+_2n_-^-_2n_-=∏_j=1^2n_-(H-E_j^-),
            where  by E^-_j 
            we denote  the shifted for the same  constant 
            2Δ
            energies
            of the physical eigenstates ψ_i_s   of H
            which are present as the seed states  in the
            minimal (β_2n_-)-scheme.
            For the permuted product of the 
            intertwining operators 
            we have  ^-_2n_-^+_2n_-=
            𝒫_(H̆+2Δ).
            Let us also note  here a useful relation 
            𝒫_(H̆)
            𝒫_(H̆+2Δ)=
            ∏_j=0^2n_-+n_+-1(H̆-2j) ,
            
            that reflects  the complementarity of the minimal 
            (α_n_+)- and (β_2n_-)-schemes.
            
            The operator 𝒜^-  annihilates
            all the n_+ physical eigenstates 
            Ψ_0,…, Ψ_n_+-1
            of the  system H̆
            whose energies  lie below the infinite equidistant
            part of the spectrum and are separated from 
            it by some gap of 2n_0 missing energy levels,
            n_0≥ 1.
            Between these n_+ 
            separated energy levels  there can appear g,
            0≤ g<n_+,
            `internal' gaps each one containing   an
            even number of missing energy levels.
            We name the g+1 sets of energy levels in the lower
            separated part of the spectrum which do not contain 
            internal gaps as valence bands.
            If g>0, we denote by 2n_1,…, 2n_g
            the number of missing energy levels in the corresponding 
            internal energy gaps assuming that the highest
             value g of index i 
            in n_i corresponds here 
            to the lowest energy gap in the spectrum.
            The total number of the missing energy levels 
            2(n_0+…+n_g) is equal to the number
            2n_-
            of the physical eigenstates ψ_i_s which participate 
            as the seed states in the minimal (β_2n_-)-scheme.
            In addition to the gapped physical eigenstates 
            Ψ_0,…,
            Ψ_n_+-1, 
            the operator 𝒜^- also annihilates 
            the lowest state 
            Ψ_n_+ of energy E_n_+=2Δ
            in 
            the infinite equidistant part
            of the spectrum due to the presence of the operator
            a^- in its  structure.
            Besides, 𝒜^- annihilates 
            some n_+ non-physical eigenstates of H̆.
            Note in particular that  
            Ψ_0=^-_n_+ψ^-_j_+=
            ^-_2n_-ψ_0 ,    Ψ_n_+-1=^-ψ^-_2n_0=
            ^-_2n_-ψ_j_+-2n_0 ,
            Ψ_n_+=^-_n_+ψ_0=
            ^-_2n_-ψ_2n_-+n_+ .
            
            Here j_+=2n_-+n_+-1 is the maximal value
            of the index j_s of the non-physical eigenstates
            ψ^-_j_s∈  (^-_n_+)
            from the minimal (α_n_+)-scheme;
            it coincides with the maximal value 
            of the  index i_s
            of the physical eigenstates
            ψ_i_s∈ (^-_2n_-) from 
            the minimal scheme (β_2n_-).
            The energy values of the indicated 2n_++1 physical 
            and non-physical eigenstates 
            from the kernel of 𝒜^- are 
            the roots of the polynomial 𝒫_𝒜(H̆)
            which appears in the first identity in (<ref>).
            The kernel of 𝒜^+ is spanned by 
            the  n_+ lowest separated  physical eigenstates 
            Ψ_0,…, Ψ_n_+-1,
            and  by some n_++1 non-physical eigenstates
            of H̆. The energy values of these
            eigenstates from (𝒜^+)
            correspond to the roots of the polynomial 𝒫_𝒜(H̆+2)
            which appears in the second
            relation  in (<ref>).
            
            The kernel of the ladder operator 
            𝒞^- is 
            spanned by n_++2n_- physical states,
            n_+ of which, Ψ_0,…, Ψ_n_+-1,  
            correspond to the lowest separated 
            (gapped) energy values. The  other  2n_- eigenstates of H̆ 
            in the kernel of  𝒞^- are 
            the physical states ^-_n_+ψ_i_s
            in a lower part
            of the infinite equidistant part of the spectrum,
            where ψ_i_s∈ (^-_2n_-).
             The number of those `supplementary' states in the lower part
            of the equidistant spectrum which are not annihilated 
            by 𝒞^- and whose eigenvalues
            lie below  the highest energy value 
            of a physical eigenstate from (𝒞^-) 
            is equal to n_+.
            The number of the sequential  lowest states 
            at the very bottom of the  equidistant 
            infinite part of the spectrum of H̆ 
            which are not annihilated by 𝒞^- 
            is equal to the number of the 
            physical states 
            in the lowest valence band of the separated part of the spectrum. 
            The kernel of 𝒞^+ will be spanned
            by some n_++2n_- non-physical eigenstates.
            The energies of the corresponding eigenstates from the 
            kernels of the ladder operators 𝒞^- and 
            𝒞^+ correspond to the roots of the polynomials 
            in H̆ which appear in the first and the second identities in
            (<ref>). 
            
            The kernel of the ladder operator ℬ^- 
            contains g+1 physical eigenstates, each one lying
             at the very  bottom of each valence band. Besides, 
            it also contains the physical eigenstate Ψ_n_+
            of energy E_n_+=2Δ which is the lowest state
            of the equidistant infinite part of the spectrum.
            In addition, (ℬ^-) 
            contains 2n_- non-physical eigenstates 
            of the form Ψ^non_i_s≡^-ψ_i_s,
            which correspond to the missing energy values 
            in the gaps. Finally, it also involves 2n_- - (g+1) 
            non-physical eigenstates
            of the form Ψ^non_i_s
            for all values of the index i_s
            except those g+1 values,  each one of which 
            corresponds to a lowest eigenstate
            in each gap.  
            The kernel of the increasing 
            ladder operator ℬ^+ 
            contains g+1 physical eigenstates
            whose eigenvalues lie at the top of
            each valence band. It also involves  4n_- - g 
            non-physical eigenstates from the gaps.
            Eigenvalues of the eigenstates  from the kernels
            of the ladder operators ℬ^- and 
            ℬ^+ correspond to the roots of the 
            polynomials which appear, respectively,  in the first and 
            the second identities
            in (<ref>). 
            
            
             The relations of the form (<ref>) and 
            (<ref>)  are valid
            for the basic intertwining operators 
            of the arbitrary REQHO system
            we consider here.
            Indeed, the operator ^-_2n_- annihilates all the 2n_-
            physical eigenstates  ψ_i_s of the QHO which
            participate as the seed states in the (β_2n_-)-scheme.
            On the other hand, when ^-_2n_-  acts on the n_+  
            `supplementary'  eigenstates  in the lower part of the 
            spectrum of the QHO, it transforms these states
            into the separated  lowest n_+ 
            physical eigenstates 
            of the REQHO system
            H̆ which constitute the kernel 
            of the intertwining operator _n_+^+. 
            In particular, as we saw ^-_2n_- maps the ground-state ψ_0 
            of the shifted QHO, Hψ_0=2Δψ_0, 
            into the zero energy ground-state Ψ_0 of H̆,
            ^-_2n_-ψ_0=Ψ_0, 
            H̆Ψ_0=0.
            We conclude then that the composite operator 
            _n_+^+^-_2n_- annihilates all the 2n_-+n_+ 
            eigenstates ψ_n
            of the QHO with n=0,…,2n_-+n_+ -1.  
            But the same job is made by the lowering 
            ladder operator 
            (a^-)^2n_-+n_+ of the QHO.
            From here we obtain the operator equalities
             ^+_n_+^-_2n_-=(-1)^n_+(a^-)^2n_-+n_+ ,    ^+_2n_-^-_n_+=(-1)^n_+(a^+)^2n_-+n_+ .
             
             The  relations (<ref>) reflect
             the complementary nature of the 
             involved minimal (α_n_+)- and (β_2n_-)- schemes.
            The identities (<ref>) are employed
            to establish the operator identities in (<ref>).
            They are also essential for  the analysis of the
            kernels of the basic ladder operators.
            
            Similarly to the simplest case of the REQHO system,
            one can construct other pairs  of secondary   
            ladder operators  different from the described 
            basic  ladder operators  
            𝒜^±,  ℬ^± and 𝒞^±.
            This can be done effectively by introducing 
            additional factors  (a^±)^n
            inside the structure of these operators :
            𝒜^±_n≡^-_n_+(a^±)^n^+_n_+, 
            ℬ^±_n≡^-_2n_-(a^±)^n^+_2n_-,
            n=1,…,
            𝒜^±_1=𝒜^±,
            ℬ^±_1=ℬ^±,
            and
            𝒞^+_n+1≡^-_n_+(a^+)^n_2n_-^+,
            𝒞^-_n+1≡_2n_-^- (a^-)^n ^+_n_+, 
            where n=0,…,
            𝒞^±_1=𝒞^±.
             One can also consider
            the operators 𝒞^-_-n=^-_2n_-(a^+)^n^+_n_+,
            𝒞^+_-n=^-_n_+ (a^-)^n ^+_2n_- with
            n=1,…,2n_-+n_+-1. In 
            𝒞^±_-n we restrict the values of  
            the  index n from above
            having in mind the identity  
            𝒞^-_-(2n_-+n_+)=
            (-1)^n_+ 𝒫_(H̆) 𝒫_(H̆+2Δ) ,
            
            see Eq. (<ref>), 
            and so, for n>2n_-+n_+-1 these operators 
            do not provide  essentially  new structures. 
            The secondary, higher-order ladder operators 
            can also be obtained 
            by taking  the compositions
             of the intertwining operators of the 
            corresponding 
            (α)-, (β)- and  (γ)-schemes.
            They also are generated via the 
            composition  of the basic ladder operators
            𝒜^±, ℬ^± 
             and 𝒞^±.
            In particular,  the quadratic compositions 
            of 𝒜^±
             and 𝒞^±
            are given by (<ref>), (<ref>),
            and by the relations 
            
            (𝒜^+)^2=𝒫_(H̆-2) 𝒜^+_2 ,
              
            (𝒞^+)^2=(-1)^n_+𝒞^+_2n_-+n_++1 ,  𝒜^+𝒞^-=(-1)^n_+(H̆-2Δ)𝒜^-_2n_-+n_+-1 ,
            𝒜^+𝒞^+=𝒫_(H̆-2) 𝒞^+_1 ,  𝒜^-𝒞^-=(-1)^n_+𝒜^-_2n_-+n_++1 ,  𝒜^-𝒞^+=𝒫_(H̆+2)·𝒞^+_-1 ,
            
            and by the relations conjugate to (<ref>) and (<ref>).
            The relation 
            
            (𝒜^-)^2n_-+n_+=
            (-1)^n_+∏_l=0^2n_-+n_+-1𝒫_(H̆+2l)·𝒞^-
            
            shows that as in the considered particular cases 
            of the REQHO systems, 
            the ladder operators 𝒞^± can be
            generated  by the operators 𝒜^±.
            Also,  the following operator identity is valid: 
            
            (𝒜^-)^2n_-+n_+-1=
            (-1)^n_+1/H̆∏_j=0^2n_-+n_+-2𝒫_(H̆+2j) 𝒞^-_-1 .
             
            Here the operator multiplier before 𝒞^-_-1
            is  the polynomial of order n_+(2n_-+n_+-1)-1
            in H̆
            since the j=0 term 𝒫_(H̆) 
            in the product is equal to 
            the factor H̆ which cancels the multiplier 
            1/H̆ before 
            the product symbol.
            We also have the  identity 
            which relates the operators 𝒜^- and 
            ℬ^-,
            
            (H̆-2Δ+2)𝒫_(H̆+2) ℬ^-=
            (H̆+2)𝒫_(H̆+2Δ) 𝒜^- .
            
            The analogous identity for 𝒜^+ and 
            ℬ^+ is obtained from 
            (<ref>) by Hermitian conjugation.
            In particular cases of the three REQHO systems considered in the previous
            two sections, relation (<ref>)
            reduces to the first relation in (<ref>) and to
            the
            identities (<ref>) and 
            (<ref>). 
            
            In conclusion of this section, let us  show that the trinity 
            (𝒜^±,ℬ^±,𝒞^±)
            of the pairs of the lowering  and raising ladder operators 
            allows us to generate an arbitrary physical eigenstate from the 
            ground state Ψ_0, and as a consequence, any two 
            physical eigenstates can be related by the appropriate consecutive 
            action of the basic ladder operators from the trinity.
            First, from the described properties of the operators 
            and commutation relations (<ref>) 
            it follows that in the equidistant infinite part of the spectrum 
            any two eigenstates can be related by
            the ladder operators 𝒜^± and ℬ^± 
            in the same way as  the ladder operators a^± relate the states in the 
            QHO system. 
            The only difference will appear in the numerical 
            coefficients which have to be included into the composition of the
            indicated basic operators when we work with the normalized
            eigenstates. 
            If a valence band contains more than one eigenstate,
            different states in this band can be connected by application to them of the
            appropriate degrees of the lowering and raising operators   ℬ^- and 
            ℬ^+. Note that within the valence band  with n_i states 
            these operators satisfy the identity (ℬ^±)^n_i=0.
            Recall also that the lowest state Ψ_n_+ in the equidistant 
            infinite part of the spectrum is related with the ground state 
            by the action of the ladder operators 𝒞^±:
             Ψ_0=𝒞^-Ψ_n_+ and 
            Ψ_n_+=𝒞^+Ψ_0.
            In the same way 
            one can relate any  state Ψ_n 
            of energy 0<E_n<2Δ with 0<n≤n_+-1 from 
            the separated part of the spectrum with the corresponding
            state Ψ_n_++n of energy E_n+2Δ from the equidistant 
            infinite part of the spectrum. Then, if a REQHO system contains more than 
            one valence band, the ground state Ψ_0 from the 
            lowest valence band can be related to some state Ψ_l with eigenvalue 
            E_l from some higher valence band, for instance,
            by the following composition of the ladder operators:
            Ψ_l=𝒞^-(𝒜^+)^r_l𝒞^+Ψ_0,
            Ψ_0=𝒞^- (𝒜^-)^r_l𝒞^+Ψ_l, 
            where r_l=E_l/2. This shows finally  
            that the ladder operators from the  trinity are the 
            spectrum-generating operators for  the REQHO system of a 
            general form. 
            
            The described properties of the REQHO systems 
            of a general form are illustrated by  Figure <ref>.
            1cm
            
            
            
            § SUMMARY AND OUTLOOK
            
            
            In conclusion, we summarize the obtained results
            and indicate some interesting problems  for 
            further investigation.
            
            
            The  REQHO  system of a general form is characterized
            by n_+≥ 1   low-lying  energy levels 
            which are separated from 
            the higher  equidistant infinite part of the spectrum by some 
            gap  of  an even number 2n_0≥ 2 of missing  levels. 
            Between  separated energy levels there can be 
            additional gaps of an  even number of missing levels 
            in each such a gap. As a result,
            in the lower  part
            of the spectrum there appears a total number 2n_-≥ 2n_0
            of missing levels, and  the separated part of the  spectrum is organized  
            into (g+1)≥ 1  `valence bands'.
            Such a peculiar structure of the spectrum of the 
            REQHO system characterized by the three integer numbers
            (n_+,2n_-,g+1) is detected and reflected 
            by the trinity (𝒜^±,ℬ^±,𝒞^±)
            of the pairs of the lowering and raising ladder operators.
            Any raising or lowering  ladder operator from the trinity when acts on a 
            physical eigenstate either transforms it into another physical 
            eigenstate with the changed value of energy, or annihilates it.
            
            The separated states 
            are detected by 
            𝒜^- and 𝒜^+ that are differential
            operators of order 2n_++1.  Each one of these two
            operators annihilates all the n_+ states
            Ψ_0,…, Ψ_n_+-1
             in the valence bands.
            The operator 𝒜^-
            annihilates in addition the lowest state 
            Ψ_n_+
            in the equidistant infinite part of the spectrum.
            In this aspect the lowering ladder operator 𝒜^-
            has  properties very similar to those
             of the Lax-Novikov integral in reflectionless
            systems we discussed in Section <ref>. 
            Besides, the operators 𝒜^-
            and 𝒜^+ annihilate some 
            n_+ and n_++1 non-physical eigenstates, respectively.
            Due to the relation
            [H̆,𝒜^±]=± 2𝒜^±,
            the 𝒜^±
            act in the equidistant infinite  part of the spectrum
            as the spectrum-generating operators 
            like the ladder operators a^± in the 
            QHO system. Namely, 
            they transform  physical eigenstates 
            Ψ_j with j≥ n_+
            into the states Ψ_j± 1 
            by shifting the energy values in 
            ± 2, and 𝒜^-Ψ_n_+=0.
            It is also worth to note here that 
            the quadratic in the
            ladder operators 𝒜^±
            relations in (<ref>) 
            are analogous to the 
            Burchnall-Chaundy polynomial 
            identity
            <cit.>
            that relates the Lax-Novikov  integral 
            with the corresponding Hamiltonian 
            of a reflectionless (or a finite-gap) system,
            and underlies the modern theory of 
            integrable systems <cit.>. 
              
            
             The ladder operators ℬ^±
             are differential operators of order 4n_-+1,
             and each of them  effectively counts
             the number g+1  of the valence bands
             and measures
              the size of each valence band.
             This is done as follows. The lowering operator
             ℬ^-  annihilates each physical eigenstate
             which lies at the very bottom of each valence band,
             whereas the raising operator  ℬ^+ 
            annihilates each state  at the very top of
            each valence band.
            So, if a valence band contains only one state,
            this state 
            is annihilated by both ℬ^- and ℬ^+,
            and in such a one-dimensional 
            valence band the action of ℬ^±
            is similar to that of 𝒜^±.
            However, if we have a valence band with more than one 
            state, 
            the operators ℬ^+ and ℬ^-,
            unlike 𝒜^±,
            act in such a band as the raising 
            and lowering operators,
            [H̆,ℬ^±]=± 2ℬ^±,
            which satisfy there the 
            relations (ℬ^±)^n_i=0,
            where n_i is the number of 
            states in the band.
            Like 𝒜^-, the  operator ℬ^- 
            annihilates the lowest state Ψ_n_+ 
            in the equidistant infinite
            part of the spectrum.
            The kernels of 
            ℬ^-  and ℬ^+ 
            also include some 4n_–g-1
            and  4n_–g non-physical eigenstates, respectively.
            Similarly to 𝒜^±, 
            in the equidistant infinite part of the spectrum
            the operators ℬ^± also act 
            as the spectrum-generating operators.
            
            Under the action of the ladder operators 𝒜^± 
            and ℬ^± the n_+ states 
            from the valence bands turn out to be completely disconnected 
            from the physical states
            in the equidistant infinite part of the spectrum.
            Both parts of the spectrum are connected
            by means of the third pair of the 
            mutually conjugate ladder operators 𝒞^+ and 
            𝒞^-, that are 
             differential 
            operators of order n_++2n_-. 
            The kernel of the lowering operator 
            𝒞^- is spanned by physical
            eigenstates, n_+ of which 
            correspond to  all the n_+ eigenstates 
            from the valence bands.
            The rest 2n_- states from  
            (𝒞^-) 
            are some eigenstates 
            with energy levels lying in the low part
            of the  equidistant  infinite part of the spectrum.
            The positions of those 2n_- energy levels 
            correspond to the missing energy levels 
            in the gaps moved up for the distance 
            2Δ=2n_++4n_-≥ 6 which is 
            exactly equal   to  the distance  between  the energy levels of the 
            ground state Ψ_0 and the lowest state Ψ_n_+ 
            in the equidistant infinite part  of the spectrum.
            The kernel of 𝒞^+ is spanned by
            some non-physical eigenstates only.
            Due to the relation 
            [H̆,𝒞^±]=± 2Δ 𝒞^±,
            the operators 𝒞^+ and 𝒞^-
            act as the 
            raising and lowering operators
            changing the energy for ± 2Δ.
            All the states from the valence bands
            are obtained by application of the lowering 
            operator 𝒞^- to those low-lying  states with energies 
            2Δ≤ E<4Δ in the 
            equidistant infinite part of the spectrum 
            which are not annihilated by it.
            In particular, the lowest state 
            Ψ_n_+ with energy E_n_+=2Δ 
            in the equidistant infinite part of the spectrum is transformed by 𝒞^-
            into the ground-state Ψ_0 of zero energy.
            The state Ψ_n_+, in turn,  
            can be obtained from Ψ_0 
            by action of the raising operator 𝒞^+, and 
            is also generated 
            from the state Ψ_n_++Δ
            of energy E_n_++Δ=4Δ by applying to the 
            latter the lowering operator   𝒞^-.
            As a consequence of the described properties, 
            any two states in the spectrum of  a REQHO system can be
            related by an appropriate consecutive 
            action of the basic ladder operators from the trinity.
            In particular, arbitrary excited state from any valence band
            or from the equidistant infinite part of the spectrum 
            can be obtained from the ground state Ψ_0.
            This means that the basic ladder operators from the trinity 
            are the spectrum-generating operators of the REQHO system.
            
            The energies of all the physical and non-physical eigenstates
            of the kernels of the lowering operators 
            𝒜^-, ℬ^- and 𝒞^-
            are   the roots of the corresponding polynomials in H̆  
            which appear on the right hand side in  the first relations from equations 
            (<ref>), (<ref>) and (<ref>), respectively.
             The eigenvalues of the physical and non-physical eigenstates 
             from the kernels of the conjugate operators
             𝒜^+, ℬ^+ and 𝒞^+
            are the roots of the corresponding polynomials 
            which appear in the second operator identity relations
            in the same equations. 
            The   basic ladder operators 𝒜^- and ℬ^- 
            satisfy the  two-term  identity relation (<ref>)  
            which is  linear in 
            both of these operators but involve the 
            coefficients that 
            are certain polynomials in the Hamiltonian H̆.
            The presence of such polynomial coefficients reflects 
            a difference in action of these operators on the states
            in a separated part of the spectrum. 
            The operators 𝒜^- and 𝒞^-
            are related by the operator identity of the form
            (<ref>). 
            Proceeding from these relations,
            one can obtain the identity  that  
            relates the operators ℬ^- 
            and 𝒞^-, and by conjugation 
            one can find the identities that relate 
            the raising operators of the trinity.
            
            The operators 𝒜^± 
            are constructed as the ladder operators a^± of the 
            QHO dressed by means of the Darboux-Crum-Krein-Adler  intertwining 
            operators ^-_n_+ and ^+_n_+ constructed on the basis of 
            the minimal set of n_+ non-physical eigenstates
            of the QHO which are used as the seed states
            in the corresponding DCKAT based on the (α_n_+)-scheme,
            see Eq. (<ref>). 
            The operators ℬ^±
            are constructed in the same way with the help of the intertwining operators
            ^-_2n_- and ^+_2n_- obtained 
            on the basis 
            of the minimal set of the 2n_- physical eigenstates 
            which are employed as the seed states in the DCKAT 
            in the corresponding (β_2n_-)-scheme,
            see Eq. (<ref>). 
            The operators 𝒞^± can be obtained 
            as the composition (<ref>) of the corresponding 
            intertwining operators from both indicated 
            schemes of the DCKATs. 
            The minimal schemes  (α_n_+) and 
            (β_2n_-) are complementary, 
            what  is reflected in particular by the relations 
            (<ref>) and
            (<ref>).
            The  secondary, higher-order ladder operators can be constructed 
            in analogous way by dressing the higher-order 
            ladder operators (a^±)^n of the QHO, or
            by the composition of the appropriate 
            intertwining operators from the non-minimal (α_n_++n)
            and (β_2n_-+n) schemes,
            or by employing the intertwining operators from the corresponding 
            intermediate (γ)-type schemes 
            which use  both physical and non-physical eigenstates
            of the QHO as the seeds states of the corresponding 
            DCKATs. The  secondary ladder operators can  also be generated
            via the appropriate composition of the 
            basic (primary)  ladder operators 𝒜^±,
             ℬ^± and 𝒞^±.
             
             It seems to be interesting to investigate the 
             quantum mechanical systems related to the
             exceptional Laguerre and Jacobi orthogonal polynomials
             in the light of the results on the ladder operators
             obtained here. 
             The  results of such an investigation will be presented elsewhere.
             
             
             The ladder operators ℬ^± 
             have a nature of the polynomially deformed 
              bosonic creation and annihilation
              operators in the equidistant infinite part 
              of the spectrum.  On the other hand, 
              these operators act trivially on
             the one-state valence bands 
            whose singleton states  are annihilated by both 
            the lowering ℬ^-
            and the raising ℬ^+ 
            operators. The same operators
            reveal the properties of 
            the deformed fermionic
            creation and annihilation operators
            in the 
            valence bands consisting from 
            two states. 
            They have the properties of the deformed para-fermion
            creation and annihilation operators of order n>2 in those 
            valence bands which contain   n>2
            eigenstates of H̆.   
            The interesting question is then if there exist some 
            concrete physical systems which would 
            reveal  the spectrum of the 
            REQHO systems. 
            If so, it seems that the trinity of the ladder operators 
            should play a fundamental  role in the physics 
            associated  with such  systems. 
            In the same direction the interesting question is whether
            the quantum mechanical REQHO systems and the    
            structures associated with them 
            can be generalized somehow for the case of  
            the quantum fields.  
            
            
            In <cit.>, the family of  the REQHO 
            systems with two separated states 
            generated by  non-physical seed states
            ψ^-_m_1 and ψ^-_m_2, 
            m_2-m_1≡ℓ=1+2r, r=0,1,…,
            m_1=2k, k=1,…,  
            was considered.
            For such class of the systems,
            there the  lowering, c,  and 
            increasing, c^†, ladder operators of the differential order
            2+ℓ were constructed
            by employing auxiliary systems  some of which 
            are singular and have a nature similar
            to that of the isotonic oscillator (<ref>).
            In the simplest case m_1=2 and m_2=3
            such a system corresponds to the REQHO system (<ref>)
            we considered in Section <ref>.
            Like our fifth order ladder operator ℬ^-,
            the third order ladder operator c from <cit.> 
            annihilates the
            ground state Ψ_0 and the lowest state Ψ_2 
            in the infinite equidistant part of the spectrum.
            The increasing operator  c^† like 
            our ℬ^+ annihilates the excited state Ψ_1
            in the separated two-state lower part of the spectrum.
            In the systems with ℓ>1,  however, the kernel of the increasing operator
            c^† still includes only one physical state which is, again, 
            the separated state Ψ_1, while our ℬ^+
            operator annihilates both separated states Ψ_0 and Ψ_1. 
            In addition to the separated ground 
            state Ψ_0 and the lowest state Ψ_2 
            in the equidistant part of the spectrum,
            the kernel of  the lowering operator c in this case includes
            also the ℓ-1 excited states 
            Ψ_3,…,Ψ_2+ℓ in the equidistant part of the spectrum.
            In this aspect, the lowering operator c from  <cit.> 
            has some similarity 
            with our operator 𝒞^-. But our  
            ladder operator 𝒞^- is
            of differential order 2m_1+2r and its kernel includes
            some 2m_1+2r-2 excited states in the equidistant 
            part of the spectrum together with
            both separated states 
            Ψ_0 and Ψ_1. 
            Thus, in the case of ℓ>1 the nature of the 
            operators c and c^†  in the sense  
            of the physical states which  they annihilate  is different from
            the nature of any of
            our lowering and increasing ladder operators 
            𝒜^-, 𝒜^+,
            ℬ^-, ℬ^+ and 𝒞^-, 𝒞^+.
            It would be interesting to investigate whether the analogs 
            of the ladder operators c and c^† from  <cit.>  
            can be 
            constructed for REQHO systems containing more than two 
             states in the lower separated part of the spectrum,
             and what is  the exact  relation of such ladder operators with
            our trinity (𝒜^±,
            ℬ^±, 𝒞^±)  of the ladder operators.
            
            
            
            0.2cm
             
              
              
             
             Acknowledgements0.4cm
            
            JFC and MSP acknowledge support from
            research projects FONDECYT 1130017 (Chile),
            Proyecto Basal USA1555 (Chile),
            MTM2015-64166-C2-1 (MINECO, Madrid) and DGA E24/1 (DGA, Zaragoza).
            MSP is grateful  for the warm hospitality at  Zaragoza University.
            JFC thanks for the kind hospitality at Universidad de Santiago de Chile.
            
            
            
            2cm
            
            
            
            
            
            
            99
            KaiMos
            I. Kay and H. E. Moses, 
             “Reflectionless transmission through dielectrics and scattering potentials,",
             http://aip.scitation.org/doi/abs/10.1063/1.1722296J. Appl. Phys.  27, 1503 (1956).
             Darb
            G. Darboux,“Sur une proposition relative aux équations 
            linéaires,"  C. R. Acad. Sci Paris  94,  1456 (1882).
            
            
            Schr
              E. Schrödinger,
              “A method of determining quantum-mechanical eigenvalues and eigenfunctions,”
              Proc. Roy. Irish Acad. (Sect. A)  46, 9  (1940).
              
              InfHull 
              L. Infeld and T. E. Hull,
            “The factorization method,”http://journals.aps.org/rmp/abstract/10.1103/RevModPhys.23.21
              Rev. Mod. Phys.  23, 21 (1951).
             
              
              Crum
            M. M. Crum, “Associated Sturm-Liouville systems,"http://qjmath.oxfordjournals.org/content/6/1/121.full.pdf+html
            Quart. J. Math. Oxford   6, 121 (1955).
            
             MatSal
            V. B. Matveev and M. A. Salle, Darboux Transformations and Solitons
            (Springer, Berlin,
            1991).
            
            CR00 J. F. Cariñena and  A. Ramos, 
             “Riccati equation, factorization method and shape invariance,"http://www.worldscientific.com/doi/abs/10.1142/S0129055X00000502
            Rev.  Math. Phys.  12, 1279  (2000)http://arxiv.org/abs/math-ph/9910020[arXiv:math-ph/9910020].
            
               CR01 J. F. Cariñena, D. J. Fernández  and  A. Ramos, 
              “Group theoretical approach to the intertwined  Hamiltonians,”http://www.sciencedirect.com/science/article/pii/S0003491601961792
            Ann.  Phys.  292, 42  (2001)http://arxiv.org/abs/math-ph/0311029[arXiv:math-ph/0311029].
             MelRos
            B. Mielnik and O. Rosas-Ortiz,
            “Factorization: little or great algorithm?,”http://iopscience.iop.org/article/10.1088/0305-4470/37/43/001/meta
            J. Phys. A   37, 10007 (2004). 
              
              
            CR08 J. F. Cariñena and  A. Ramos, 
              “Generalized Bäcklund-Darboux transformations in one-dimensional 
            quantum mechanics,”http://www.worldscientific.com/doi/abs/10.1142/S0219887808002989
            Int. J. Geom. Methods Mod.  Phys.  5, 605 (2008). 
              
            GroNev
              D. J. Gross and A. Neveu,
              “Dynamical symmetry breaking in asymptotically free field theories,”http://journals.aps.org/prd/abstract/10.1103/PhysRevD.10.3235Phys. Rev. D  10, 3235 (1974).
             
              DasHasNev
            R. F.  Dashen, B. Hasslacher and A. Neveu,
              “Semiclassical bound states in an asymptotically free theory,"http://journals.aps.org/prd/abstract/10.1103/PhysRevD.12.2443Phys. Rev. D  12, 2443 (1975).
            
              CamBis
            D. K.  Campbell and  A. R. Bishop,
              “Soliton excitations in polyacetylene and relativistic field theory models,"http://www.sciencedirect.com/science/article/pii/055032138290089XNucl. Phys. B  200, 297 (1982).
            
              CoJaPl 
             F. Correa, V. Jakubsky and M. S. Plyushchay,
               “Aharonov-Bohm effect on AdS(2) and nonlinear supersymmetry of reflectionless Pöschl-Teller system,”http://www.sciencedirect.com/science/article/pii/S000349160900044X Annals Phys.  324, 1078 (2009)
             https://arxiv.org/abs/0809.2854 [arXiv:0809.2854 [hep-th]].
              LaxNov 
              A. Arancibia, J. Mateos Guilarte and M. S. Plyushchay,
             “Effect of scalings and translations on the supersymmetric 
             quantum mechanical structure of soliton systems,”http://journals.aps.org/prd/abstract/10.1103/PhysRevD.87.045009
              Phys. Rev. D  87,  045009 (2013) 
             http://arxiv.org/abs/1210.3666  [arXiv:1210.3666 [math-ph]].
              BBEIM
            E. D. Belokolos, A. I. Bobenko, V. Z. Enolskii, A. R. Its, and V. B. Matveev, 
              Algebro-Geometric Approach to Nonlinear Integrable Equations,
            (Springer, Berlin, 1994).
            
            
            GesHol
             F. Gesztesy and H. Holden,
               Soliton Equations and their Algebro-Geometric Solutions, (Cambridge Univ. Press, 2003). 
               
             exosusy 
             F. Correa, V. Jakubsky, L. M. Nieto and M. S. Plyushchay,
              “Self-isospectrality, special supersymmetry, and their effect 
              on the band structure,”http://journals.aps.org/prl/abstract/10.1103/PhysRevLett.101.030403Phys. Rev. Lett.  101, 030403 (2008) 
              http://arxiv.org/abs/0801.1671[arXiv:0801.1671 [hep-th]].
              
               
            
            Transmut 
              A. Arancibia and M. S. Plyushchay,
               “Transmutations of supersymmetry through soliton scattering, and self-consistent condensates,”http://journals.aps.org/prd/abstract/10.1103/PhysRevD.90.025008Phys. Rev. D  90, no. 2, 025008 (2014)
               http://arxiv.org/abs/1401.6709 [arXiv:1401.6709 [hep-th]].
              Dub
            S. Yu. Dubov, V. M. Eleonskii and N. E. Kulagin, 
             “Equidistant spectra of anharmonic oscillators,"http://www.jetp.ac.ru/cgi-bin/e/index/r/102/3/p814?a=list
            Zh. Eksp. Teor. Fiz.  102, 814 (1992);
            
            http://scitation.aip.org/content/aip/journal/chaos/4/1/10.1063/1.166056
            Chaos  4, 47 (1994).
            
             Adler
             V. E. Adler, “A modification of Crum's method,"http://link.springer.com/article/10.1007
             Theor. Math. Phys.  101, 1381 (1994).
            
            SamOvc
            B. F. Samsonov and I. N. Ovcharov, 
            “Darboux transformation and exactly solvable potentials with quasi-equidistant spectrum,"
            https://link.springer.com/article/10.1007/BF00559274Russ. Phys. J.  38, 765 (1995).
            
            BagSam
            V.G.Bagrov and B.F.Samsonov, 
            “Darboux transformation, factorization, and supersymmetry 
            in one-dimensional quantum mechanics,"https://link.springer.com/article/10.1007
            Theor. Math. Phys.   104,  1051 (1995).Spir
            V. Spiridonov,
             “Universal superpositions of coherent states and self-similar potentials,"http://journals.aps.org.ezproxy.usach.cl/pra/abstract/10.1103/PhysRevA.52.1909 Phys. 
            Rev. A  52, 1909 (1995)https://arxiv.org/abs/quant-ph/9601030[arXiv:quant-ph/9601030].
            
             JunRoy
            G. Junker  and  P. Roy, 
            “Conditionally exactly solvable problems and non-linear algebras,"
            http://www.sciencedirect.com/science/article/pii/S0375960197004222Phys. 
            Lett. A  232, 155 (1997).CPRS 
            J. F. Cariñena, A. M. Perelomov, M. F. Rañada and M. Santander, 
             “A quantum exactly solvable nonlinear oscillator related
            to the isotonic oscillator,"http://iopscience.iop.org/article/10.1088/1751-8113/41/8/085301/meta
             J. Phys. A: Math. Theor.  41, 085301 (2008) http://arxiv.org/abs/0711.4899  [arXiv:0711.4899 [quant-ph]].
            
            
            
            
            FellSmi
            J. M. Fellows and R. A. Smith, 
             “Factorization solution of a family of quantum nonlinear oscillators,"http://iopscience.iop.org/article/10.1088/1751-8113/42/33/335303/meta
              J. Phys. A  42, 335303 (2009) .
              
              
            Sesma
            J. Sesma, 
             “The generalized quantum isotonic oscillator,"http://iopscience.iop.org/article/10.1088/1751-8113/43/18/185303/meta
            J. Phys.  A
             43, 185303 (2010)  
             http://arxiv.org/abs/1005.1227[arXiv:1005.1227 [quant-ph]].
              
             
             GGM
            D. Gómez-Ullate, Y. Grandati, and R. Milson, 
             “Rational extensions of
            the quantum harmonic oscillator and exceptional Hermite polynomials,"http://iopscience.iop.org/article/10.1088/1751-8113/47/1/015203/metaJ.
             Phys. A  47, 015203 (2014) http://arxiv.org/abs/1306.5143[arXiv:1306.5143 [maph-ph]]. 
            
             
             Pupas
             A. M. Pupasov-Maksimov,
             “Propagators of isochronous an-harmonic oscillators 
             and Mehler formula for the exceptional Hermite polynomials,"http://www.sciencedirect.com.ezproxy.usach.cl/science/article/pii/S0003491615003565
            Annals of Physics  363, 122 (2015)https://arxiv.org/abs/1502.01778[arXiv:1502.01778 [maph-ph]].  
              
            Krein
            M. G. Krein,
            “On a continuous analogue of a Christoffel 
             formula from the theory of orthogonal polynomials," 
             Dokl. Akad. Nauk SSSR  113, 970 (1957).
             
             GomGrMil
             D. Gomez-Ullate, Y. Grandati, and R. Milson,
             “Extended Krein-Adler theorem for the translationally shape 
             invariant potentials,"http://aip.scitation.org/doi/10.1063/1.4871443J. Math. Phys.  55, 040201 (2014)http://arxiv.org/abs/1309.3756[arXiv:1309.3756 [nlin.SI]].
            
              
              Exc1
              D. Gómez-Ullate, N. Kamran, and  R. Milson,
             “An extended class of orthogonal polynomials defined by a Sturm-Liouville problem,"http://www.sciencedirect.com/science/article/pii/S0022247X09004569
             J. Math. Anal. Appl.  359
            (2009) 352https://arxiv.org/abs/0807.3939[arXiv:0807.3939 [math-ph]].
            
              Exc1+
              D. Gómez-Ullate, N. Kamran, and  R. Milson,
             “An extension of Bochners problem: exceptional invariant subspaces,"http://www.sciencedirect.com/science/article/pii/S0021904509001853
            J. Approx. Theory  162, 897 (2010)https://arxiv.org/abs/0805.3376[arXiv:0805.3376 [math-ph]].
            
            OdSas
             S. Odake and  R. Sasaki, 
             “Infinitely many shape invariant potentials and new 
             orthogonal polynomials, http://www.sciencedirect.com/science/article/pii/S0370269309009186
             Phys. Lett. B  679, 414 (2009)https://arxiv.org/abs/0906.0142[arXiv:0906.0142 [math-ph]].
              
             
             OdSas+
             S. Odake and R. Sasaki, 
             “Another set of infinitely many exceptional (X) 
             Laguerre polynomials, http://www.sciencedirect.com/science/article/pii/S0370269310000158
            Phys. Lett. B  684, 173 (2010)https://arxiv.org/abs/0911.3442[arXiv:0911.3442 [math-ph]].
            Ques+
            C. Quesne, 
             “Exceptional orthogonal polynomials, exactly 
             solvable potentials and supersymmetry,"http://iopscience.iop.org/article/10.1088/1751-8113/41/39/392001/meta
            J. Phys. A  41, 392001 (2008) http://arxiv.org/abs/0807.4087[arXiv:0807.4087 [quant-ph]].
            
            Ques++
            C. Quesne, 
             “Solvable rational potentials and exceptional 
            orthogonal polynomials in supersymmetric quantum mechanics,"http://www.emis.de/journals/SIGMA/2009/084/
            SIGMA  5, 084  (2009)https://arxiv.org/abs/0906.2331[arXiv:0906.2331 [math-ph]].
            
            
            SasTsZh
            R. Sasaki, S. Tsujimoto and  A. Zhedanov, 
             “Exceptional Laguerre and Jacobi polynomials 
            and the corresponding potentials through Darboux-Crum
             transformations,"http://iopscience.iop.org/article/10.1088/1751-8113/43/31/315204/meta
            J. Phys. A  43, 315204 (2010) https://arxiv.org/abs/1004.4711[arXiv:1004.4711 [math-ph]].
             Grand
            Y. Grandati, 
             “Solvable rational extensions of the isotonic oscillator,"http://www.sciencedirect.com/science/article/pii/S000349161100039X
            Ann. Phys.   326, 2074 (2011)  https://arxiv.org/abs/1101.0055[arXiv:1101.0055 [math-ph]].
            
            Sasaki 
              R. Sasaki,
              
              http://www.apcospa.org/vol2/2-2.htmlThe Universe  2, no. 2, 2 (2014)https://arxiv.org/abs/1411.2703[arXiv:1411.2703 [math-ph]].
              VesSha
            A. P. Veselov  and A. B. Shabat,   
            “Dressing chains and the spectral theory of the Schrödinger operator,"http://link.springer.com/article/10.1007
            Funct. Anal. Appl.  27, 81 (1993). 
            
              
            AdlPai
             V. E. Adler,  “Nonlinear chains and Painlevé equations,"http://www.sciencedirect.com/science/article/pii/016727899490104X
             Physica D  73, 335 (1994).
              
            Hiet
            R. Willox and J. Hietarinta, “Painlevé equations from Darboux chains: I.  PIII - PV,"http://iopscience.iop.org/article/10.1088/0305-4470/36/42/014
            J. Phys. A  36, 10615 (2003).  
            
            oblom
            A. Oblomkov, 
            “Monodromy-free Schrödinger operators with quadratically
             increasing potentials,"http://link.springer.com/article/10.1007/BF02557204Theor. Math. Phys.  121, 1574 (1999).
             
            
            MarQue1+
            I. Marquette and  C. Quesne,
              “Two-step rational extensions of the harmonic oscillator: exceptional orthogonal 
              polynomials and ladder operators," 
              http://iopscience.iop.org/article/10.1088/1751-8113/46/15/155201/metaJ. 
              Phys. A: Math. Theor.  46,  155201 (2013)https://arxiv.org/abs/1212.3474 [arXiv:1212.3474[math-ph]].MarQue1
              I. Marquette and  C. Quesne,
             “New ladder operators for a rational extension of the 
              harmonic oscillator and superintegrability of some two-dimensional systems,"http://aip.scitation.org/doi/10.1063/1.4823771J. Math. Phys.  54, 102102 (2013)https://arxiv.org/abs/1303.7150[arXiv:1303.7150 [math-ph]].
            
            
            MarQue1++
             I. Marquette and C. Quesne,
              “Combined state-adding and state-deleting approaches
               to type III multi-step rationally-extended potentials: applications to ladder operators and 
               superintegrability,”http://aip.scitation.org/doi/10.1063/1.4901006J. Math. Phys.  55, 112103 (2014)
             https://arxiv.org/abs/1402.6380 [arXiv:1402.6380 [math-ph]].
              MarI
            I. Marquette,
              “New families of superintegrable systems from k-step rational 
            extensions, polynomial algebras and degeneracies,"http://iopscience.iop.org/article/10.1088/1742-6596/597/1/012057/metaJ. Phys.: 
            Conf. Ser.  597, 012057 (2015)https://arxiv.org/abs/1412.0312 [arXiv:1412.0312 [math-ph]].
            
            Ques
            C. Quesne,
             “Ladder operators for solvable potentials connected with 
             exceptional orthogonal polynomials,"http://iopscience.iop.org/article/10.1088/1742-6596/597/1/012064/metaJ. Phys. 
             Conf. Ser.  597, 012064 (2015)https://arxiv.org/abs/1412.5874 [arXiv:1412.5874 [math-ph]].
            
            CarPly 
              J. F. Cariñena and M. S. Plyushchay,
               “Ground-state isolation and discrete flows in a rationally extended quantum harmonic oscillator,"http://journals.aps.org/prd/abstract/10.1103/PhysRevD.94.105022Phys. Rev. D 
             94, no. 10, 105022 (2016)
              https://arxiv.org/abs/1611.08051[arXiv:1611.08051 [hep-th]].
              SchwAn 
              M. S. Plyushchay,
              “Schwarzian derivative treatment of the quantum second-order supersymmetry anomaly, and coupling-constant metamorphosis,”http://www.sciencedirect.com/science/article/pii/S0003491616302743Ann. Phys., in press
             https://arxiv.org/abs/1602.02179[arXiv:1602.02179 [hep-th]].
              BurCha
            J.L. Burchnall, T.W. Chaundy,
             Commutative ordinary differential operators,
            http://onlinelibrary.wiley.com/doi/10.1112/plms/s2-21.1.420/abstract;jsessionid=0DB9784C385AC0EDE2AAD1190695C1AF.f02t01Proc. London Math. Soc. Ser. 2,  21, 420 (1923);
            
            
             Commutative ordinary differential operators,
            http://rspa.royalsocietypublishing.org/content/118/780/557Proc.  Royal Soc. London A  118, 557  (1928).
            
            Krich
            I. M. Krichever, 
              “Integration of nonlinear equations by the methods of 
              algebraic geometry,http://link.springer.com/article/10.1007
            
              “Baker-Akhiezer functions and integrable systems," 
              in:   “Integrability. The Seiberg-Witten and Whitham Equations", Edited by
              H. W. Braden and I. M. Krichever, p. 1 (Gordon and Breach Science Publishers, Amsterdam, 
              2000). ]
