Stability boundaries of a Mathieu equation having PT symmetry

TL;DR

This paper derives analytic formulas for the stability boundaries of a PT-symmetric Mathieu equation, showing how non-Hermitian parameters influence solution stability, filling a gap in existing theoretical understanding.

Contribution

It provides the first analytic expressions for stability boundaries in a PT-symmetric Mathieu equation, enhancing understanding of parameter control over solution behavior.

Findings

Non-Hermitian parameter controls boundary shape and curvature

Analytic formulas for stability boundaries are derived

Bounded and unbounded solution regions are characterized

Abstract

I have applied multiple-scale perturbation theory to a generalized complex $PT$-symmetric Mathieu equation in order to find the stability boundaries between bounded and unbounded solutions. The analysis suggests that the non-Hermitian parameter present in the equation can be used to control the shape and curvature of these boundaries. Although this was suggested earlier by several authors, analytic formulas for the boundary curves were not given. This paper is a first attempt to fill this gap in the theory.

Equations

$$\frac{d^2\psi (x)}{d{x}^2}+(a+2\epsilon \cos x)\psi (x)=0,$$

$$\frac{d^2\psi (x)}{d{x}^2}+[a+2\epsilon V(x)]\psi (x)=0,$$

$$\psi (x,\tau )=\sum _{n=0}^{\infty }{\epsilon }^n{Y}_n(x,\tau )=Y_0(x,\tau )+\epsilon Y_1(x,\tau )+{\epsilon }^2{Y}_2(x,\tau )+...$$

$$\frac{d^2\psi }{d{x}^2}=\sum _{n=0}^{\infty }{\epsilon }^n\left(\frac{{\partial }^2{Y}_n}{\partial x^2}+2\epsilon \frac{{\partial }^2{Y}_n}{\partial \tau \partial x}+{\epsilon }^2\frac{{\partial }^2{Y}_n}{\partial {\tau }^2}\right)$$

$$a(\epsilon )=\sum _{n=0}^{\infty }{a}_n{\epsilon }^n=a_0+a_1\epsilon +a_2{\epsilon }^2+...$$

$${\epsilon }^0:\frac{{\partial }^2{Y}_0}{\partial x^2}+a_0{Y}_0$$

$${\epsilon }^1:\frac{{\partial }^2{Y}_1}{\partial x^2}+a_0{Y}_1$$

$${\epsilon }^n:\frac{{\partial }^2{Y}_n}{\partial x^2}+a_0{Y}_n$$

$$\frac{{\partial }^2{Y}_0}{\partial x^2}=0$$

$$\frac{{\partial }^2{Y}_1}{\partial x^2}=-h(\tau )[a_1+2(\cos x+i\beta \sin x)].$$

$$Y_1=2h(\tau )(\cos x+i\beta \sin x)+g(\tau ),$$

$$\frac{{\partial }^2{Y}_2}{\partial x^2}$$

$$-4ih(\tau )\beta \sin 2x-2[g(\tau )+2i\beta h^{\prime }(\tau )]\cos x-2[ig(\tau )\beta -2h^{\prime }(\tau )]\sin x.$$

$$h^{\prime \prime }(\tau )+[a_2+2(1-{\beta }^2)]h(\tau )=0$$

$$a_2$$

$$a_2=-2(1-{\beta }^2)$$

$$a(\epsilon )=-2(1-{\beta }^2){\epsilon }^2+O({\epsilon }^3),$$

$${\kappa }_1=\left|\frac{d^{2}a(0)}{d{\epsilon }^2}\right|=4(1-{\beta }^2)$$

$$Y_2=\frac{h_0(1+{\beta }^2)}{2}\cos 2x+i\beta h_0\sin 2x+2g(\tau )\cos x+2ig\beta (\tau )\sin x+f(\tau ),$$

$$\frac{{\partial }^2{Y}_0}{\partial x^2}+\frac{Y_0}{4}=0.$$

$$Y_0=A_1(\tau )e^{ix/2}+A_2(\tau )e^{-ix/2}$$

$$\begin{array}{rl}\frac{{\partial }^2{Y}_1}{\partial x^2}+\frac{Y_1}{4}& =[-a_1{A}_1-i{A}_1^{\prime }-(1+\beta )A_2]e^{ix/2}\\ & +[-a_1{A}_2+i{A}_2^{\prime }-(1-\beta )A_1]e^{-ix/2}\\ & +-(1+\beta )A_1{e}^{3ix/2}-(1-\beta )A_2{e}^{-3ix/2}\end{array}$$

$$a_1{A}_1+i{A}_1^{\prime }+(1+\beta )A_2$$

$$a_1{A}_2-i{A}_2^{\prime }+(1-\beta )A_1$$

$$A_j^{\prime \prime }=[(1-{\beta }^2)-a_1^2]A_j,$$

$$a_1^2$$

$$a(\epsilon )=\frac{1}{4}\pm \sqrt{1-{\beta }^2}\epsilon .$$

$$Y_1=\frac{(1+\beta )}{2}{K}_{+}{e}^{3ix/2}+\frac{(1-\beta )}{2}{K}_{+}{e}^{-3ix/2}+B_1(\tau )e^{ix/2}+B_2(\tau )e^{-ix/2},$$

$$\frac{{\partial }^2{Y}_2}{\partial x^2}+\frac{Y_2}{4}$$

$$+e^{-ix/2}\left[-\sqrt{1-{\beta }^2}{B}_2-a_2{K}_{-}+i{B}_2^{\prime }-\frac{(1-{\beta }^2)}{2}{K}_{-}-(1-\beta )B_1\right]$$

$$+e^{3ix/2}\left[-\frac{\sqrt{1-{\beta }^2}(1+\beta )}{2}{K}_{+}-(1+\beta )B_1\right]$$

$$+e^{-3ix/2}\left[-\frac{\sqrt{1-{\beta }^2}(1-\beta )}{2}{K}_{-}-(1-\beta )B_2\right]$$

$$-\frac{(1+{\beta }^2{)}^2}{2}{K}_{+}{e}^{7ix/2}-\frac{(1-\beta {)}^2}{2}{K}_{-}{e}^{-5ix/2},$$

Full text

††journal: Physics Letters A

1 Introduction

Non-Hermitian Hamiltonians having $PT$-symmetry were introduced in 1998 by Bender and Boettcher in a first attempt to generalize the standard postulates of Hermitian quantum mechanics to include operators invariant under parity ($\mathcal{P}$) and time reversal ($\mathcal{T}$) symmetries [1]. The authors showed that if a Hamiltonian commutes with the $PT$ operator, the energy eigenvalues could consist of entirely real, discrete and positive numbers. However, this is only a necessary but not a sufficient condition for the reality of the spectrum. If a $PT$-symmetric Hamiltonian is found to possess complex eigenvalues, we say that it has a broken symmetric phase [2, 3, 4]. If a $PT$-symmetric Hamiltonian is found to possess real spectra, then it is possible to introduce a new symmetry, described by an operator $\mathcal{C}$, such that the evolution of the entire system is unitary, consistent with physical requirements of probability conservation [3]. The condition imposed on the potential function $V(x)$ of the Schrödinger equation to guarantee $PT$ symmetry is that $V(x)=V(-x)^{*}$, as can be easily verified. This, in turn, implies that the real (imaginary) part of the potential must be even (odd) under $x\rightarrow-x$. Reviews on the subject of $PT$ symmetry in quantum and classical mechanics are already available [5, 6].

The Schrödinger equation with a periodic potential $V(x)=V(x+P)$ ($P$ being the period) having parity and time-reversal symmetry was initially considered by Bender and coauthors [7] and also by Z. Ahmed [8]. The authors found that real energy bands can be created even though the potential is a complex function of position. Moreover, the bands could be controlled by tuning the non-Hermitian parameter present in the Hamiltonian. There is also a vast literature dealing with the dynamics of wavepackets propagating in periodic potentials with $PT$ symmetry [9, 10, 11, 12, 13, 14, 15]. Apart from a purely mathematical context, the Mathieu equation is one of the most studied second order linear differential equation in applied mathematics [16, 17]. Its solutions, known as Mathieu functions, play a key role in physics with special interest in the propagation of electrons in crystal lattices, where the potential function of the Schrödinger equation reduces to a Mathieu type. The Mathieu equation is a special case of a more general second order linear differential equation known as a Hill’s equation [18]. Given a differential equation of a Mathieu type,

[TABLE]

where $a$ and $\varepsilon$ are real and positive constants, it is well known that, depending on the numerical values of $(a,\varepsilon)$, the solutions can be bounded or unbounded with respect to $x$. One of the main objectives is to find the boundary lines in the $(a,\varepsilon)$-plane which separates bounded from unbounded solutions. The problem of finding the stability boundaries of a Mathieu equation is well documented [19]. The fact that one can control the band diagram by tuning the non-Hermitian character of the lattice is also not a new result and it was suggested earlier by several authors [7, 9, 11, 20, 21]. There are also some very interesting connections between the Mathieu equation and the E2 algebra in non-Hermitian models [22, 23, 24]. However, to my knowledge, a detailed analysis aimed to obtain analytic expressions for the stability boundaries in the $(a,\varepsilon)$ plane for a $PT$ symmetric lattice has not been considered. Therefore, the objective of this paper is to close this gap by deriving the first set of approximate expressions that can be used to understand in a clear way the shape and curvature of the stability boundaries of a generalized Mathieu equation with $PT$ symmetry.

2 Multiple-Scale Analysis: Setting up the problem

The primary goal is to find the boundaries between stable and unstable solutions in the $(a,\varepsilon)$-plane of the generalized $PT$-symmetric Mathieu equation and to understand their behavior as the non-Hermitian parameter $\beta\in[0,1]$ varies. The choice for this range of $\beta$ is dictated by the fact that when $\beta>1$ the Schrödinger equation is no longer normalizable under a $CPT$ inner product and therefore it is not possible to construct a physical theory for values outside this range [3]. In other words, the Hamiltonian has a broken phase when $\beta>1$ [25]. To this end, define the $PT$-symmetric Mathieu equation

[TABLE]

where the periodic potential $V(x)=\cos x+i\beta\sin x$ has $PT$ symmetry since it satisfies $V(x)=V(-x)^{*}$ as discussed previously. If $\beta=0$, the usual Mathieu equation is recovered, whose solutions converge for $a_{n}\neq n^{2}/4$ ($n=0,1,2,...$ and $\varepsilon\rightarrow 0$), as can be easily verified by using regular perturbation expansions around small values of $\varepsilon$ [19]. Multiple-scale perturbation analysis starts by assuming that the solution can be expanded in a perturbation series of the form

[TABLE]

where $\tau=\varepsilon x$ is a long spatial scale (treated as an independent variable) since $\tau$ is not negligible when $x$ is of order $1/\varepsilon$ or larger. By taking the total second order derivative of (3) with respect to $x$, we obtain the expansion

[TABLE]

We also expand the parameter $a$ in a perturbation series in powers of $\varepsilon$:

[TABLE]

where $a_{j}$ are real numbers. It should be mentioned that I am assuming that expansions (3) and (5) are a good representation of the true solution at least as long as $\varepsilon\rightarrow 0$. Since multiple-scale perturbation performs a global analysis, we expect the results to be valid at every point on the $x$ axis. In general, it is only through numerical methods applied to Equation (2) that we can really obtain any sort of information about the errors in the expansion [19]. However, since my objective is to obtain analytical formulas for the stability boundaries, numerical methods are usually not enlightening in this matter and that is why I chose to work with perturbation expansions. After substituting (3), (4) and (5) into (2), and collecting powers of $\varepsilon$, we obtain the following system of coupled linear partial differential equations:

[TABLE]

Equations (6), (7) and (8) provide the starting point for the analysis. In what follows I am going to consider the expansion around two particular values of $a_{0}$: {0, $\frac{1}{4}$} in order to find the stability boundaries between stable and unstable solutions in the vicinity of these specific values. As I pointed out at the beginning of this section, the solutions are stable if $a_{n}\neq n^{2}/4$ ($n=0,1,2,...$) in the vicinity of $\varepsilon\rightarrow 0$ and therefore the stability lines cut the $a$-axis at these points. The analysis presented here can be easily applied to any value of $a_{0}$ and the above two values are considered because the formulas for the stability lines around [math] and $1/4$ are already known (for the Hermitian case) and so we can compare the results obtained here with the literature [19].

3 Perturbative analysis for $a_{0}=0$

The structure of the set of coupled equations (6)-(8) is such that if $Y_{0}$ is found from (6), $Y_{1}$ can be calculated from (7) and all the successive orders can be found from (8) where the right hand side depends on previous calculated values. Therefore, we first set $a_{0}=0$ in (6) such that the function $Y_{0}$ satisfies

[TABLE]

and the bounded solution is $Y_{0}(x,\tau)=h(\tau)$ with $h(\tau)$ an arbitrary bounded function of $\tau$. To find $a_{1}$ we substitute the above $Y_{0}$ into (7) to obtain

[TABLE]

If we impose the condition that $Y_{1}$ be a bounded function, we must take $a_{1}=0$. Therefore, we have the $a$ expansion up to second order in $\varepsilon$: $a=a_{2}\varepsilon^{2}$. To find $a_{2}$ we need the solution of the differential equation (10) with $a_{1}=0$. This is easily done and the result is given by

[TABLE]

where $g(\tau)$ is an arbitrary bounded function of $\tau$. The differential equation (8) satisfied by $Y_{2}$ can be written as

[TABLE]

If we impose the condition that $Y_{2}$ remain bounded, we must have

[TABLE]

so that

[TABLE]

and we can take

[TABLE]

as the boundary stability line with $h(\tau)=h_{0}$ a constant value. This gives us the second order contribution to $a$:

[TABLE]

so that the curvature $\kappa$ of the boundary is

[TABLE]

and the function $Y_{2}$ (which is $PT$-symmetric if we take $g$ as a real function of $\tau$) is given by

[TABLE]

where $f(\tau)$ is an arbitrary bounded function of $\tau$. This process can be repeated if one wishes to find the third order correction $a_{3}$ and so on 111It can be shown that $a_{3}$ is, in fact, zero.. It is seen from (18) and (19) that the non-Hermitian parameter $\beta$ quadratically controls the curvature of the stability boundary for $\varepsilon$ near zero. In particular, if $\beta$ assumes the unit value, which is the critical number in which the symmetry of the lattice is spontaneously broken, the curvature vanishes. Figure 1 shows a plot of the stability line for three values of $\beta$ in the plane $(a,\varepsilon)$. When $\beta=0$ we recover the usual approximation found in textbooks that $a=-2\varepsilon^{2}$, shown as a continuous line in the figure. As soon as $\beta\neq 0$ the stability line deforms in such a way that it approximates to a straight line when $\beta$ is close to 1. The introduction of the non-Hermitian parameter therefore has a very important consequence on the position in the $(a,\varepsilon)$-plane of stable and unstable solutions.

4 Perturbative analysis for $a_{0}=1/4$

Let us now consider a slightly more complicated situation. If we choose $a_{0}=\frac{1}{4}$ then, according to (6), the function $Y_{0}$ satisfies

[TABLE]

The general solution of this linear partial differential equation is

[TABLE]

with $A_{1}$ and $A_{2}$ arbitrary bounded functions of the parameter $\tau$. According to (7), the function $Y_{1}$ satisfies the differential equation

[TABLE]

In order to obtain bounded solutions, the secular terms must vanish. This implies the following coupled second order differential equations for $A_{1}$ and $A_{2}$:

[TABLE]

Both, $A_{1}$ and $A_{2}$, satisfy the same differential equation, namely

[TABLE]

and it is easy to see that we have two conditions on the boundness of $A_{j}$:

[TABLE]

so that the stability boundary line in the $(a,\varepsilon)$ plane (up to first order in $\varepsilon$) is given by

[TABLE]

To find the second order correction, $a_{2}$, we notice that $Y_{1}$ is given by

[TABLE]

where $B_{j}$ are arbitrary bounded functions of the parameter $\tau$. After substituting (30) into (8) we get

[TABLE]

and imposing the condition that the secular terms vanish we have

[TABLE]

and

[TABLE]

which implies that $B_{j}$ with $j=1,2$ both satisfy the same differential equation

[TABLE]

For a bounded solution to exist, we must impose that

[TABLE]

with $B_{j}$ a constant independent of $\tau$. Equation (35) is consistent with the literature when $\beta=0$ [19]. Therefore, we can write the boundary line up to second order in $\varepsilon$ as

[TABLE]

The curvature $\kappa_{2}$ of this stability line is now given by

[TABLE]

which is always less than the curvature of the stability line around $a_{0}=0$. Figure 2 shows the stability lines as given by (36). One can clearly see from the figure that the boundaries can be tuned depending on $\beta$. Once this non-Hermitian parameter acquires physical interpretation (for example, it can be considered as the physical parameter which controls loss and gain in waveguide arrays in optical systems), the above approximated equations for the stability boundaries can be used to engineer a specific function on the Bloch waves. It should be noted that these lines are a global property of the system, implying that Figures 1 and 2 are composing the same system for the same values of $\beta$.

5 Conclusions

The introduction of a non-Hermitian term into the Mathieu equation (1) [giving a generalized complex Mathieu equation (2)] opens new possibilities in the exploration of the band structure of physical systems. It is well known that the Schrödinger equation with a periodic potential is the main tool behind solid state physics, and, within the independent electron approximation, the band structure of the material determines, for example, its conductivity and optical properties. Therefore, it is very important to understand the topological aspects of the band structure and we analytically provide the first set of equations that may be used to characterize the physical material [see Equations (18) and (36)]. Although this result was verified numerically by a number of authors, the exact form of these stability lines were not known. I have also shown that the curvature of these lines (as $\varepsilon\rightarrow 0$) are related to the non-Hermitian parameter $\beta$. Thus, in principle, if one can measure the curvatures in some experimental arrangement, the determination of $\beta$ could be obtained in laboratory conditions. Future work will be aimed at the deformation of the stability lines as the Mathieu equation becomes nonlinear and non-Hermitian. These results will be published in near future.

Acknowledgements

The author would like to acknowledge the Brazilian Agencies CNPq, CAPES and FAPEAL for partial financial support.