Amplitude-phase method for solving Floquet-type problems

The understanding of one-dimensional quantal and dynamical problems is fundamental to problems involving stability of material constructions, quantum tunneling in solids, chemical selection of gas components, electronic properties of material nanostructures. An amplitude-phase method for one-dimensional Schrödinger/Hill-type equations with periodic potentials/coefficients is shown to provide detailed insights into Floquet-type quantal and dynamical problems. Method-independent real and periodic local amplitude and wave-number (or angular frequency) functions are found for bounded Floquet solutions. Regular and weakly singular periodic potentials/coefficients of analytic forms apply.


Introduction
Second-order differential equations in one space dimension with periodic coefficient/potential behaviors are relevant in many mathematical, physical and chemical disciplines. Such equations similar to a quantum mechanical Schrödinger equation with periodic potentials appear in: stability analysis of periodic motions obtained in nonlinear dynamical systems [1][2][3][4][5], chemical problems of rotating molecules subject to periodic angular potentials [6], studies of electronic properties of solid structures [7][8][9][10][11][12]. Studies aiming to understanding relativistic electronic corrections of realistic arrays of atoms and molecules are in progress [13,14].
The present topic has a well established mathematical foundation based on the Floquet theory [15,16] and there are numerous physical and chemical applications [6]. Several standard numerical and various approximate methods are known. Approximations are often based on perturbation theory [4] and semiclassical techniques [6,17].
A new approach based on an old method is presented here. The present idea of a decomposition of waves in terms of an amplitude and a phase dates back to Milne [18] and Wheeler [19], and even earlier in the development of Wentzel-Kramer-Brillouin-type approximations [17]. The amplitude-phase approach replaces a real-valued wave function with two functions; a local phase function and a local amplitude function. Numerical computations are efficient and results are expressed in terms of a phase and an amplitude. These two quantities contain the essence of an oscillation or a wave caused by the physical parameters. The original wave is reproduced without approximations other than the limitation of the numerical accuracy [20,21]. One must realize that real-valued, local phase and amplitude functions are not unique quantities in general, unless the wave function is harmonic. However, It will be pointed out that they can be uniquely (independent of methods) defined for all solutions that are bounded.
The amplitude-phase method [18,19] applied to Floquettype problems [15,16] has not been found in the research literature by the author. With the amplitude-phase method, physical behaviors are understood in terms of amplitude and phase functions. Its relation to Floquet theory is close but not perfect; an arbitrary Floquet solution function is described by m P z z exp ( ) ( ), with a complex periodic quantity P(z) and a complex z-independent constant μ, where z is a variable. A corresponding amplitude-phase solution is expressed by A z p z exp i ( ) ( ( )) with a real-valued, positive amplitude function A(z) and a real phase function p(z). Physically acceptable solutions with Re(μ)=0 make these solutions more similar, although μ is constant and p(z) is a function of z. This similarity is explored in the present investigation. Section 2 presents the second-order Schrödinger-like differential equation for z-symmetric periodic potentials that may be weakly singular. In section 3 the basic amplitude-phase approach is presented and fundamental as well as principal fundamental solutions are constructed. Section 4 analyses transitions of principal fundamental solutions from one period to the next in terms of a period map. Characteristic phase propagations per period and Lyapunov exponents are defined and illustrated using periodic amplitude functions. Conditions for bounded solutions and periodic amplitude functions are specified and calculated. Section 5 discusses an application to weakly singular potentials. Conclusions are in section 6.

Second-order differential equation
The dynamical equations are formulated with an independent coordinate z, as if being a space coordinate. As a result of this choice many quantities become closely related to those of quantum mechanics. Interpretations in terms of quantities related to time-dependent dynamical systems should still be straight forward.
Many important equations in mathematics, physics and chemistry may be reduced to a dimensionless equation where V(z) may in quantum mechanics correspond to a potential energy function and E may correspond to a total, constant energy. Differentiation in (1) is performed with respect to a dimensionless z-axis. In cases V(z) is periodic, equation (1) is a special case of a so called Hill equation [22]. Most solution properties of (1) can be analyzed by studying a principal fundamental solution in the interval 0za, where a is such that V(a)=V(0) [16]. The z-axis is arbitrarily defined with the origin z=0 located at an overall minimum of V(z) in the range 0z<a.
From Floquet theory [15,16] it is known that (1) has two independent solutions  F z ( ) ( ) of the form where P(z) is a complex-valued periodic function and μ is a characteristic exponent. P(z) and μ are often computed analytically by introducing non-analytic potential models like the Kronig-Penney ones [12].
A solution of the form (2) satisfies for any z on the real axis.

Model potentials
Numerical illustrations of the amplitude-phase method, to be introduced in section III, are computed with potentials V(z) chosen to be of two types. One type is of the harmonic type, i.e.
with v 0 =1 and 8. This type of potential is seen in many publications and is related to the Mathieu equation [6,16]. Note that z=0 is a minimum of V(z) in (4) and reflection symmetries of periodic solutions are thereby opposite to those with the origin defined at a maximum of the potential. The proposed method applies also to a potential model with periodically separated weak singularities. The potential (5) is a simple model of a light particle (electron) moving along an array of attracting heavy particles (ions). Note that z=0 is at a maximum of V(z) in (5). Both potential models are illustrated in figure 1. 3. Amplitude-phase method

Principal fundamental solutions
A principal fundamental solution matrix is defined as a unit matrix at some boundary point [16], a point which is chosen at the origin = z 0. It is convenient to introduce the amplitude-phase solutions in complex form and later construct the real-valued linear combinations of them.
Two independent amplitude-phase solutions of equation (1) for real E and real V(z) on the real z-axis are written as where A(z) and p′(z) are positive scalar functions [20,21]. The amplitude function A(z) is the key quantity in the present amplitude-phase analysis. This is so because a given amplitude function defines a corresponding phase function p(z) by the differential relation except for an integration constant. The Wronskian of F (±) (z) is constant [20]. Unphysical solutions are characterized by amplitude functions A(z) diverging as  ¥ z | | . This follows from the definition and requirement of a real-valued phase function p (z) in (6) and (7).
If boundary conditions are such that A(0) is positive, then A(z) stays positive for all z. The reason is the 'repulsive' character (with respect to A(z)) of the nonlinear component The boundary conditions for A(z) are taken at z=0. The potentials used are symmetric with respect to z=0. Special care is taken to guarantee a real positive value A(0) and a positive function p′(z) in (7). For this purpose two different choices of boundary conditions are suggested, depending on the sign of the coefficient value E−V(0) in (8). From semiclassical quantum theories one knows that for E−V > 0 the corresponding z-region is classically allowed and quantal waves are locally oscillating.
For E−V < 0 the region is classically forbidden and quantal waves have exponential behaviors as functions of z.
For singular attractive, periodic potentials, like in (5), the origin z=0 is defined at a position of a main potential maximum (within a period range 0za). When the energy difference E−V(0) is negative, parts of the region 0za is classically forbidden. Then a general boundary condition can be taken as which allows the (generally non-periodic) amplitude function to be symmetric with respect to = z 0. When the energy difference E−V(0) is positive, then the entire region 0za is classically allowed. In that case an adiabatic approximation of the nonlinear equation (8), i.e.
can be utilized. This approximation is periodic, but the exact amplitude functions will only be approximately periodic if (10) is a good approximation. Suggested boundary conditions at z=0 for (8) are taken as The condition (11) can be used for all relevant energies if the periodic potential is finite like in (4) and the origin is defined at a potential minimum. Note also that condition (11) imply A′(0)=0 in (8). Condition (11) merely select reasonably slowly varying amplitude functions, being non-periodic in general.
To form a principal fundamental solution matrix one can normalize any real-valued fundamental solution matrix [16]. A real-valued fundamental solution matrix, denoted g(z), is chosen as It consists of real and imaginary parts of the exponential representation (6). The derivatives of the matrix elements in the first row in (12) constitute the second row, and are obtained by using (7). This fundamental solution matrix provides a matrix g(0) for the boundary values at z=0, given by (9) and (11), and p(0)=0 as in (6). It also provides two matrices g(±a) at the neighboring equivalent positions z=±a in the potential: By normalizing the fundamental solution matrix in (12) to satisfy a unit matrix at z=0 instead of satisfying (13), one obtains a representation of the principal fundamental solution M(z) satisfying g(z)=M(z)g(0) and given by For arbitrary values of z in (15), the first row are the principal fundamental solutions with their derivatives expressed in the second row. The amplitude function A(z) is a symmetric function with respect to z=0 for the potentials (4) and (5)  The main assumption in the derivation of (15) is that the amplitude functions are consistent with A′(0)=0 and p (0)=0, which is always possible. The value A(0) is not essential in theory, only suggested above. Any value leads to exact wave functions and to a unique principal fundamental solutions matrix M(z). The arbitrariness of A(0) does not change any of the matrix elements of M(z), which implies that relations between representations due to specific choices of A (0) can be obtained. In the subsequent analysis a relation to representations based on possible periodic amplitude functions is investigated.

Periodic amplitude functions
Condition (9) or (11) do not provide boundary conditions for exactly periodic amplitude functions. Exact periodic amplitude functions A p (z) must be obtained by other algorithms or formulas. It is shown below that periodic amplitudes A p (z) exist only when a general solution is bounded.
The reason for amplitude functions A p (z) not being periodic when m ¹ Im 0 has to do with the definition of the phase function p p (z) related to A p (z) by equation (7). Any exponential behaviors of amplitude-phase waves must be due to A p (z). Note that the periodic amplitude functions are always positive, real-valued functions, in contrast to the periodic functions P(z) of the Floquet theorem [16] in (2). By comparing a scalar Floquet solution with a corresponding amplitude-phase solution, for example  (15)) cannot be periodic [16], and A p (z) cannot exist.
If periodic amplitude functions exist, they define bounded solutions equivalent to those defined by amplitude functions specified by boundary conditions (9) and/or (11). Figure 2 illustrates for the potential (4) with v 0 =1, that A(z) satisfying (11) and A p (z) are different (right subplots), but provide the same (principal) fundamental waves in (15) (left subplots). The upper subplots are obtained for an energy near the lower limit of an energy interval of unphysical solutions, a so-called gap. A gap is here associated with unbounded solutions (M 11 and/or M 12 ) and their corresponding energies. The lower subplots are obtained for an energy near the upper limit of the same gap. Both energies almost satisfy conditions for period-2a solutions of (1). At each energy only one of the fundamental solutions is periodic with period 2a (a=π). By comparing the upper and lower subplots, the periodic amplitude function A p (z) is seen to change character from one side of the gap to the other; an indication of its singular behavior.

Period map and primary resonances
A period map is here the special value a M( ) of the principal fundamental solution, or possibly the value M(−a). The behavior of solutions as  ¥ z | | can be analyzed in terms of a period map M(a) (and/or M(−a)) in at least two ways with the present approach; by using the amplitude functions obtained from boundary condition (9) or (11), or by using periodic amplitude functions. Both representations of the period map must be equivalent when they exist.

Representation with periodic amplitude functions
In special cases of existing periodic amplitude functions, one obtains from (15) 0 cos 0 sin 0 cos 0 sin 0 sin 0 cos . 1 5 which follows by induction based on trigonometric formulas for cos/sin-functions of two angles. One realizes that = = np a p na n , 1, 2, 3, , 20

Exact representation
To generalize the applicability of the period map M p (a) in (18) one can relate it to the principal fundamental solution M (a), obtained from the general formula (15). The matrix (15) is valid for arbitrary amplitude functions satisfying A′(0)=0 as in (9) or (11) with p(0)=0, and does not rely on further requirements. This means, that with z=a in (15) one can compare the matrix elements obtained from amplitude functions satisfying (9) or (11) with those based on periodic amplitudes in (18). All corresponding elements must be equal whenever A p (z) exists. This matrix equivalence allows a generalization of p p (a) and a formula for calculating A p (0). can be used as a generalization of p p (a), to include imaginary values. The real part of p p (a) is chosen to lie in the same angular branch as that of p(a). This means that ' p a a i p ( ) ' becomes the general Floquet exponent 'μ' (see equation (3)). The quantity p p (a) can then be used for calculating μ for all types of solutions satisfying (25) and (22), independent of an existing periodic amplitude function. In the general case the rate of phase propagation per period, a 'period phase constant' α, can be defined as a = > p a a Re 0.
The quantity α corresponds to the period wave number (k z ) of wave propagation. Note that the amplitude-phase relation (7) also suggests a real-valued local wave number (or angular frequency, when z is a time variable) p′ p (z) for bounded solutions. The period phase constant α affects both solutions F (±) (z), propagating their phases in different z-directions. In solid-state physics one often discusses the energy dependence α(E) by an inverse function E(α), a so-called dispersion relation; see right subplot in figure 3.
In case of unbounded solutions, one can identify the imaginary part p a a Im p ( ) as the so-called maximal Lyaponov exponent for the period map, i.e. l = p a a Im . 28 The quantity λ is non-zero only within the energy gaps; see left subplot of figure 3. For the potential model V(z) in (4) with v 0 =8 in figure 3, bounded solutions are possible in narrow regions of low energies. This is seen in the left subplot, where λ=0.
The phase quantity α (solid curve) in the middle subplot is constant for energies within the gaps, where the Lyaponov exponent of the left subplot is non-vanishing. This behavior causes characteristic energy jumps in the piecewise curve of the right subplot. The quantity α increases smoothly between its constant values separated by units of a (=π) in the middle subplot. α-increments are steep at low energies, but not infinitely steep. The phase p(a) defined by the non-periodic amplitude function divided by a (=π) is always real valued and appears as a broken curve in the middle subplot. The lines of α and p(a)/a cross or touch at distinct energies where their values are in the middle between two integers, consistent with equation (26). At integer values the curves touch on one side of the corresponding energy gap.
The two energies defining each gap allow periodic solutions; the solution z M 11 ( ) in (15) on one side, and M 12 (z) on the other side (see figure 1). The gaps correspond to α with integer values j=0, 1, 2, L in (25). Each value of j has an even and an odd periodic solution with respect to z=0, M 11 (z) respectively M 12 (z), which will be denoted '(+)' respectively '(−)' in the tables. The gap referred to in figure 2 corresponds to j=1 in (25).
It is clarified below that the odd fundamental solution M 12 (z) is periodic on the side where the curves of α and p(a)/ a touch.

Primary resonance energies defining the gap regions.
The generalization of p p (a) allows a determination of primary resonance energies defining the gaps. The energy gaps are related to the existens of periodic amplitude functions A p (z). In fact, the knowledge of the generalized p p (a) from (26) implies a value A p (0) from A A a p a p a 0 0 sin sin must be 0 . 29 One can see that A p (0) has a singular behavior as function of E. Branch points at = p a sin 0 p ( ) and = p a sin 0 ( ) seen in (29) imply that continuous ranges of E do not allow realvalued, positive amplitude functions.
An alternative way to express A p (0), without the use of an intermediate phase p p (a), is to use the formula A a A a p a 0 0 sin sin cos must be 0 .
This formula follows from the identification of the offdiagonal elements of (15) at z=a with those of (18), and by eliminating p a sin p ( ). Instead of a search algorithm to find periodic amplitude functions, one can from (30) calculate the value A p (0) directly with the basic (non-periodic) amplitude functions defined by (9) or (11).
It is seen in (30) that the periodic amplitude function has singular points whenever = -¢ = p a p a A a A a p a sin 0 and sin cos 0. 31 The energy limits of each gap (gap edges) can be computed from (26) or from the singularities of (30). Numerical results for the V(z) in (4) defined by v 0 = 8 are given in table 1. This potential is thus significantly stronger than that corresponding to figures 2 and 4, but is the one corresponding to figure 3. Table 1 shows energies consistent with those published by Connor et al [6]. The energies corresponding to a given gap are numbered by j=0, 1, 2, L followed by '±', indicating the existing even and odd periodic solutions M 11 (z) respectively M 12 (z). The For +1 these energies correspond to periodic solutions of period-a of either odd or even symmetry with respect to z=0. For −1 these energies correspond to periodic solutions of period-2a of either odd or even symmetry with respect to z=0. The broken curves corresponding to (33) respectively (34) allow supplementary conditions at zero crossings for close-lying energy roots of (32). Table 1. Primary resonance energies defining the gaps in figure 3. The potential (4) with v 0 =8 is used. The notation 'j±' indicates the gap number j=0, 1, 2, L and the corresponding even/odd symmetry of the existing periodic solution, M 11 (z) (+) respectively M 12 (z) (−). The entry '0−' corresponds to an energy E=−8, which is the energy limit below which solutions cannot have oscillations. A numerical tolerance 'tol' of about 10 −12 is used in a MatLab © environment. The present notation differs from that used in reference [6], and in references specifically using the Mathieu equation notation.

Application to weakly singular periodic potentials
Periodic potentials with singular points of attraction need special care. The singularity cannot be too strong without essential modifications of the above formulas. For the potential (5) the origin is chosen at a potential maximum. For energies lower than the potential maximum one cannot use the adiabatic expression (10). Instead, one can use the boundary conditions in equation (9). Note that the value of A (0) is not essential for the validity of the fundamental matrix solution (15) (only for the numerical efficiency), as long as it is positive. Note also that periodic amplitude functions are not needed in computations of the quantities α, λ.
Potential (5) results in energy behaviors of α and λ illustrated in figure 5. Such plots are useful guides for precise computations of energy gaps such as those given in table 2. A notable difference compared to entries in table 1 is the order of±-states for j>0.

Conclusions
The present application of the amplitude phase method to Schrödinger-like equations with periodic potentials provides consistent results with existing publications. The amplitude function is a key quantity of the method, which in principal is an arbitrary function of a nonlinear (Milne-type) equation. The existence of positive periodic amplitude functions of the present method is shown being equivalent to the existence of bounded solutions. Such amplitudes allow definitions of a real periodic wave number (or angular frequency) function and a period phase constant α. These quantities specify method-independent aspects of the fundamental solutions within energy bands. Unphysical (unbounded) solutions belong to so-called energy gaps. Each gap is specified by two energy limits satisfying alternative conditions expressed in terms of phase and amplitude values. These conditions differ for gap edges allowing even respectively odd periodic solutions with respect to the origin.