São Paulo Journal of Mathematical Sciences https://doi.org/10.1007/s40863-020-00186-0 SPECIAL SECTION: GEOMETRY IN ALGEBRA AND ALGEBRA IN GEOMETRY Liouvillian solutions for second order linear differential equations with polynomial coefficients Primitivo B. Acosta‑Humánez1,2 · David Blázquez‑Sanz3  · Henock Venegas‑Gómez3 Accepted: 10 September 2020 © Instituto de Matemática e Estatística da Universidade de São Paulo 2020 Abstract In this paper we present an algebraic study concerning the general second order linear differential equation with polynomial coefficients. By means of Kovacic’s algorithm and asymptotic iteration method we find a degree independent algebraic description of the spectral set: the subset, in the parameter space, of Liouville inte- grable differential equations. For each fixed degree, we prove that the spectral set is a countable union of non accumulating algebraic varieties. This algebraic descrip- tion of the spectral set allow us to bound the number of eigenvalues for algebraically quasi-solvable potentials in the Schrödinger equation. Keywords Anharmonic oscillators · Asymptotic iteration method · Kovacic algorithm · Liouvillian solutions · Parameter space · Quasi-solvable model · Schrödinger equation · Spectral varieties Mathematics Subject Classification Primary 34M15 · Secondary 81Q35 1 Introduction Let us consider the family of second order linear differential equations, Communicated by Pedro Hernandez Rizzo. * David Blázquez-Sanz dblazquezs@unal.edu.co 1 Instituto Superior de Formación Docente Salomé Ureña, Santiago, Dominican Republic 2 Fac. Ciencias Básicas y Biomédicas, Universidad Simón Bolívar, Barranquilla, Colombia 3 Escuela de Matemáticas, Facultad de Ciencias, Universidad Nacional de Colombia - Sede Medellín, Medellín, Colombia Vol.:(012134 56789) São Paulo Journal of Mathematical Sciences u�� + P(x)u� + Q(x)u = 0, (1) with polynomial coefficients of bounded degree. This family is parameterized by the coefficients of P and Q and therefore endowed of an structure of affine algebraic variety. We are interested in characterizing the moduli of Liouville integrable dif- ferentential equations in (1) and describing how the Liouvillian solutions of those integrable equations depend on the coefficients. From a result of Singer [17], we expect that this moduli to be enumerable union of constructible set corresponding to possible choices of local exponents at infinity of Liouvillian solutions. With this purpose we explore the application of Kovacic’s algorithm (see [10, 12]) to the family (1). Some steps of the algorithm, dealing with polynomial solu- tions of auxiliar equations, are very sensitive to changes of the parameters. However, the Asymptotic Iteration Method (see [7]) allows us to describe the algebraic condi- tions on the parameters giving rise to the existence of Liouvillian solutions. The structure of the paper is as follows. Section 2 is devoted to the definitions of parameter space ℙ2n , spectral set 2n , spectral varieties 2n,d and the statement of our first main result, Theorem 2.3. Section 3 is devoted to the definition of pol- ynomial-hyperexponential solutions, the reduction of the parameter space through D’Alembert transformation, and Kovacic’s algorithm. The analysic of equation y�� = (x2n + xn−1)y allows us to prove the non-emptyness of the spectral varie- ties 2n,kn and 2n,kn+1 (Corolary  3.6). Section  4 contains the results of this paper related to Asymptotic Iteration Method. We find a sequence of differential polyno- mials Δd(a, b) in two variables that codify the equations of the spectral varieties 2n,d inpendently of n (Theorem 4.3). The proof of Theorem 2.3 is included at the end of the section. Section 5 is devoted to the Liovullian solutions of Schrödinger equations with polynomial potentials. We proof that the number of the values of the energy parameter allowing a Liouvillian eigenfunction is bounded by the arithmetic condi- tion which is a simple function of the coefficients of the potential (Theorem 5.2). 2 P arameter space ∑ Let us c∑onsider Eq. (1) with polynomial coefficients n P(x) = pjx j ∈ ℂ[x]≤n and j=1 2n Q(x) = q xjj ∈ ℂ[x]≤2n . We also take into account a non-degeneracy condition p2 − q2n ≠j=10 , which implies that the equation can be reduced to trace free form with n a polynomial coefficient of degree 2n. Thus, the parameter space corresponding to the family of Eq. (1) is, ℙ2n = ( 2 ℂ[x]≤n× ∈ ℂ[x]≤2n) − {p − q2n = 0},n that we consider an an affine algebraic variety of dimension 3n + 2 with affine coor- dinates p0,… , pn, q0,… , q2n . Our purpose is to describe algebraically the spectral set 𝕃2n ⊆ ℙ2n . That is, the set of equations in the family (1) admitting a Liouvillian solution. An important class of Liouvillian functions, specifically relevant for the integrability of Eq. (1) is the following. 1 3 São Paulo Journal of Mathematical Sciences Definition 2.1 A polynomial-hyperexponential function of polynomial degree d and exponential degree k is a function of the form u(x) = P (x)e∫ Ak(x)dxd , (2) with Pd(x) and Ak(x) polynomials of degree d and k respectively. Definition 2.2 The spectral subvariety 2n,d is the subset of 2n corresponding to equations in the family (1) having a polynomial-hyperexponentional solution of pol- ynomial degree d. Theorem 2.3 Let 𝕃2n ⊂ ℙ2n be the set of equations in family (1) having a Liouvil- lian solution, and 2n,d be the set of equations in family (1) having a polynomial- hyperexponential solution of polynomial degree d . The following statements hold: (a) For any fixed n ∈ ℕ there is an infinite set of values of d such that 2n,d is not empty. (b) If not empty, 2n,d is an algebraic variety of codimension ≤ n in ℙ . (c) For any d ≠ 2nk the algebraic varieties 2n,d and 2n,k are disjoint in ℙ2n. (d) Any compact subset of ℙ2n intersects only a finite number of algebraic varieties of the family {𝕃2n,d}d∈ℕ. Furthermore, ⋃∞ 2n = 2n,d. d=0 Therefore we conclude that 2n is a singular analytic submanifold of ℙ2n consisting in the enumerable union of pairwise disjoint algebraic varieties of codimension ≤ n in ℙ2n. In what follows we will deal with the proof of Theorem 2.3 and the calculation of the equations of the spectral subvarieties 2n,d in suitable coordinates. 3 L iouvillian solutions 3.1 Reduction of the parameter space As it is well known, Eq. (1) can be reduced to trace free form y�� = R(x)y (3) ( ) by means of D’Alembert transform 1u = exp − ∫ P(x)dx y , where 2 P(x)2 P�(x) R(x) = + − Q(x) . Note that the degree of R(x) is not greater than 4 2 max{deg(Q(x)), 2deg(P(x))} . Note that the family of equations of the form (3) with 1 3 São Paulo Journal of Mathematical Sciences R(x) of fixed degree 2n are parameterized by the space 𝕂2n = ℂ[x]2n of polynomials of degree 2n that we see as an affine algebraic variety of dimension 2n + 1 , parame- terized by the coefficients of R(x) and thus isomorphic to ℂ∗ × ℂ2n . Note that the family (3) is included in (1), where R(x) ∈ 2n corresponds to (0,−R(x)) ∈ ℙ2n . The D’Alembert transformation is a polynomial map in the coefficients of P(x) and Q(x) and it can be seen as a retract, P(x)2 P�(x) dal2n ∶ ℙ2n → 𝕂2n ⊂ ℙ2n, (P(x),Q(x)) ↦ R(x) = + − Q(x) ↦ (0,−R(x))4 2 of the natural inclusion 𝕂2n ⊂ ℙ2n . Taking into account that the ratio between u and y is the exponential of a polynomial, we obtain that (P(x),Q(x)) ∈ 2n,d if and only if R(x) ∈ 2n,d ∩ 2n . Therefore, the analysis of polynomial-hyperexponential solu- tions of a given polynom∑ial degree can be restricted to the trace free family 2n. Let us write 2nR(x) = rjxj . Equation (3) can be red√uced to the case of monic j=0 polynomial coefficient by the change of variables 2n+2 1x ↦ x which lead us to the r2n equation ⎜⎛ ⎞2�n−1 r ⎟ y�� = ⎜x2n j+ � xj⎟y, a ∗k ∈ ℂ , ai ∈ ℂ. (4) ⎝⎜ j=0 2n+2 1 ⎟ r2n ⎠ For the next step, let us consider 2n ⊂ 2n the family of Eq. (3) with monic polyno- mial coefficient. It is an an algebraic variety isomorphic to ℂ2n . Since the (2n + 2)-th root of r2n is an algebraic multivalued function of r2n , any equation in 2n has 2n + 2 different equivalent reductions in 2n . This can be seen as an algebraic correspond- ence C2n ⊂ 2n ×2n . This algebraic correspondence is a (2n + 2)-fold covering space of 2n by the first projection, 1 and the (2n + 2) monic reductions of the equa- tion of coefficient R(x) are given by 2(−1({R(x)})) . Note that R(x) is in 2n,d if and 1 only if so are any of its monic reductions. Therefore, if suffices to focus our analysis to equations in the family 2n. 3.2 K ovacic’s algorithm and adapted coordinates in  2n From now on let � = 2n ∩2n be the reduced spectral set consisting of equations 2n in  �2n having a Liouvillian solution, and let  = 2n,d 2n,d ∩2n be the reduced spec- tral variety consisting of equations in 2n having a polynomial-hyperexponential solution of polynomial degree d. Note that, since D’Alembert reduction does not affect the polynomial degree of polynomial-hyperexponential solutions then a differential equation in the family (1) has a polynomial-hyperexponential solution of polynomial degree d if and only if so has any of its monic D’Alembert reductions. Therefore, if 2n,d is a subvariety of ℙ2n then codim(𝕃2n,d,ℙ2n) = codim(𝕃� ,𝕄2n,d 2n). Here we will analyze the existence of Liovillian solutions of equations in the family 2n . This is done in terms of some known theoretical results obtained by 1 3 São Paulo Journal of Mathematical Sciences application of Kovacic’s algorithm [12]. A first step is to introduce a system of coordinates in 2n that fits our analysis of Eq. (3) better than the coefficients of R(x). The following Lemma that can be traced back to [15, p. 474], allows to decompose the monic polynomial R(x) in a suitable form for the application of the algorithm. Lemma 3.1 Every monic polynomial M(x) of even degree 2n can be written in one only way completing squares, that is, M(x) = A(x)2 + B(x), (5) ∑ with n−1 ∑ A(x) = xn + ajx j is a monic polynomial of degree n and n−1B(x) = b xj j=0 j=0 j is a polynomial of degree at most n − 1. According to the proof given in [1, Lemma 2.4, p. 275] it also clear that the decomposition map 𝕄2n → ℂ2n , R(x) ↦ (a0,… , an−1, b0,… , bn−1) where R(x) = A(x)2 + B(x) is a regular invertible polynomial map. Therefore, we may consider the coefficients of A(x) and B(x) as a system of regular coordinates is 2n . The following results gives us precise information about the sets ′ and 2n ′ . 2n,d Theorem 3.2 [1, Theorem 2.5, pp. 276] Let us consider the differential equation, y�� = M(x)y, (6) with M(x) ∈ ℂ[x] a monic polynomial of degree k > 0. Then its differential Galois Group G with coefficients in ℂ(x) falls in one of the following cases: 1. G = SL2(ℂ) (non-abelian, non-solvable, connected group). 2. G = ℂ∗ ⋉ ℂ (non-abelian, solvable, connected group). Furthermore, the second case is given if and only if the following conditions holds: 3. M(x) has even degree k = 2n, 4. Writing M(x) = A(x)2 + B(x) as in Lemma 3.1, the quantity ±bn−1 − n is a non- negative even integer 2d, d ∈ ℤ≥0. 5. There exist a monic polynomial Pd of degree d satisfying at least one of the fol- lowing differential equations, P�� + 2AP� − (B − A�)Pd = 0,d d (7) P�� − 2AP� − (B + A�)Pd = 0.d d (8) In such case, Liouvillian solutions are given by ∫ −2 ∫ Adxy = P e Adx e1 d , y2 = y1 � , or2 (9)P d 1 3 São Paulo Journal of Mathematical Sciences e2 ∫ Adx y1 = Pde − ∫ Adx, y2 = y1 � .P2 (10) d A careful read of Theorem 3.2 gives us the following. Corollary 3.3 The sets ′ and ′ in 2n satisfy the following. 2n 2n,d ⋃ 1. ∞� = � . 2n d=0 2n,d 2. 2n,d is contained in the hypersurface of 2n of equation b2 − (n + 2d)2 = 0.n−1 ⋃ Therefore, the sets ∞2n and 2n,d in ℙ2n satisfy 2n = d=0 2n,d. Proof 1. It is a consequence of the dichotomy of the Galois group. In case the group is not SL2(ℂ) it leads to a polynomial-hyperexponential solution. 2. It is a direct consequence of point 2 in the second part of Theorem 3.2, The last statement of the corollary is a consequence of the point 1. and the fact the the reductions process from ℙ2n to 2n preserves polynomial-hyperexponential solu- tions. ◻ 3.3 Canonical equation The following example: y�� = (x2n + xn−1)y,  ∈ ℂ. (11) that we refer to as canonical equation gives us some information about the non emptiness of the sets L2n,d for large d. Due to theorem 3.2, if (11) has a Liouvillian solution, the parameter  in the canonical coefficient x2n + xn−1 is forced to be a discrete parameter that can be  = 2d + n or either  = −2d − n , where d is a non- negative integer, which lead us to deal with two different equations, y�� =(x2n + (2d + n)xn−1)y, or (12) y�� =(x2n − (2d + n)xn−1)y. (13) Proposition 3.4 The differential equation (12) is integrable in the liouvillian sense if and only if, d = (n + 1)k or d = (n + 1)k + 1 where k is a non-negative integer. Proof The differential equation (12), is transformed into the Whittaker differential equation, 1 3 São Paulo Journal of Mathematical Sciences ( −2d−n 1 ) 1 4( ) 2 − 1 �� 2n+2 2n+2 W = − + W, (14) 4 z 4z2 n through the change of variables 2z = xn+1 , −y = z 2n+2W . Applying Martinet-Ramis n+1 theorem, see [13], we have that −2d − n 1 1 ± ± = + k, k ∈ ℤ , 2n + 2 2n + 2 2 ≥0 which left only two posibilities, d = (n + 1)k or d = (n + 1)k + 1 . ◻ It is easy to see that the change of variables made in above proof also trans- form the Eq. (13) into a Whittaker equation. Nevertheless this new equation will have parameters 2d+n 1 = and  = . So via Martinet–Ramis theorem we can 2n+2 2n+2 enunciate the following result analogous to the previous proposition. Proposition 3.5 The differential equation (13) is integrable in the liouvillian sense if and only if, d = (n + 1)k or d = (n + 1)k + 1 where k is a non-negative integer. Moreover, the solutions to the Eq. (11) can be explicitly written as xn+1 y (x) = P n+1d,n d,n(x)e , if = 2d + n, or n+1 (15)x − yd,n(x) = Pd,n(x)e n+1 , if = −(2d + n), where the polynomials Pd,n can be find by a Frobenius-like method. Having said that, it is a tedious process. In any case, for d = (n + 1)k we have that ∑d P (x) = xdd,n +  x d−j j ,wherej = 0forj ≠ (n + 1)m. (16) j=n+1 On the other hand, for d = (n + 1)k + 1 ∑d Pd,n(x) = x d + jx d+1−j,wherej = 0forj ≠ (n + 1)m + 1. (17) j=n+2 j j = (n + 1)m , j = (n + 1)m + 1∏  = 2d + 2 m (d+2−r(n+1))(d+3−r(n+1))− ∏r=1 −2(d+2−r(n+1)−n)−2d  = −2d − 2 m (d+2−r(n+1))(d+3−r(n+1))− r=1 2(d+2−r(n+1)−n)+2d Corollary 3.6 For any pair (n, d) of degrees with d ≡ 0 or d ≡ 1 mod (n + 1) there exist a monic polynomial M(x) of degree 2n such that the equation 1 3 São Paulo Journal of Mathematical Sciences y�� = M(x)y (18) has a polynomial-hyperexponential solution of exponential degree n + 1 and polyno- mial degree d; therefore 2n,d is non-empty. 4 A nalysis of auxiliary equations We refer to Eqs. (7) and (8) as auxiliary equations for Eq. (6). As it is stated in The- orem 3.2 the existence of a Liouvillian solution of Eq. (6) depends of the existence of a polynomial solution of the auxiliary equations. In what follows we will show that conditions for the existence of a polynomial solution Pd of given degree is alge- braic in the coefficients of A(x) and B(x), and therefore in the coefficients of M(x). 4.1 A symptotic iteration method The asymptotic iteration method or AIM was introduced by Ciftci et al in [7] as a tool to solve homogeneous differential equations of the form y�� = �0y + r0y (19) where 0 and r0 are smooth functions defined on a real interval. Nevertheless, the method is purely differential algebraic, so we can extent the result to differential rings of characteristic zero. By derivation of Eq. (19) we obtain a sequence of dif- ferential equations, y(j+2) = jy � + rjy (20) where the sequences {j}j∈ℕ and {rj}j∈ℕ are defined by the recurrence, = �j+1  + rj + 0j, r = r � j+1 + r0j.j j (21) and the sequence of obstructions, j = rjj−1 − jrj−1. We say that the AIM stabilizes at p > 0 if p = 0 . The following statement is a dif- ferential algebraic translation of [16, Theorem 1]. Theorem  4.1 Let 0 and r0 be elements of a differential field R of characteristic zero. If there exist p > 0 such that rp rp−1 = ∶= , (22) p p−1 then differential equation (19) has general solution, 1 3 São Paulo Journal of Mathematical Sciences y = u−1(c2 + c1), c1, c2 arbitrary constants, (23) in the extension R⟨, u, u−1, v, ⟩ where u, v,  are solutions of u� = u , v� = 0v y � = u2v, respectively. Proof By derivation of Eq. (19) we obtain, y(p+2) = �py + rpy, and from there ( )  y� rp p + y log(y(p+1) � p ) = ( ) . r  y� p−1 p−1 + yp−1 If condition (22) is satisfied, then we have  (y(p+1))� p = y(p+1). p−1 On the other hand, from the recurrence, we have, p = log(p−1) � +  + 0, p−1 and replacing into the above equation we obtain, (y(p+1))� = (log( )� +  + )y(p+1)p−1 0 . We have that y(p+1) = c1p−1uv is a general solution for this equation and finally we obtain y� + y = c1uv that yields the general solution of the statement. ◻ The AIM method tests whether the auxiliary equations have polynomial solu- tion. The following statement is a differential algebraic translation of [16, Theo- rem 2]. There is no difference in the proof, so we refer the reader to the original text. Theorem 4.2 Let 0, r0 be elements in a differential field R of characteristic zero that contains ℂ[x] . (i) If (19) has a polynomial solution of degree p, then p = 0 (ii) If pp−1 ≠ 0 and p = 0, then the differential equation (19) has a polynomial solution of degree at most p. 1 3 São Paulo Journal of Mathematical Sciences 4.2 Liouvillian solutions by means of AIM Let us proceed to the AIM of auxiliary equations (7) and (8). For Eq.  (7) we should start with + = −2A(x) and r+ = B(x) − A�(x) . By the recurrence law (21) 0 0 we have a sequence: [ + ] [ ]� [ ][ ] + + p+1 p −2A(x) 1  + = + + p r r B(x) − A�(x) 0 r+ p+1 p p A condition for the existence of a polynomial solution of degree at most p of (7) is the vanishing of the polynomial + = r++ − r+ + . We proceed analogously p p p−1 p−1 p with Eq. (8) obtaining sequences of polynomials r− , − and −. p p p In order to model this process, let us consider ℚ{a, b} the ring of differential polynomials in two differential variables a, b. We may consider the following ℚ -linear differential operator in the space of 2 by 2 matrices (Table 1). [ ] � −2a 1 ∶ Mat2×2(ℚ{a, b}) → Mat2×2(ℚ{a, b}), C ↦ (C) = C + Cb − a� 0 We consider the iterations of this map. If we give to the differential variables a, b the values of the polynomials A(x) and B(x) we obtain: ( ([ ])) [ ] p+1 1 0   (A(x),B(x)) = p p−1 . 0 1 rp rp−1 Let us define the sequence of universal differential polynomials, ( ([ ])) Δ = − det p+1 1 0 p ∈ ℚ{a, b}.0 1 As we will see this sequence {Δd}d∈ℕ of differential polynomials governs the Liouvillian integrability of Eq. (6) for any even degree 2n of M(x). Table 1 First values of the universal differential polynomials Δ n Δ0 b − � ( a ) Δ1 2a a �� − b� + 4ba� − 3(a�)2 − b2 Δ2 −9b 2a� − 15(a�)3 − 2(−3a��b� + 2(a2(a(3) − b��) + (a��)2) + (b�)2) +b(−a(3) + 2a(3b� − 5a��) + 23(a�)2 + b��) + a�(3a(3) − 2a(7b� − 9a��) − 3b��) + b3 Δ3 2(2a �� + 2a(b − 4a�) + 4a3 − b�)(a(4) − 4a2(b� − a��) + b(10a�� − 4b�) + 10a�b� + 8a3(b − a�) +2a(−a(3) − 14ba� + 12(a�)2 + b�� + 2b2) − 16a�a�� − b(3)) − (−5a(3) + a(26a�� − 10b�) +12a2(b − 5a�) − 10ba� + 21(a�)2 + 16a4 + 3b�� + b2)(−a(3) − 2a(b� − a��) + 4a2(b − a�) −6b(a�) + 5a�2 + b�� + b2) 1 3 São Paulo Journal of Mathematical Sciences Theorem 4.3 Equation (6) with M(x) = A(x)2 + B(x) has a polynomial-hyperexpo- nential solution of polynomial degree d if and only if, b2 = (n + 2d)2 and Δd(A(x),B(x))Δd(−A(x),B(x)) = 0.n−1 Therefore ′ is an algebraic subvariety of  contained in the union of the irre- 2n,d 2n ducible hypersurfaces of equations: bn−1 = 2d + n, −bn−1 = 2d + n. Proof Note that, by definition of the sequence Δd we have + = Δd(A(x),B(x)) . d Analogously, the application of the AIM to Eq.  (8) produces an obstruction, − . d Note that, because of the symmetry between Eqs. (7) and (8) − = Δd(−A(x),B(x)).d We need only to check that ± ±  ≠ 0 for the auxiliar equations. This comes d−1 d easily from the fact that ± = ±2A(x) is of bigger degree than r0 = B ± A� . Note 0 that Δd(A(x),B(x))Δd(−A(x),B(x)) is a polynomial in x, a0,… , an−1, b0,… , bn−1 . Its coefficients as a polynomial in x are the algebraic equations of the restricted spectral variety ′ in  2n,d 2n . ◻ Example 4.4 As a first example of AIM applications, let us consider an equation on 2 y�� = ((x + a )20 + b0)y. (24) An elementary traslation as x ↦ x + a0 reduces the determination of ′2 structure to an analysis of liouvillian-integrability conditions for quantum harmonic oscillator y�� = (x2 + b0)y. (25) These conditions are b2 = (2d + 1)∏2 and Δd(x, b0)Δd(−x, b0) = 0 . It is easy to verify 0 that Δ (x, b ) = Δ (−x, b ) = 2d+1 dd 0 d 0 d − k . Therefore,k=0 �  2,d ∶ (b0 + 2d + 1)(b0 − 2d − 1) = 0 (26) Let us note that for a given equation 2n , conditions bn+1 = 2d + n and bn+1 = −2d − n are mutually exclusive. In the first case, auxiliary equation (7) may have a polynomial solution but not (8). The opposite occurs in the second case. Therefore, we decompose the spectral variety ′ as the disjoint union of two com- 2n,d ponents � + +2n,d =  ∪ − . The first component  correspond to equations 2n,d 2n,d 2n,d whose auxiliary equation (7) has a polynomial solution of degree d and the second component − correspond to equations whose auxiliary equation (8) has a polyno- 2n,d mial solution of degree d. Definition 4.5 � + −2n,d =  ∪  , where2n,d 2n{,d + bn−1 = 2d + n = 2n,d Δd(A(x),B(x)) = 0, and 1 3 São Paulo Journal of Mathematical Sciences { − b  = n−1 = −2d − n 2n,d Δd(−A(x),B(x)) = 0. As in the previous Example 4.4 it is always possible to get rid of the coefficient an−1 by means of a translation in the x axis. Therefore is convenient to consider the sets, V± = {an−1 = 0} ∩ ±  . 2n,d 2n,d whose equations are easier to describe. For instance, in 4 and 6 we restrict our analysis to equations of the forms: y�� = ((x2 + a 20) + b1x + b0)y. (27) and y�� = ((x3 + a x + a )21 0 + b 2 2x + b1x + b0)y (28) respectively. The following calculations of the equations of V+ for n = 2, 3 and 2n,d small values of d, in Tables 2 and 3 is performed by means of the universal differen- tial polynomials Δd. 4.3 C odimension of the spectral variety As the degree in x of the polynomials Δp(A(x),B(x)) grows quickly with p and the degree of M(x) = A(x)2 + B(x) it seems that the sets 2n,p are smaller as p grows. However, a direct analysis of the auxiliary equations allows us to bound the codi- mension of the spectral varieties 2n,d in 2n . As we have seen before the algebraic Table 2 Algebraic equations of restricted spectral varieties + = +V  ∩ {a1 = 0} for small values of d 4,d 4,d { V+ 4,0 b1 = 2 b0 = 0 { V+ b = 4 4,1 1 b2 + 4a0 = 0 { 0 V+ b 4,2 1 = 6 b3 + 16a 0 0 b0 − 16 = 0 { V+ b = 8 4,3 1 b4 + 40a 2 2 0 0 b − 96b 0 0 + 144a = 0 { 0 V+ b 4,4 1 = 10 b5 + 80a b3 − 336b2 + 1024a2b − 3072a = 0 { 0 0 0 0 0 0 0 V+ b1 = 124,5 b6 − 140a b4 + 896b3 − 4144a2b20 + 28160a0b0 − 14400a 3 − 25600 = 0 { 0 0 0 0 0 0 V+ b = 14 4,6 1 b7 + 224a b50 − 2016b 4 + 12544a2b3 − 142848a0b 2 + 147456a3b + 288000b − 884736a2 = 0 0 0 0 0 0 0 0 0 0 0 1 3 São Paulo Journal of Mathematical Sciences Table 3 Algebraic equations of restricted spectral varieties +V = + ∩ {a1 = 0} for small values of d 6,d 6,d V+ 6,0 ⎧⎪ b2 = 3⎨ b = 0⎪ 1⎩ b0 − a1 = 0 V+ b = 5 6,1 ⎪⎧ 2⎨ 2a1b1 − 8a0 − 2b0b1 = 0 ⎩⎪ −6a1 − b 2 + 2b0 = 01 4a1b0 − 2a b − 3a 2 0 1 − b 2 = 0 1 0 V+ 6,2 ⎧⎪ b2 = 7⎪ 23a2b0 − 9a b21 − 14a0a1b1 + 6a0b 30b1 − 15a − 24a1 + 32a2 + b3 − 2b2 + 8b0 = 01 0 1 0 0 1⎨ 9a2b1 − 12a1b0b1 + 6a0b2 + 24a 2⎪ 1 1 0 b0 − 24a0a1 + 3b b1 − 12b1 = 00 −3a b2⎪ 1 + 36a1b0 + 24a0b1 − 30a 2 − 6b2 + 3b0b 2 − 48 = 0 1 1 0 1 ⎩ 22a b 31 1 − 32a0 + b − 6b b1 0 1 = 0 V+ 6,3 ⎧⎪ b2 = 9176a3b − 86a2b2⎪ 0 − 116a a 2 0 b1 + 16a1b 3 − 20a1b 2 + 264a1b0 + 80a0a1b b1 1 0 1 0 1 0 1 ⎪ − 12a2b2 − 144a2⎪ b0 − 12a b 2 0 b1 + 120a b − 105a 4 − 372a20 1 + 432a 2a 0 1 0 0 1 1 0 1 ⎪ − b4 − 36b2 + 8b b20 − 288 = 0⎪ 0 0 1⎪ 60a 3b − 92a2b b + 56a a b2 + 192a a b + 36a b2b + 96a b − 48a b2 1 1 1 0 1 0 1 1 0 1 0 1 0 1 1 1 0 0 ⎨⎪ − 24a b 2 0 0b − 288a 2b1 − 144a a 2 0 + 576a0 + 8b 3 − 4b3b = 0 1 0 1 1 0 1 ⎪ − 18a 2b2 + 372a2b0 − 132a1b 2 + 24a 2 3 1 1 1 0 1 b0b − 72a0a1b1 − 12a0b − 72a1 1 0b0b1 ⎪ − 252a3 − 288a + 12b3 − 6b2 2 2⎪ 1 b + 48b + 288b0 = 01 0 0 1 1⎪ 4a b31 − 48a 20b + 136a2b1 − 160a1b0b1 + 288a0b0 − 864a0a1 − 4b 30b + 24b2b1 1 1 1 0 1 ⎪⎪ + 192b1 = 0⎩ −52a 21b + 144a0b + 120a b − 252a2 − b41 1 0 + 12b b2 − 12b2 − 288 = 01 1 1 0 1 0 equations for 2n,0 are well expressed by the obstruction Δ0(a, b) = b − a� , so henceforth we will consider d > 0. Proposition 4.6 If ′2n,d is not empty, then codim(�2n,d,2n) ≤ n. ∑ Proof Now, let us suppose that dPd = pkxk is a solution to one of the following k=0 algebraic equations P�� ± 2AP� − (B ∓ A�)P = 0 d d d (29) ∑n ∑where nA = xn + a n−kn−kx and B = b n−kn−kx . Hence the coefficients of the k=1 k=1 polynomial 1 3 São Paulo Journal of Mathematical Sciences ∑ ( ∑ )( )d n ∑d k(k − 1)p xk−2 n n−k k−1k ± 2x + 2an−kx kpkx k=2 ( k=1 k=1∑ )(n ∑n−1 ∑ ) (30)d − b xn−kn−k ∓ nx n−1 + (n − k)a xn−1−kn−k p x k k = 0 k=1 k=1 k=0 in ℂ[x] vanish. This give place to a system of equations which are sufficient condi- tions for the existence of Pd, ⎡ a1 − b0 2a2 2 0 ⋯ 0 0⎢ ⎤ ⎡p ⎤⎢2a2 − b1 3a1 − b0 4a 6 ⋯ 0 0 ⎥ ⎢ 0 0 p ⎥ ⎢ ⋱ ∗ ∗ 1 ⎢ ⎥ ⎥ ⎢ ⎥ = 0.⋮ 0 0 0 0 ⋯ 2(d − 1) + n − bn−1 (2d + n − 1)a⎣ n−1 − bn−2⎥ ⎣⎢ ⎦⎥p 0 0 0 0 ⋯ 0 2d + n − b ⎦ dn−1 (31) We will denote the coefficient matrix of the system (31) by M± (A,B) . Note this d,n matrix has size (d + n) × (d + 1) and it also has the property ⎡ 0 ± ⎢ ⎤ ⎢ M ± (A,B) ⋮ ⎥ M (A,B) = d−1,n . d,n (32)⎢ ∗ ⎥⎣ ⎥0 ⋯ 0 2d + n ± bn ⎦−1 Remark 4.7 As there is no solution P of degree less than d, then rank(M± (A,B)) = d. d−1,n In order to determinate the codimension of ′ around a point (A0,B0) we shall 2n,d choose a suitable d × d submatrix D of M± (A0,B0) such that its determinant is dif-d−1,n ferent from zero. In addition, the vanishing of the determinants of the matrices set by adding one of the remaining n rows of M± (A0,B0) to D, generates n conditional d−1,n equations which guarantees the existence of a non-trivial solution to (31). ◻ 4.4 An example: case n = 3 As an useful example in order to illustrate further computes, specially for look- ing accurate spectral values on Schrödinger type problems, let us assume that A(x) = x3 + a 23,1x + a 2 3,2x + a3,3 and B(x) = b2,0x + b2,1x + b2,2 . So, the analysis on previous Sect. 4.3 for case n = 3 can be summarized with the following proposition. Proposition 4.8 A necessary condition for equation y�� − (2x3 + 2a 2 �3,1x + 2a3,2x + 2a3,3)y − ((b 2 2,0 + 3)x + (2a3,1 + b2,1)x + a3,2 + b2,2)y = 0 (33) in order to have a polynomial solution of degree d is b2,0 + 3 = −2d for d = 0, 1, 2,… 1 3 São Paulo Journal of Mathematical Sciences On the other hand, sufficient conditions are coded by the solutions of the linear system associated to the matrix ⎡0 0 0 ⎤ ⎢⎢1 1 1 1 ⎥⎢2 2 2 2  ⎥2 − ⎢ ⎥⋱ ⋱ ⋱ ⋱ ⋱M (A,B) = ⎥ d ⎢ (34)⎢ d−2 d−2 d−2 d−2 d−2 ⎥ ⎢ d−1 d−1  ⎥d−1 d−1⎢ ⎥d ⎣ d d ⎥ 0  ⎦d+1 d+1 where k = −a3,2(2k + 1) − b2,2, k = 0, 1, 2,… k = −2a3,3(k + 1), k = 0, 1, 2,… k = (k + 2)(k + 1), k = 0, 1, 2,… k = −2a3,1k − b2,1, k = 1, 2, 3… k = −2k − b2,0 + 1, k = 2, 3, 4… It creates a set of at most two polynomial equations in the variables a3,0 , a3,1 , a3,2 , a3,3 , a2,0 , a2,1 , a2,2 which guarantees likewise a non-trivial solution to the associated system to M−(A,B) and a polynomial solution to Eq. (33). d Proof This is a restriction of the analysis developed on Sect. 4.3 to n = 3 . ◻ Generically we can suppose that our d × d principal minor is different from zero, so the equations given by Proposition 4.8 are the determinants || ||0 0 0 || |1 1 1 1 | |||     | 2 2 2 2 2 || Δ1 (−A,B) = | ⋱ ⋱ ⋱ ⋱ ⋱ | = 0 d,3 | | (35)||  || d−2 d−2 d−2 d−2 d−2 | | || d−1 d−1 d−1 d−1 || d d d || and | | || |0 0 0 | |1 1  | 1 1 | || | 2 2 2 2 2 | 2 | |Δ (−A,B) = | ⋱ ⋱ ⋱ ⋱ ⋱ | = 0d,3 | (36) || | d−2 d−2 d−2 d−2 d 2 ||− | || d−1 d−1 d−1 d−1 || 0  |d+1 d+1 | 1 3 São Paulo Journal of Mathematical Sciences Several detailed examples of this equations can be found on [8, 19]. 4.5 Proof of Theorem 2.3 We can now state the proof, which follows easily from the other results. Statement (a) is a direct consequence of Proposition  3.5. Statement (b) is a consequence of Theorems  4.3 and Proposition  4.6. Note that, from d’Alembert reduction, the codimension of ′ in 2n coincide with that of 2n,d in ℙ2n . Statement (c) and 2n,d (d) are also clear, as ′ is contained in the union of hyperplanes of equations 2n,d bn+1 = 2d + n and bn−1 = −2d − n . ◻ 5 S chrödinger equation Let us summarize briefly the known results about explicit solutions for the one dimen- sional stationary Schrödinger equation. We start mentioning that Natanzon in 1971, see [14], introduced exactly solvable potentials, which today are known as Natanzon potentials. The seminal work of Natanzon inspired further researchers about exactly solvable potentials, although in the sense of Natanzon exactly solvable potentials also include potentials in where Schrödinger equations have eigenfunctions of hypergeo- metric type, not necessarily Liouvillian functions. The exactly solvable potentials, also known as solvable potentials, we extended to Schrödinger equations with explicit eigenfunctions. In this sense, solvable potentials are related to Schrödinger equations with eigenfunctions belonging to the set of special functions (Airy, Bessel, Error, Ei, Hypergeometric, Whittaker, Heun), not necessarily Liouvillian! Moreover, in case of Coulomb and 3D harmonic oscillator potentials correspond to Schrödinger equa- tions which are transformed into Whittaker differential equations, Martinet-Ramis in [13] established the necessary and sufficient conditions to determine the obtaining of Liouvillian solutions of the Whittaker differential equations. Recently Combot in [9] developed another method to obtain exactly solvable potentials, in the sense of Natanzon, involving rigid functions in the sense of Katz. To avoid confusion between explicit and Liouvillian solutions it was introduced the concept of algebraic spectrum in [2]. Also known as Liouvillian spectral set it is the set of eigenvalues for which the Schrödinger equation has Liouvillian eigen- functions, see also [3, 4]. In some scenarios it is known that bounded eigenfunctions of Schrödinger operator are necessarily Liouvillian, see [6]. Potentials with infinite countable algebraic spectrum are called algebraically solvable potentials and those with finite algebraic spectrum algebraically quasi-solvable potentials, for complete details see [4, §3.1, pp. 316] and see also [2, 3]. On the other hand, Turbiner in 1988, see [18], following the same philosophy of Natanzon, introduced quasi-solvable potentials. The seminal paper of Turbiner leaded to the seminal paper of Bender and Dunne in 1996, see [5], in where they obtain a family of orthogonal polynomials in the energy values of the Schrödinger equation with sextic anharmonic potentials, see also [11] for the study of more general sextic anharmonic oscillators. Due to Schrödinger equation with quartic anharmonic oscillator potential 1 3 São Paulo Journal of Mathematical Sciences falls in triconfluent Heun equation, see [10], it is in some sense a generalized Natan- zon potential (exactly solvable) although there no exist Liouvillian eigenfunctions. In a similar way for algebraically solvable potentials, in [4, §3.1, pp. 316] also was introduced the concept of algebraically quasi-solvable potential as those finite non empty algebraic spectrum, see also [2, 3]. Examples of algebraically solvable potentials and algebraically quasi-solvable potentials (quartic and sextic oscillators) were presented in [2–4] using [1, Theorem 2.5, pp. 276], which corresponds to the application of Kovacic algorithm for reduced second order linear differential equation with polynomial coefficients. Let us consider the one dimensional stationary Schrödinger equation  �� = ( − U(x)) (37) with a polynomial potential U(x). It is clear that the potential U(x) is algebraically quasi-exactly solvable if there are some values of  for wich equation (37) has a Liouvillian solution. This is equivalent to say that the line, {𝜆 − U(x) ∶ 𝜆 ∈ ℂ} ⊆ 𝕄2n paremeterized by  , intersects the spectral set 2n. As it is well know, and we examined in Example 4.4, any quadratic potential is quasi-exactly solvable (and more over, exactly solvable). It is also clear that any quasi-exactly solvable potential is of even degree. Let us assume from now on that U(x) is of degree 2n ≥ 4. We consider the decomposition −U(x) = A(x)2 + B(x) as in Theorem  3.2. We define the arithmetic condition of U(x) as the complex number, |bn−1| − n d = 2 where bn−1 is the coefficient of xn−1 the polynomial B(x) appearing in the unique decomposition −U(x) = A(x)2 + B(x) . Note that a necessary condition for U(x) to be quasi exactly solvable is its arithmetic condition to be a non-negative integer. In such case the intersection between the line: {𝜆 − U(x) ∶ 𝜆 ∈ ℂ} ⊆ 𝕄2n and 2n is confined to the spectral variety 2n,d. Let us consider the universal sequence of differential polynomials Δd ∈ ℚ{a, b} as in Theorem 4.3. The following lemma allows us to bound the number of admissible values of energy (for which the Schrödinger equation admits a Liouvillian solution) of any quasi-exactly solvable polynomial potential. Let us make clear that by the degree of a differential polynomial Δd in the variable b we mean its ordinary degree: that is we consider a, a�, a��,… , b, b�, b��,… as an infinite set of independent variables. Lemma 5.1 The degree of Δd in the variable b is at most d + 1. Proof Let us recall the differential polynomials d and rd appearing in the definition of Δd . Let us prove first: 1 3 São Paulo Journal of Mathematical Sciences (a) The degree of d in the variable b is small or equal to d+1. 2 (b) The degree of rd in the variable b is small or equal to d+2. 2 The degree of 0 = −2a in the variable b is 0 an the degree of r0 = b − a� in the vari- able b is 1. Therefore (a) and (b) hold for d = 0 . Now, from the recurrence law (21) we have that the degree in b of j+1 is at most that of rj and that the degree in b of rj+1 is at most a unit bigger that the degree of j . This proves (a) and (b). The degree of d is at most the maximum between the sum of the degrees of d and rd−1 and the sum of the degrees of d−1 and rd ; which is at most d + 1 . ◻ Theorem  5.2 Let U(x) be an algebraically quasi-solvable polynomial potential, and let d be its arithmetic condition. The number of values of the energy parameter  such that Eq. (37) has a Liouvillian solution is at most d + 1. Proof Generically, we may consider that U(x) has no independent term. Then the condition on  for the existence of a Liouvillian solution is the vanishing of Δd(A(x),B(x) + ) which is a polynomial in x of  . The number of values of  for which this polynomial vanish can not be greater than its degree in  . Clearly, the degree in  of Δd(A(x),B(x) + ) can not exceed the degree in b of Δd(a, b) which is bounded by d + 1 by Lemma 5.1. ◻ Example 5.3 In order to illustrate the procedures developed here let us consider the non-singular Turbiner potential Table 4 Spectral system of d  Schrödinger equation associated d to (38) {1 J = 1  = 0 { 3 J = 2 2 − 24 = 0 { 5 J = 3 3 − 128 = 0 { 7 J = 4 4 − 4002 + 12096 = 0 { 9 J = 5 − 5 + 9603 − 129024 = 0 { 11 J = 6 6 − 19604 + 7292802 − 26611200 = 0 { 13 J = 7 − 7 + 35845 − 29347843 + 438829056 = 0 1 3 São Paulo Journal of Mathematical Sciences U(x) = x6 − (4J + 1)x2 (38) where J is a non-negative integer. This potential has been studied in several papers, including [5]. Let ′d ⊂  be the set consisting of all possible values for J and 6,d  with polynomial hyperexponential solutions of polynomial degree d. In virtue of Theorem 4.3 it is a subvariety of V(2J − d − 1) . So, d shall only take non-negative odd values (Table 4). On the other hand, we can easily compute the equations of d through the univer- sal differential polynomial Δd(x3,−(4J + 1)x2 − ) for the auxiliary equation P�� − 2x3P� − (3x2 − (4J + 1)x2 − )Pd = 0.d d (39) For the case d = 1 we get the following equations ⎧ 2J − 2 = 0 ⎨⎪− 2 = 0 (40) ⎪ 2(−4J − 1) + 12 = 0⎩−(−4J − 1)2 − 8(−4J − 1) − 15 = 0. Taking into account above consideration we compute the first seven equations for d 6 F inal remarks In this paper we developed a technique to obtain Liouvillian solutions for parameter- ized second order linear differential equations with polynomial coefficients. In par- ticular case, we study the set of possible values of energy to get Liouvillian solutions of Schrödinger equations with anharmonic potentials. We adapted asymptotic itera- tion method, Kovacic’s Algorithm and previous results provided in [1–4] in terms of algebraic varieties extending slightly the known results about polynomial quasi- solvable potentials. Acknowledgements This research has been partially funded by Colciencias Project “Estructuras lineales en geometría y topología” 776-2017 code 57708 (Hermes UN 38300). We also thank the anonymous ref- eree for his/her valuable comments and suggestions. Compliance with ethical standards Conflict of interest On behalf of all authors, the corresponding author states that there is no conflict of interest. References 1. Acosta-Humánez, P., Blázquez-Sanz, D.: Non-integrability of some Hamiltonian systems with rational potential. Discrete Contin. Dyn. Syst. Ser. B 10(2–3), 265–293 (2008) 1 3 São Paulo Journal of Mathematical Sciences 2. Acosta-Humánez, P.B.: Galoisian Approach to Supersymmetric Quantum Mechanics. PhD thesis, Universitat Politècnica de Catalunya (2009). https ://www.tdx.cat/handle /10803 /22723 3. Acosta-Humánez, P.B.: Galoisian Approach to Supersymmetric Quantum Mechanics. The Integra- bility Analysis of the Schrödinger Equation by Means of Differential Galois Theory. VDM Verlag, Dr Müller, Saarbrücken, Deutschland (2010) 4. Acosta-Humánez, P.B., Morales-Ruiz, J.J., Weil, J.-A.: Galoisian approach to integrability of schrödinger equation. Rep. Math. Phys. 67(3), 305–374 (2011) 5. Bender, C., Dunne, G.: Quasi-exactly solvable systems and orthogonal polynomials. J. Math. Phys. 37(1), 6–11 (1996) 6. Blázquez-Sanz, David, Yagasaki, Kazuyuki: Galoisian approach for a Sturm–Liouville problem on the infinite interval. Methods Appl. Anal. 19(3), 267–288 (2012) 7. Ciftci, H., Hall, R., Saad, N.: Asymptotic iteration method for eigenvalue problems. J. Phys. A Math. Gen. 36(47), 11807–11816 (2003) 8. Ciftci, H., Hall, R., Saad, N., Dogu, E.: Physical applications of second-order linear differential equations that admit polynomial solutions. J. Phys. A Math. Theor. 43(41), 415206–415219 (2010) 9. Combot, T.: Integrability of the one dimensional Schrödinger equation. J. Math. Phys. 59(2), 022105 (2018) 10. Duval, A., Loday-Richaud, M.: Kovacic’s algorithm and its application to some families of special functions. Appl. Algebra Eng. Commun. Comput. 3(3), 211–246 (1992) 11. Hall, R., Saad, N., Ciftci, H.: Sextic harmonic oscillators and orthogonal polynomials. J. Phys. A Math. Gen. 39(26), 8477–8486 (2006) 12. Kovacic, J.: An algorithm for solving second order linear homogeneous differential equations. J. Symb. Comput. 2(1), 3–43 (1986) 13. Martinet, J., Ramis, J.-P.: Theorie de galois differentielle et resommation. In: Tournier, E. (ed.) Computer Algebra and Differential Equations, pp. 117–214. Academic Press, London (1989) 14. Natanzon, G.A.: Investigation of a one dimensional Schrödinger equation that is generated by a hypergeometric equation. Vestnik Leningrad. Univ 10, 22–28 (1971). in Russian 1 5. Rainville, E.D.: Necessary conditions for polynomial solutions of certain Riccati equations. Am. Math. Mon. 43(8), 473–476 (1936) 16. Saad, N., Hall, R., Ciftci, H.: Criterion for polynomial solutions to a class of linear differential equa- tions of second order. J. Phys. A Math. Gen. 39(43), 13445–13454 (2006) 1 7. Singer, Michael F.: Moduli of linear differential equations on the Riemann sphere with fixed Galois groups. Pac. J. Math. 160(2), 343–395 (1993) 18. Turbiner, A.V.: Quantum mechanics: problems intermediate between exactly solvable and com- pletely unsolvable. Soviet Phys. JETP 10(2), 230–236 (1988) 19. Venegas-Gómez, H.: Enfoque galoisiano de la ecuación de schrödinger con potenciales polinomi- ales y polinomios de laurent. Master’s thesis, Universidad Nacional de Colombia, sede Medellín (2018). http://bdigit al.unal.edu.co/71580 / Publisher’s Note Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations. 1 3