Perturbation of nonlinear operators in the theory of nonlinear multifrequency electromagnetic wave propagation

Commun Nonlinear Sci Numer Simulat 75 (2019) 76–93

Contents lists available at ScienceDirect

Commun Nonlinear Sci Numer Simulat journal homepage: www.elsevier.com/locate/cnsns

Research paper

Perturbation of nonlinear operators in the theory of nonlinear multifrequency electromagnetic wave propagation S.V. Tikhov, D.V. Valovik∗ Department of Mathematics and Supercomputing, Penza State University, 40, Krasnaya street, Penza, Russia

a r t i c l e

i n f o

Article history: Received 24 February 2019 Accepted 18 March 2019 Available online 23 March 2019 MSC: 47H14 35P30 35Q61 Keywords: Maxwell’S equations Nonlinear multiparameter eigenvalue problem Nonlinear perturbation Nonpertubative solutions

a b s t r a c t The paper develops an original approach to study nonlinear multiparameter eigenvalue problems arising in the theory of nonlinear multifrequency electromagnetic wave propagation. The problem under consideration is a multiparameter eigenvalue problem that under some conditions degenerates into n nonlinear one-parameter eigenvalue problems. Further simplification reduces the one-parameter nonlinear problems to linear (one-parameter) eigenvalue problems. Each of the linear problems has a finite number of positive eigenvalues, whereas each of the nonlinear (one-parameter) problems has an infinite number of positive eigenvalues. Using the nonlinear one-parameter problems as ’nonperturbed’ ones, one can prove existence of eigentuples of the multiparameter problem that have no connections with solutions to the linear (one-parameter) problems even if the nonlinear terms have small factors. © 2019 Elsevier B.V. All rights reserved.

1. Statement of the main problem and introductory remarks Recently there was introduced a large class of nonlinear multifrequency electromagnetic wave propagation phenomena described by nonlinear multiparameter eigenvalue problems (NMEPs) [1]. In spite of the fact that different types of NMEPs are known long ago and there are several approaches to study them, e.g., variational methods [2–4], bifurcation theory methods [5–8], perturbations of linear problems [9], etc., the problems introduced in [1] cannot be studied with the above mentioned methods. In fact, only a perturbation approach based on linear ’nonperturbed’ problems can be applied [10,11]. Indeed, NMEPs come from electromagnetics, on the one hand, have only local boundary conditions, on the other hand, they have linearised as well as nonlinearised solutions. The former property does not allow one to apply variational methods, the latter property essentially limits the use of perturbation methods. This paper focuses on a particular NMEP described propagation of a sum of transverse-electric (TE) guided waves [1]. To study the problem we will develop the integral dispersion equation method (IDEM), which is completely different from the methods mentioned in the previous paragraph [12]. Let us briefly discuss perturbation approaches. Classical perturbation theory methods are based on the following idea [13]. Let R(λ, α) be a nonlinear problem, where λ = (λ1 , . . . , λn ) and α = (α1 , . . . , αl ) are n-tuple and l-tuple of parameters, respectively. Such tuples can be considered as vectors. Let R be an operator-function corresponding to the problem R(λ, α). Thus the problem R(λ, α) can be written in an operator form R(u; λ, α ) = 0, where R depends on u nonlinearly. Let ∗

Corresponding author. E-mail addresses: [email protected] (S.V. Tikhov), [email protected] (D.V. Valovik).

https://doi.org/10.1016/j.cnsns.2019.03.020 1007-5704/© 2019 Elsevier B.V. All rights reserved.

S.V. Tikhov and D.V. Valovik / Commun Nonlinear Sci Numer Simulat 75 (2019) 76–93

77

there be a solution u ≡ u(x; λ) to the problem R(λ, 0), where 0 = (0, . . . , 0 ). If the problem R(λ, 0) is linear with respect to u, then under wide assumptions on R one can prove that the problem R(λ, α) also has a solution u ≡ u(x; λ, α) as soon as |α| is sufficiently small and u(x; λ, α) → u(x; λ) as |α| → 0. In the case of linear dependence of the problem R(λ, 0) on u, there are powerful methods to study its solvability [14,15]. For the case of two-parameter eigenvalue problems, where the nonperturbed problems are linear, such an approach was used, for example, in [10,11]. Nonlinear multiparameter problems, however, admit a more complicated but more fruitful approach. Let the problem R(λ, α ), where α = (α1 , . . . , αl  , 0, . . . , 0 ) is another l-tuple with α1 , . . . , αl  = 0 (l < l), be nonlinear but (in some sense) simpler than the problem R(λ, α). We assume that solvability of the problem R(λ, α ) can be established and let u ≡ u(x; λ, α ) be its solution. Then under certain conditions imposed on R, one can prove that the problem R(λ, α) is also solvable as soon as |α − α | is sufficiently small and u(x; λ, α) → u(x; λ, α ) as |α − α | → 0. In this case |α | is not necessarily small. The nonlinear approach discussed above can be helpful if: (i) there is an l-tuple α such that the problem R(λ, α ), being a simpler nonlinear problem, can effectively be studied; (ii) solutions to the problem R(λ, α ) exist not only for ’small’ |α |; (iii) the problem R(λ, α ) has solutions without linear counterparts. Clearly, using the problem R(λ, 0) as ’nonperturbed,’ one can get very scanty information about solutions to the problem R(λ, α) in the cases (ii), (iii). Moreover, the problem R(λ, 0) can have no solutions. Theory of nonlinear electromagnetic wave propagation is a rich source of problems satisfying conditions (i)–(iii), see, for example, [16–21] for scalar α and [1] for multidimensional α. The paper is organised as follows. Below in this section, statements of the multiparapameter nonlinear problem and auxiliary one-parameter nonlinear and linear problems are given; formulation of the results are given in Section 2; numerical simulations are presented in Sections 3, 4 provides a physical background for applications of the found results and shows where the studied problems come from; proofs are given in Section 5. First of all, introduce the notation we use. Everywhere below integer indexes i, j vary from 1 to n ≥ 2 and often we do not indicate this explicitly. Let ε l,i , ε s,i , ε c,i , Ai be 4n positive constants such that the inequalities ε l,i > ε s,i ≥ ε c,i > 0 hold, α ij be n2 nonnegative constants, and λi be n positive parameters. In addition, we use n-tuple λ and (n × n)-tuples α, α , and 0. The tuple λ consists of n parameters λi and can be considered as an n-dimensional vector λ = (λ1 , . . . , λn ). The tuple α consists of n2 parameters α ij and can be considered as a (n × n)-matrix α = (αi j ). The tuple α consists of n2 parameters α ij , where α ii > 0 and αi j = 0 for i = j; the tuple 0 can be considered as a zero (n × n)-matrix. We define sets i = (εs,i , +∞ ), ∗i = (εs,i , λ∗i ), where λ∗i are positive sufficiently big constants. The choice of λ∗i will be clear from Corollary 1, Theorems 4 and 6. It is assumed that αi j ∈ A∗i j , where A∗i j = (0, αi∗j ). In this notation αii∗ are arbitrary but fixed positive constants and αi∗j for i = j are positive constants that depend on αii∗ and λ∗i . The choice of αi∗j for i = j will be clear from statement 4, Theorems 4, 6.  The notation l Cl as well as C1 × . . . × Ck are used to define a (finite) Cartesian product of sets Cl . Below we use the   ∗ following Cartesian products  = i ∗i , A∗ = i, j A∗i j . The notation λ ∈ ∗ and α ∈ A∗ mean that λi ∈ ∗i and αi j ∈ A∗i j , respectively. As usually, R+ = (0, +∞ ). The (real) interval (0, h) and segment [0, h] are denoted by I and I¯, respectively. Now let us consider the system of n coupled equations

⎧    2 2 ⎪ ⎨ u1 = −(εl,1 − λ1 )u1 − α11 u1 + . . . + α1n un u1 , ⎪ ⎩

· · · · · · · · · · · · · · · · · · · · · · · · · · · · · · · · · · · · · · · · · · · · · · · · · · · · · ·,

un = −(εl,n − λn )un −



(1)

 αn1 u21 + . . . + αnn u2n un ,

∗

where the prime marks denote differentiation with respect to x; here it is assumed that (x, λ, α ) ∈ I¯ ×  × A∗ . Solutions to (1) are denoted by ui , ui (x), or ui (x; λ, α). The problem P(α) consists in finding n-tuples λ for which there exist solutions

u 1 ≡ u 1 ( x ; λ, α ) , . . . , u n ≡ u n ( x ; λ, α )

(2)

to system (1) that satisfy boundary conditions

ui ( 0 ) = Ai = 0,

κs,i ui (0 ) − ui (0 ) = 0,

κc,i ui (h ) + ui (h ) = 0, where κs,i =



λi − εs,i > 0, κc,i =

u1 , . . . , un ∈ C (I¯ ). 2



(3) (4)

λi − εc,i > 0, and such that (5)

The quantities κ s,i , κ c,i are assumed to be positive due to special condition imposed to solutions to the Maxwell equations, see Section 4. Thus λi > ε s,i .

78

S.V. Tikhov and D.V. Valovik / Commun Nonlinear Sci Numer Simulat 75 (2019) 76–93

Definition 1. If λ ∈ ∗ solves the problem P(α), then it is called an eigentuple and the corresponding solution (2) are called eigenfunctions of the problem P(α).  2 2 Remark 1. It is also important to study system (1) for (x, λ, α ) ∈ I¯ ×  × Rn+ , where  = i i and Rn+ is a Cartesian 2 product of n copies of R+ . However, this is also a mathematically challenging problem, see a discussion after Theorem 6. If α → α for i = j, then the problem P(α) degenerates into P(α ), which consists of n independent nonlinear problems Pi satisfying conditions (ii), (iii). In order to formulate the problems Pi rigorously let us consider the equation

vi = −(εl,i − λi )vi − αii v3i ,

(6)

where the prime marks denote differentiation with respect to x; here it is assumed that (x, λi , αii ) ∈ I¯ × i × R+ . Solutions to (6) are denoted by vi , vi (x), or vi (x; λi , α ii ). Every problem Pi consists in finding values λi for which there exist solutions

vi ≡ vi (x; λi , αii )

(7)

to Eq. (6) that satisfy boundary conditions

vi (0 ) = Ai = 0, κs,i vi (0 ) − vi (0 ) = 0,

(8)

κc,i vi (h ) + vi (h ) = 0

(9)

and such that

vi ∈ C 2 (I¯ ).

(10)

Definition 2. If λi ∈ i solves the problem Pi , then it is called an eigenvalue and the corresponding solutions (7) are called eigenfunctions of the problem Pi . Setting αi j = 0 for all i, j in (1), one arrives at the linear problem P(0). The problem P(0) consists of n independent linear problems, which are denoted by P0,i . To be more precise, every problem P0,i consists in finding values of λi for which there exist nontrivial twice-differentiable solutions wi ≡ wi (x; λi ) to equation wi = −(εl,i − λi )wi satisfying boundary conditions

κs,i wi (0; λi ) − wi (0; λi ) = 0, κc,i wi (h; λi ) + wi (h; λi ) = 0. Remark 2. The eigenvalues of the problems Pi and P0,i are squared propagation constants (PCs) of the corresponding waveguide systems. Since the PCs are real, then the eigenvalues are positive. Entries of the eigentuples of the problem P(α) are also squares of real numbers and, therefore, the eigentuples consist of positive entries, see Section 4 for details. Well, the problems P0,i are solved elementary. Using the problems P0,i as nonperturbed ones and applying, for example, inversion of the linear parts of the differential operators defined by the problems Pi and P(α), one can prove existence of eigenvalues and eigentuples of the problems Pi and P(α), respectively, under conditions that all α ij are sufficiently small, see, for example, [11,18]. Such an approach looks as evident as boring. Below we show that every problem Pi has infinitely many positive solutions for any α ii ( > 0), where infinitely many of the positive solutions do not have linear counterparts even if αii → +0. At the same time every problem P0,i has only a finite number of positive eigenvalues. Taking this into account and bearing in mind properties (ii) and (iii) given in the beginning of this section, one can see that perturbation methods based on using linear (nonperturbed) problems have very limited applicability to the problems Pi and P(α). In this paper, we are going to develop a more complicated but more delicate approach. Using the problems Pi as nonperturbed, one can prove existence of solutions to the problem P(α) that are close to solutions to the problems Pi . Going this way, one also proves existence of those solutions to the problem P(α), which have no connections with solutions to the problem P0,i even if αi j → +0 for all i, j. Smallness of α ii is not assumed in this case. In contrast to [2,22–25], we consider a multiparameter problem for which the number of spectral parameters coincides with the number of sought-for functions. It is worth saying that linear and nonlinear multiparameter eigenvalue problems attract considerable attention, see [26–30] and [5,6,9,31,32] for the linear and nonlinear cases, respectively. 2. Results Below we use additional notation for the eigentuples and eigenvalues. Eigentuples λ of the problem P(α) are denoted by λ¯ k1 ...kn = (λ¯ 1,k1 , . . . , λ¯ n,kn ), where k1 , . . . , kn are nonnegative integer indexes. Eigenvalues λi of the problems Pi and P0,i are denoted by

λi,k and  λi,k , respectively, where ki are nonnegative integer indexes. It is assumed that λ¯ i,ki ,

λi,k , and  λi,k are i

i

i

i

arranged in the ascending order. In some cases we also use the above notation without subindexes ki , ki . Below we equivalently reduce the eigenvalue problems to a transcendental equation or system of such equations w.r.t. spectral parameter(s) called the dispersion equation (DE) [12]; the same term is used in the multiparameter case.

S.V. Tikhov and D.V. Valovik / Commun Nonlinear Sci Numer Simulat 75 (2019) 76–93

79

2.1. Problems P0,i Due to its simplicity we immediately start to the following Statement 1. There exists a constant hmin ≥ 0 such that for any h > hmin the problem P0,i has a finite number of simple positive eigenvalues  λi,0 , . . . ,  λi,k−1 and  λi,k ∈ (εs,i , εl,i ). If ε l,i ≤ ε s,i , then the problem P0,i does not have (positive) solutions. Statement 1 results from the formula

tan(κl,i h ) =



κl,i (κs,i + κc,i ) , κl,i2 − κs,i κc,i

(11)

where κl,i = εl,i − λi > 0. Exact value for hmin can also be extracted from (11). Eq. (11) is well-known in (linear) optics, see, for example, [33]. One can use the problems P0,i as nonperturbed. In this case, assuming that all α ij are sufficiently small, it is possible to prove existence of eigenvalues

λi,k that are close to  λi,k (for finite k), see, for example, [18]. In the same way one can prove   ¯ ¯ ¯ existence of eigentuples λk1 ,...,kn = (λ 1,k1 , . . . , λn,kn ) that are close to n-tuples (λ1,k1 , . . . , λn,kn ) (for finite k1 , . . . , kn ), where  λ are solutions to the problem P , see, for example, [1,11]. It is clear that if ε ≤ ε , then perturbation methods based 0,i

i,ki

l,i

on the problems P0,i are not applicable to find solutions of the problems Pi and P(α).

s,i

2.2. Problems Pi Under condition (8) the first integral of (6) has the form

vi2 = Ci − (εl,i − λi )v2i − 12 αii v4i , 

(12)



+ εl,i − εs,i > 0. Let vi ≡ vi (x; λi , α ii ) be a solution to the Cauchy problem (6), (8) defined globally on I¯. This assumption, proved below, allows us to derive preparatory results that shed light on the idea of the method developed in this paper. Well, let us consider functions where Ci =

A2i 12 αii A2i

θi = v2i , μi = vi /vi ,

(13)

where vi is the above-mentioned solution to the Cauchy problem. Function θ i ≡ θ i (x) is defined and continuous for x ∈ I¯; function μi ≡ μi (x) is defined and continuous everywhere in I¯ where vi does not vanish. Functions θ i and μi are depend on λi , α ii , however, often we will not show this explicitly in order to shorten the notation. Obviously, vi and vi cannot take zero values at the same point. By virtue of (6), functions θ i (x) and μi (x) satisfy the following system

θi = 2θi μi ,

(14)

μi = −(μ2i + εl,i − λi + αii θi ). Taking into account (8) and (12), the first integral of system (14) takes the form 1 2

αii θi2 + (μ2i + εl,i − λi )θi = Ci .

(15)

It is easy to check that μi < 0 for all λi ∈ i . Indeed, assuming that μi vanishes, equating the right-hand side of the second equation in (14) to zero, and using (15), one comes to the equality αii θi2 = −2Ci . Since α ii , Ci > 0, then this equality cannot be true. The contradiction found proves the validity of the initial claim. Taking into account (15), one rewrites the second equation of (14) in the form

μi = −χi (μi ; λi ),

(16)

μ2i

where χi (μi ; λi ) = + εl,i − λi + αii θi > 0 and θ i ≡ θ i (μi ; λi ) is found from (15). Expressing θ i from (15), taking into account that θi = v2i  0, and then substituting the found expression into the second equation in (14), one arrives at the



formula χi (μi ; λi ) = (μ2i + εl,i − λi )2 + 2αiiCi

1/2

> 0, where min χi (μi ; λi ) =



2αiiCi .

Let vi have mi  0 zeros xi,r ∈ I, where r = 0, mi . If mi = 0, then vi does not vanish for x ∈ I¯. We remind that if vi ≡ 0, then vi (xi,r ) = 0. We use the index r (not ri ) in order to shorten the notation; it is not assumed that vi and vj for i = j necessarily have the same number of zeros. Let us define the intervals

Ii,0 = (0, xi,1 ),

Ii,r = (xi,r , xi,r+1 ),

where r = 1, mi − 1; if mi = 0, then Ii,0 = I.

Ii,m = (xi,m , h ), i

i

(17)

80

S.V. Tikhov and D.V. Valovik / Commun Nonlinear Sci Numer Simulat 75 (2019) 76–93

It follows from conditions (8)–(9) and inequality μi < 0 that

lim

θi (x ) = A2i ,

lim

θi ( x ) = 0 ,

x→+0

x→xi,r

lim

x→h−0

θi (x ) = Bi2 ,

lim

x→+0

μi (x ) = κs,i ,

(18)

lim

μi (x ) = ±∞,

lim

μi (x ) = −κc,i ,

x→xi,r ±0

x→h−0

r = 1, mi

(19)

(20)

where Bi are real unknown constants. Thus functions θ i (x) and μi (x) can be defined on each interval Ii,r as solutions to Eq. (14) with conditions appropriately chosen from (18)–(19). Integrating Eq. (16) on every interval Ii,r and using (18)–(19), one gets Statement 2. The Cauchy problem (6), (8) is globally unique solvable for x ∈ I¯ and its (classical) solution vi ≡ vi (x; λi , α ii ) continuously depends on (x, λi , αii ) ∈ I¯ × i × R+ and satisfies (10). In addition, the following formula

κs,i −∞

ds + (mi − 1 ) χi (s; λi )



+∞ ds + χi (s; λi ) μi ( h )

+∞

−∞

ds

χi (s; λi )

=h

(21)

is valid, where, in general, μi (h) is not prescribed. The DE results from formula (21) under condition (20).

i  0 such Theorem 1 (of equivalence). The value

λi is a solution to the problem Pi if and only if there exists an integer mi = m that λi =

λi is a solution to DE

i (λi ; mi ) ≡

κs,i −∞

ds + (mi − 1 ) χi (s; λi )

+∞

−∞

ds + χi (s; λi )

+∞ −κc,i

ds

χi (s; λi )

=h

(22)

i ; the corresponding eigenfunction vi ≡ vi (x;

i (simple) zeros for mi = m λi , αii ) has m

xi,r =

κs,i −∞

ds + (r − 1 ) χi (s;

λi )

+∞

−∞

ds

χi (s;

λi )

i . for r = 1, m

The equivalency between the problem Pi and DE (22) motivates the following

λi is a root of DE (22) of multiplicity pi , then

λi is a multiple eigenvalue of the problem Pi of the same Definition 3. If λi =

multiplicity. Below we use the notation

Ti (λi ) =

+∞ −∞

ds = χi (s; λi )

+∞ −∞

ds , s2 + εl,i − λi + αii θi (s )

(23)

where θ i (s) is defined from (15) with μi = s and χ i (s; λi ) is defined in (16). Since the function i (λi ; mi ) is estimated as

m Ti (λi ) < i (λi ; mi ) < (m + 1 )Ti (λi ),

(24)

then the properties Ti (λi ) w.r.t. λi play a crucial role. This is the reason to formulate Statement 3. Function Ti is positive and depends continuously on (λi , αii ) ∈ i × R+ ; for big λi it is true that − 12

Ti (λi ) = 2λi

− 12

ln λi + O(λi

).

(25)

Solvability of the problem Pi is established below. Theorem 2. There exists an integer mi  0 such that for every integer m  mi Eq. (22) has at least one (positive) solution

λi =

λi,m , where

λi,m → +∞ as m → +∞, and, therefore, the problem Pi has infinitely many (positive) eigenvalues

λi,m with an accumulation point at infinity. Furthermore, (1) if the problem P0,i has p (positive) solutions  λi,0 <  λi,1 < . . . <  λi,p−1 , then there exists a constant αii > 0 such that for any (positive) αii = αii < αii it is true that

λ

i,m

∈ (εs,i , εl,i ) and

lim

λi,m =  λi,m for m = 0, p − 1,

αii →+0

where

λi,0 , . . . ,

λi,p−1 are first p solutions to the problem Pi with αii = αii ;

S.V. Tikhov and D.V. Valovik / Commun Nonlinear Sci Numer Simulat 75 (2019) 76–93

81

(2) if m ≥ p, then

λi,m has no linear counterpart and limαii →+0

λi,m = +∞; (3) for big λi and any fixed i > 0 the asymptotical inequality

(1 − i )λi (m ) <

λi,m < (1 + i )λi (m + 1 ) is valid, where λi (m ) = g−1 ( 2hm ), g−1 is the inversion of g(t ) = t −1/2 ln t; /2 (4) maxx∈(0,h ) |vi (x;

λi,m , αii )| = O(

λ1i,m ) as m → ∞. In the subsequent analysis we will use items 1–2 of Theorem  ∞ 2 . Well, every problem Pi has an unbounded sequence

λi,m  of eigenvalues. Using statement 3 and Theorem 2, one m=mi

arrives at



Corollary 1. There is a subsequence

λi,mr properties:

∞

r =ri



with limmr →∞

λi,mr = ∞ of the sequence

λi,m

∞

m=mi

that has the following

(1) for any mr , the value λi =

λi,mr satisfies (22) with mi = mr ; (2) for any mr there exists an interval Ui,mr = (

λi,mr − δi,mr ,

λi,mr + δi,mr ), where δi,mr > 0 is a constant, containing no other  ∞ elements of

λ and such that the function h − (λ ; m ) takes values of different signs at the opposite endpoints i,mr r =r  i

i

i

r

of the interval Ui,mr .

λi,mr of the subseRemark 3. Values λ∗i introduced in Section 1 are chosen such that every set (εl,i , λ∗i ) contains elements

 ∞ quence

λi,mr from Corollary 1 .  r =ri

2.3. Problem P(α) We begin this subsection with the study of global unique solvability (GUS) of the Cauchy problem (1), (3), which we need to derive the DEs for the problem P(α). The GUS will be proved under restrictions λ ∈ ∗ and sufficiently small α ij with i = j. In this case, α ii are not assumed to be small. Following [34], let us consider a normal system of equations

ui = fi (x, u1 , . . . , un , γ1 , . . . , γl ),

i = 1, n,

(26)

where the right-hand sides depend on the parameters γ1 , . . . , γl . System (26) can be written in the vectorial form

u = f(x, u, γ ).

(27)

We assume that the right-hand sides of system (26) are defined and continuous together with their partial derivatives ∂    ∂ u f j (x, u, γ ) in a domain  ⊂ R, where R is a (1 + n + l )-dimensional space of variables x, u1 , . . . , un , γ1 , . . . , γl . i

The following result is valid.

 and u = ϕ (x, γ ) be a solution to Eq. (27) satisfying Theorem 3 (Pontrjagin [34]). Let (x0 , u0 , γ 0 ) be a point of the domain  initial condition ϕ (x0 , γ ) = u0 . If the solution u = ϕ (x, γ 0 ) is defined globally for x ∈ I¯, then there exists a constant γ 0 > 0 such that the solution u = ϕ (x, γ ) is defined on the same segment I¯ as soon as |γ − γ 0 | < γ0 ; furthermore, the function ϕ(x, γ ) is continuous with respect to variables x, γ for x ∈ I¯ and |γ − γ 0 | < γ0 . In addition, ϕ(x, γ ) → ϕ(x, γ ∗ ) uniformly with respect to x ∈ I¯ as γ → γ ∗ and |γ − γ ∗ | < γ0 . Remark 4. Theorem 3 is given in [34] without the claim about uniformity of limit ϕ(x, γ ) → ϕ(x, γ ∗ ); this fact, however, results from the proof given in [34]. The following theorem results from statement 2 and Theorem 3. Theorem 4. Let αii∗ , λ∗i be fixed constants, then there exist constants αi∗j > 0 with i = j such that for all 0 < αi j < αi∗j the Cauchy problem (1), (3) is globally unique solvable for x ∈ I¯ and its (classical) solution

u = ( u1 , . . . , un ),

(28) ∗

where ui ≡ ui (x; λ, α), continuously depends on (x, λ, α ) ∈ I¯ ×  and satisfies (5). In addition, ui → vi , where vi ≡ vi (x; λi , α ii ) is a solution to the Cauchy problem (6), (8), uniformly with respect to x ∈ I¯ as α → α . × A∗

∗ Let (x, λ, α ) ∈ I¯ ×  × A∗ , where λ∗i , αii∗ > 0 are fixed and λ∗i can be sufficiently big, constants αi∗j with i = j are chosen ∗ ∗ by λi , αii in accordance with Theorem 4. We assume that functions ui defined by (28) have mi ≥ 0 zeros xi,r ∈ I each, where r = 0, mi . If mi = 0, then ui does not vanish for x ∈ I. Clearly, if ui ≡ 0, then ui (xi,r ) = 0. Here as well as in the one-parameter cases (problems Pi , see p. 7) we use the index r (not ri ) in order to shorten the notation; it is not assumed that ui and uj for i = j necessarily have the same number of zeros (and of course, the index r does not coincide with the index of the same name introduced for the function vi ).

82

S.V. Tikhov and D.V. Valovik / Commun Nonlinear Sci Numer Simulat 75 (2019) 76–93

Let us define the intervals

Ii,0 = (0, xi,1 ),

Ii,r = (xi,r , xi,r+1 ),

Ii,mi = (xi,mi , h ),

(29)

where r = 1, mi − 1; if mi = 0, then Ii,0 = I. We also introduce functions

τi = u2i , ηi = ui /ui ,

(30)

where ui are components of the solution (28). Obviously, functions τ i ≡ τ i (x) are defined and continuous for x ∈ I¯; every function ηi ≡ ηi (x) is defined and continuous for x ∈ Ii,r , r = 0, mi . We do not show explicitly dependence of τ i and ηi on λ, α in order to shorten the notation. By virtue of (1) functions τ i (x) and ηi (x) satisfy the system of 2n equations

 τi = 2 τi ηi ,    ηi = − ηi2 + εl,i − λi + nj=1 αi j τ j ;

(31)

the right-hand sides of these equations are defined in that domains where ηi exist. Using Theorem 4, one gets Statement 4. Let functions ui ≡ ui (x; λ, α) be defined by (28). Then there exist constants αi j = αi0j > 0 with i = j such that

ηi2 + εl,i − λi +

n 

αi j τ j 



2αiiCi − δi > 0

(32)

j=1

for all (x, λ, α) ∈ xi × ∗ × A0 , respectively, where the left-hand sides in (32) are defined in (31), xi = with

A0ii

=

A∗ii ,

A0i j

=

(0, αi0j )

mi

I , r=0 i,r

A0 =

 ij

A0i j

for i = j, and quantities δ i can be chosen such that |δ i | are smaller than any (fixed) constant δ > 0

as soon as |α − α | is sufficiently small; Ci is defined in (12).

From this step we assume that (positive) parameters αi∗j in A∗i j for i = j are chosen such that αi∗j  αi0j and quantities δ i are chosen such that |δ i | ≤ δ , where δ > 0 is a fixed sufficiently small constant. It follows from statement 4 that ηi < 0 for all (x, λ, α) ∈ xi × ∗ × A∗ . This means that ηi monotonically decreases from

κ s,i to −∞ for x ∈ Ii,0 , from +∞ to −∞ for x ∈ Ii,r , where r = 1, mi − 1, and from +∞ to −κc,i for x ∈ Ii,mi . In other words, there exist continuous one-to-one mappings:

gi,0 : (−∞, κs,i ) → Ii,0 : ηi → x, gi,r : (−∞, +∞ ) → Ii,r : ηi → x,

where r = 1, mi − 1,

gi,mi : (−κc,i , +∞ ) → Ii,mi : ηi → x for i = 1, n. Using these mappings, one defines continuous functions τ i (x), where x ≡ gi,r (ηi ), on every interval Ii,r , see Fig. 1.   We stress that every function τi gi,r (ηi ) is defined for all r = 0, mi . This consideration allows one to impart sense to the

Fig. 1. Function ηi (x; λ) and mapping gi,r .

S.V. Tikhov and D.V. Valovik / Commun Nonlinear Sci Numer Simulat 75 (2019) 76–93

expressions

T (λ; Ii,0 ) =

T (λ; Ii,r ) =

κs,i s2 + εl,i − λi +

−∞

+∞

s2 + εl,i − λi +

−∞

T (λ; Ii,mi ) =

ds n j=1

ds n

+∞ −κc,i

s2 + εl,i − λi +

j=1

83

 , αi j τ j gi,0 (s )

(33)

  , r = 1, mi − 1, αi j τ j gi,r (s )

(34)

ds n j=1

 . αi j τ j gi,mi (s )

(35)

Taking the computations made into account, one can see that analogous one-to-one mappings gi,r can be constructed with the help of the second equation in (14) and functions μi . The following result takes place. ∗

Theorem 5 (of equivalence). The n-tuple λ¯ ∈  is an eigentuple of the problem P(α) if and only if there exist integers mi = ¯ i  0 such that λ = λ¯ is a solution to DEs m mi 

T (λ; Ii,r ) = h

(36)

r=0

¯ i ; furthermore, every eigenfunction ui ≡ ui (x; λ¯ , α ) has m ¯ i (simple) zeros xi,r = for mi = m

r−1

k=0

¯ i. T (λ; Ii,k ) for r = 1, m

The equivalency between problem P(α) and DEs (36) motivates the following Definition 4. If a n-tuple λ = λ¯ is a root of DEs (36) of multiplicity p, then λ¯ is a multiple eigentuple of the problem P(α) of the same multiplicity. Using Theorem 4 and above introduced mappings gi,r , one arrives at Corollary 2. If λ ∈ ∗ and α → α , then mi 

T (λ; Ii,r ) →

r=0

κs,i −∞

ds + (mi − 1 )Ti (λi ) + χ (s; λi )

+∞

−κc,i

ds

χ (s; λi )

.

(37)

Now subtracting from the left- and right-hand sides of Eq. (36) the left-hand side of formula (22) for mi = mi , one arrives at the system of equations mi 

T (λ; Ii,r ) − i (λi ; mi ) = h − i (λi ; mi ).

(38)

r=0

The left-hand side of (38) is sufficiently small as soon as |α − α | is small. Zeros of the right-hand sides of (38) are eigenvalues of the problems Pi . It is known from Corollary 1 that every problem Pi has infinitely many eigenvlaues

λi,r such that for any

λi,r there exists an interval Ui,r containing

λi,r for which the right-hand side corresponding to

λi,r in (38) takes values of different signs on opposite endpoints of this interval. Taking into account this consideration, one arrives at the main Theorem 6. Let every problem Pi have mi eigenvalues

λi,1 , . . . ,

λi,mi ∈ ∗i ⊂ i defined in Corollary 1. Then there exist positive ∗ ∗ constants αi j for i = j such that for any 0 < αi j < αi j (i = j) the problem P(α) has at least m1 × . . . × mn eigentuples λ¯ k1 ,k2 ,...,kn = ¯ ¯ (λ¯ ,λ ,...,λ ), where ki = 1, mi ; furthermore, every λ¯ ∈U × ...× U , where

λ ∈U . 1,k1

2,k2

n,kn

k1 ,k2 ,...,kn

1,k1

n,kn

i,ki

i,ki

Corollary 3. If any of the eigenvalues

λ1,k1 , . . . ,

λn,kn has no linear counterpart, (

λ1,k1 , . . . ,

λn,kn ) ∈ ∗ , and

λ1,k1 , . . . ,

λn,kn are

¯ chosen in accordance with Corollary 1, then there exists at least one eigentuple λ in the vicinity of (λ1,k1 , . . . ,

λn,kn ) as soon as |α − α | is sufficiently small and this eigentuple has no linear counterpart. Remark 5. Constants λ∗i are chosen sufficiently big in order that every interval (εl,i , λ∗i ) contains eigenvalues

λi,k of the problem Pi that satisfy Corollary 1; constants αii∗ > 0 can be chosen without additional restrictions. Intervals (εl,i , λ∗i ) do not contain solutions to the problems P0,i . After fixing parameters λ∗i and αii∗ in accordance with Theorem 4 and statement 4, one fixes (in general, sufficiently small) parameters αi∗j for i = j. In view of Theorem 6, it is worth giving a few comments. In the first place, since values λ∗i in i can be chosen as ¯ ¯ big as necessary, then Theorem 2 states existence of eigentuples λ¯ k1 ,...,kn = (λ 1,k1 , . . . , λn,kn ) that, in particular, belong to the domain where there are no solutions to the problems P0,i . As far as we know, such a result in the multiparameter problem considered in this paper has not been known before. In the second place, Theorem 6 also states existence of eigentuples λ¯ k1 ,...,kn that can be considered as perturbations of n-tuples ( λ1,k1 , . . . ,  λn,kn ), where  λ1,k1 , . . . ,  λn,kn are solutions

84

S.V. Tikhov and D.V. Valovik / Commun Nonlinear Sci Numer Simulat 75 (2019) 76–93

to the problems P0,i . Similar results can be obtained using integral equation approach, see, for example, [1,10,11] and the bibliography there in. Further, only a finite number of eigentuples λ¯ k1 ,...,kn is stated in Theorem 6. Since every problem Pi has infinitely many of positive eigenvalues, then it is reasonable to assume that the problem P(α) also has infinitely many eigentuples with positive entries. Proving this result (at least for small αi∗j ) will be the next essential advance. The first and evident obstacle to 2 derive such a result is the absence of GUS of the Cauchy problem (1), (3) for (x, λ, α ) ∈ I¯ ×  × Rn+ or (x, λ, α ) ∈ I¯ ×  × A∗ .

The second obstacle is an unbounded growth of max of eigenfunctions vi as λi → +∞, see item 4 of Theorem 2. The boundedness of vi is essential for Corollary 2. Possibly, the latter obstacle can be got over using appropriate ’normalization’ of the eigenfunctions as it is done in the proof of statement 2.1 in [17]. As far as we know, DEs (36) derived in this paper is a new tool to study nonlinear eigenvalue problems. Indeed, the derivation of these equations can be made as soon as one can derive GUS of the corresponding Cauchy problem and monotonicity property for some functions related to the main problem (here such functions are ηi ). Moreover, in many cases if the monotonicity property is proved, then it will not be necessary to separately derive GUS of the corresponding Cauchy problem. For example, if in the problem P(α) one can prove that ηi < 0 for every (λ, α) ∈  × α , where  , α are continu ous sets, then one will get the required GUS for (x, λ, α ) ∈ I¯ ×  × α and derive the system of DEs, which is equivalent to

the original problem for (λ, α) ∈  × α . Probably, the main idea of the developed approach may be applicable for a wide range of nonlinear multiparameter problems. In addition, DEs (22) and (36) explicitly contain number of zeros of the corresponding eigenfunctions; exact localisation of the zeros are given in Theorems 1 and 5; distributions of the eigenvalue and eigentuples are also extracted from DEs (22) and (36), etc. 3. Numerical simulations

In this section we consider a two-parameter case. Such a case, on the one hand, is the first nontrivial multidimensional case; on the other hand, it can be illustrated sufficiently good. The following parameters are used without change in each simulation: n = 2, εl,1 = 3, εl,2 = 5, εc,1 = εs,1 = 1, εc,2 = εs,2 = 1.2, A1 = A2 = 1, h = 4. Parameters α ij are chosen differently for different cases. In spite of the fact that solutions to problems P0,1 , P0,2 and P1 , P2 are scalar eigenvalues, we will consider them as 2tuples  λk,k = ( λi,k ,  λ j,k ) and

λk,k = (

λi,k ,

λ j,k ), where i, j = {1, 2} and k, k are nonnegative integer indexes. Uniting the eigenvalues into the pairs, we plot them as points in the plane Oλ1 λ2 . We do not overload the plots with too many alphabetic notation (many of them are put into the captions). In each plot below, the 2-tuples  λk,k = ( λi,k ,  λ j,k ), which represent solutions to problems P0,1 , P0,2 , can exist only inside the rectangular with (green) dashed border. In Figs. 2 and 3 eigentuples  λk,k = ( λi,k ,  λ j,k ) of problems P0,1 , P0,2 are plotted with empty circles; eigentuples

λk,k =



(λ , λ  ) of problems P1 , P2 are plotted with (red) filled circles. In Fig. 3 eigentuples λ¯  = (λ¯ , λ¯  ) of problem P(α) i,k

j,k

k,k

i,k

j,k

are plotted with (blue) diagonal crosses. For the parameters chosen, there exist 2 and 3 eigenvalues in the problems P0,1 and P0,2 , respectively. This means that there are 6 eigentuples. To be more precise, in the linear case we get 6 eigetuples  λk,k = ( λi,k ,  λ j,k ), where  λ1,0 ≈ 1.754,  λ ≈ 2.669 and  λ ≈ 1.828,  λ ≈ 3.488,  λ ≈ 4.612. 1,1

2,0

2,1

2,2

In Fig. 2 eigenvalues of problems P1 , P2 are plotted as 2-tuples together with eigenvalues of problems P0,1 , P0,2 in order one can compare them. Problems P1 , P2 have an infinite number of positive eigenvalues each; here we present only few first eigenvalues in each case. One can see that subfigures (a) and (b) of Fig. 2 are in a good agreement with the item 1) of Theorem 2. Indeed, making α 11 and α 22 less, one can see that there exists at least one 2-tuple

λk,k = (

λi,k ,

λ j,k ) in the vicinity of each 2  

tuple λ  = (λ , λ  ) and, therefore, there exists at least one eigenvalue λ in the vicinity of each eigenvalue  λ as k,k

i,k

j,k

i,k

i,k

soon as α ii is sufficiently small. As above, one can see that subfigures (a) and (b) of Fig. 3 are in a good agreement with Theorem 6. Indeed, one can see that there exists at least one 2-tuple

λk,k = (

λi,k ,

λ j,k ) in the vicinity of each 2-tuple  λk,k = ( λi,k ,  λ j,k ) as soon as α 12 and α 21 are sufficiently small. 4. Applications In this section we discuss some physical applications of the found results. Our main interest is focused on a class of fundamental nonlinear optics problems that lead to the problems P0,i , Pi , P(α), and similar ones. However, some other applications will be pointed out in the end of this section. Let  = {(x, y, z ) ∈ R3 : 0  x  h, (y, z ) ∈ R2 } be a layer filled with nonlinear dielectric and located between two halfspaces x < 0 and x > h in the Cartesian coordinates Oxyz. The half-spaces are filled with nonmagnetic media characterised by real constant permittivities  = s  0 > 0 and  = c  0 > 0, respectively, where  0 is the permittivity of vacuum. The permittivity  l of the layer  is described below; the permeability μ in the whole space is a positive constant.

S.V. Tikhov and D.V. Valovik / Commun Nonlinear Sci Numer Simulat 75 (2019) 76–93

85

Fig. 2. In subfig. (a) α11 = 0.045, α22 = 0.055, and α12 = α21 = 0; here

λ1,0 ≈ 1.798,

λ1,1 ≈ 2.904,

λ1,2 ≈ 4.727 and

λ2,0 ≈ 1.851,

λ2,1 ≈ 3.581,

λ2,2 ≈ 5.433,

λ2,3 ≈ 5.738. In subfig. (b) α11 = 0.019, α22 = 0.022, and α12 = α21 = 0; here

λ1,0 ≈ 1.774,

λ1,1 ≈ 2.760,

λ1,2 ≈ 5.902 and

λ2,0 ≈ 1.843,

λ2,1 ≈ 3.527,

λ2,2 ≈

4.797, λ2,3 ≈ 7.407.

Fig. 3. In subfig. (a) α11 = 0.019, α22 = 0.022, α12 = 0.01, α21 = 0.02; here λ¯ 00 ≈ (1.778, 1.857 ), λ¯ 01 ≈ (1.789, 3.552 ), λ¯ 02 ≈ (1.815, 4.810 ), λ¯ 03 ≈ (1.041, 7.429 ), λ¯ 10 ≈ (2.766, 1.896 ), λ¯ 11 ≈ (2.773, 3.595 ), λ¯ 12 ≈ (2.858, 4.929 ), λ¯ 13 ≈ (2.046, 7.465 ), λ¯ 20 ≈ (5.912, 2.822 ), λ¯ 21 ≈ (5.947, 4.522 ). In subfig. (b) α11 = 0.019, α22 = 0.022, α12 = 0.001, α21 = 0.002; here λ¯ 00 ≈ (1.774, 1.843 ), λ¯ 01 ≈ (1.775, 3.530 ), λ¯ 02 ≈ (1.777, 4.799 ), λ¯ 03 ≈ (1.797, 7.409 ), λ¯ 10 ≈ (2.761, 1.848 ), λ¯ 11 ≈ (2.761, 3.534 ), λ¯ 12 ≈ (2.768, 4.809 ), λ¯ 20 ≈ (5.901, 1.977 ), λ¯ 21 ≈ (5.903, 3.602 ), λ¯ 22 ≈ (5.908, 5.176 ), λ¯ 23 ≈ (5.858, 7.190 ).

In accordance with [1], we introduce the multi-frequency field

Eω =

n  j=1

E j e−iω j t ,

Hω =

n 

H j e−iω j t ,

(39)

j=1

˜ω, H ˜ ω has the form where E j = E+j + iE−j , H j = H+j + iH−j are the complex amplitudes [35]. The real (physical) field E ˜ ω (x, y, z, t ) = Re Eω , H ˜ ω (x, y, z, t ) = Re Hω . Frequencies ωj are different but their choice submits to special terms related E to a particular nonlinear law chosen for  l [1,16,36].

86

S.V. Tikhov and D.V. Valovik / Commun Nonlinear Sci Numer Simulat 75 (2019) 76–93

We assume that the permittivity  l of the layer  is a diagonal (3 × 3)-tensor that depends on the field by the Kerr law, that is,



l (E˜ ω ) ≡

x + fx

0 y + fy 0

0 0



0 0 , z + fz

where  x ,  y ,  z are real positive constants such that min { x ,  y ,  z } >  s ≥  c > 0 and fr ≡

βy, j,r |(E j , ey

 2

)|2

(40) n  j=1

βx, j,r |(E j , ex )|2 +

+ βz, j,r |(E j , ez )| ; here βr, j,r are real constants, (·,·) is the euclidian scalar product, er is a unit vector in r-direction, r, r ∈ {x, y, z}. The permittivity in the form (40) is in agreement with important real situations [16,19,21,35,37–42]; in addition, it is sufficient to study various types of waves, for example, TE, TM, and, so called, coupled TE-TE and TE-TM waves in the Kerr case. Two-frequency guided waves (GWs) have been found recently: (coupled) TE-TM and TE-TE GWs in a layer with Kerr nonlinearity are introduced and studied in [10] and [11], respectively. The case of coupled TE-TM GWs in a layer with Kerr nonlinearity propagating at one frequency, however, has been introduced theoretically long ago in [37]; then it took attention later in [16,38,39]. The case of multifrequency GWs in nonlinear photonic crystals is reported in [41]. Quite general formulation of multiparameter problems for different field configuration in plane as well as circle cylindrical waveguides and different nonlinearities are presented in [1]. Substituting (39) into Maxwell’s equations and taking into account that the operator rot is linear, one derives that Ek , Hk satisfy the following (coupled) equations

⎧n n  −iω t  ⎪ j rot H = −i ω j E j e−iω j t , ⎪ j ⎨ e j=1

j=1

n n   ⎪ ⎪ ⎩ e−iω j t rot E j = iμ ω j H j e−iω j t , j=1

(41)

j=1

where i is the imaginary unit. Since the derived system must be fulfilled for all t, then one arrives at the following system of n (coupled) systems



rot H j = −iω j E j ,

(42)

rot E j = iμω j H j ,

where j = 1, n. Thus Ej , Hj satisfy Eq. (42) and decay as O(|x|−1 ) when |x| → ∞. Tangential components of the fields Ej , Hj are continuous at the interfaces x = 0, x = h. In addition, the fields Ej , Hj have prescribed values at the boundary x = 0. Now let us consider a particular configuration of the filed (39) that results in the problem studied in sections 1, 2. Let an integer index j be such that 1 ≤ j ≤ n. We consider the fields Ej , Hj to be of the form

E j = (0, ey( j ) , 0 ) eiγ j z ,

H j = (hx( j ) , 0, hz( j ) ) eiγ j z for 1  j  j ,

E j = (0, 0, ez( j ) ) eiγ j y ,

H j = (hx( j ) , hy( j ) , 0 ) eiγ j y for j  j  n,

( j)

( j)

( j)

( j)

(43)

( j)

where components ey , ez , hx , hy , hz depend on spatial variable x only (these quantities, as solutions to Maxwell’s equations, also depend on other parameters of the problem) and γ j are unknown real constants. In other words, we consider a sum of TE fields propagating in directions Oz and Oy, respectively. ( j) ( j) Substituting (43) into (42) and using the notation u j := ey for j = 1, j and u j := ez for j = j + 1, n, after some algebra one arrives at the following system

⎧  u1 = −(μω12  − λ1 )u1 , ⎪ ⎪ ⎪ ⎪ ····················· , ⎪ ⎪ ⎪ ⎪ ⎨ u = −(μω2  − λ j )u j , j j

⎪ uj +1 = −(μω2j +1  − λ j +1 )u j +1 , ⎪ ⎪ ⎪ ⎪ ⎪ ····················· , ⎪ ⎪ ⎩  un = −(μωn2  − λn )un ,

(44)

where λ j = μω2j γ j2 > 0 and

⎧ ⎨ s , x < 0,  = l , 0  x  h, ⎩ c , x > h.

(45)

S.V. Tikhov and D.V. Valovik / Commun Nonlinear Sci Numer Simulat 75 (2019) 76–93

87

In accordance with the condition at infinity, solutions to system (44) in the half-spaces have the form

u j (x ) = A j eκs, j x , u j (x ) = B j e

−κc, j (x−h )

for x < 0, ,

for x > h;

here Aj , Bj are real nonzero values and (without loss of generality) Aj is positive, κs, j = where εs, j = μω2j s , εc, j = μω2j c . Inside the layer  system (44) takes the form

  ⎧  u1 = −(1,1 − λ1 )u1 − β1,1 u21 + . . . + β1,n u2n u1 , ⎪ ⎪ ⎪ ⎪ ⎪ ·········································· , ⎪ ⎪ ⎪ ⎨ u = −(1, j − λ j )u j − β1,1 u2 + . . . + β1,n u2 u j , n 1 j    2 ⎪    u = − (  − λ ) u − u + . . . + β2,n u2n u j +1 , β  +1 2,1 1 2, j +1 j +1 j +1 ⎪ j ⎪ ⎪ ⎪ ⎪ ·········································· , ⎪ ⎪   ⎩  un = −(2,n − λn )un − β2,1 u21 + . . . + β2,n u2n un ,

(46)



λ j − εs, j > 0, κc, j =



λ j − εc, j > 0,

(47)

where 1, j = y μω j for j = 1, j , 2, j = z μω j for j = j + 1, n, β1, j = μω2j βy, j,y for j = 1, j , β1, j = μω2j βz, j,y for j = j + 1, n,

β2, j = μω2j βy, j,z for j = 1, j , β2, j = μω2j βz, j,z for j = j + 1, n.

Clearly, system (1) is more general than system (47) due to the fact that the coefficients in front of the nonlinear terms in (1) are all different. ( j) ( j) ( j) ( j) Tangential components of the studied field are ey and hz for j = 1, j and ez and hy for j = j + 1, n. Thus the functions uj and uj are continuous at x = 0 and x = h. Taking this into account and using solutions (46), one finds

ui |x=0 = Ai ,

κs,i ui |x=0 − ui |x=0 = 0,

−κc,i ui |x=h − ui |x=h = 0.

(48) (49)

The conditions listed in this paragraph result in conditions (3), (4) if n ≥ 2 and conditions (8), (9) if n = 1. The field (43) propagates in the layer  only for special values of γ j . These values are called propagation constants (PCs). From the mathematical standpoint, the above formulated problem is a nonlinear multiparameter eigenvalue problem for system (47) with the above listed boundary conditions. The eigentuples (or eigenvalues in the one-parameter case) are vectorial (or scalar) PCs. Since γ j in (43) are real, then μω2j γ j2 are positive and for this reason electromagnetic applications require only positive

λj in the problems P(α) and Pj . If n = 1 and, therefore, j = 0 or j = 1, then one comes to one of the problems Pj . If all β 1,j and β 2,j are zeros, then one arrives at n linear problems that arise when one needs to determine linear guided TE waves propagating in the layer  with linear permittivity. These linear problems are equivalent to the problems P0,j formulated in Section 1. Nonlinear laws that are used in the nonlinear optics of waveguides have small factors; these factors are usually small parameters (this is true for the Kerr nonlinearity) [16,20,36]. This allows one to apply perturbation methods based on linear problems and prove existence of solutions to the nonlinear problem that are close to solutions of the used linear problems (see, for example, [10,11]). As is well known, in linear eigenvalue problems eigenfunctions are determined up to a constant factor; the eigenvalues are uniquely determined [43]. For nonlinear eigenvalue problems (when equations depend nonlinearly on the searched for functions), the same boundary conditions, which are enough in the linear case, do not allow one to determine the eigenvalues uniquely. At the same time, linear electromagnetic wave propagation problems in a plane layer have discrete sets of PCs, see [33]. If one generalises a linear problem to the nonlinear situation, then it is natural to formulate the nonlinear problem in such a way that solutions to the nonlinear problem have linear counterparts at least for ’small’ nonlinearities. Thus, the necessity of an additional conditions is clear. As an additional condition, one can fix (or prescribe) value of the field components (or their derivatives) at one of the boundaries, for example at x = 0. Fixing norms (in an appropriate func( j) ( j) tion space) of ey and ez , one gets another variant of the additional condition. We should stress however that in an open waveguide, a natural additional condition is the former one; moreover, the latter condition is not suitable from the physical point of view, see also [1]. In fact, results and methods presented in this paper can be helpful in that fields where nonlinear multiparameter phenomena arise. In addition to nonlinear optics, the other wide branch of nonlinear multiparameter problems is the coupled oscillator theory with nonlinear interaction, see, for example, [44–48] and the bibliography therein; experimental observations are given in [46,48]. The coupling functions between oscillators can be different from the Kerr law used in this paper but the main idea remains the same. One should reduce (if possible) the main nonlinear multiparameter problem to several one-parameter nonlinear problems. If these nonlinear problems can effectively be studied and their behaviour is different from the behaviour of the corresponding linear problems, then one has a fair chance to extend results found for the one-parameter (nonlinear) problems to the multiparameter case. The above mentioned references also contain additional information about applications.

88

S.V. Tikhov and D.V. Valovik / Commun Nonlinear Sci Numer Simulat 75 (2019) 76–93

5. Proofs v Proof of statement 2. Introducing variables v¯ i = √i , it follows from (12) that v¯ 2i is bounded for all λi ∈ i , v¯ i 2 is bounded λi

for finite λi and v¯ i 2 = 2α1 λi + O(1 ) for big λi . This boundedness and the subsequent analysis results in the required GUS. ii Well, it follows from formulas (14)–(15) as well as (16) that μi (x) monotonically decreases on every interval Ii,r , r = 0, mi , see formula (17). If mi = 0, then Ii,0 = I. The behaviour of μi at the endpoints of every interval Ii,r is known from formulas (18)–(20). Let

μ−i,r = μi (xi,r − 0 ) = −∞,

μ+i,0 = μi (0 + 0 ) = κs,i ,

(50)

μ−i,m +1 = μi (h − 0 ) = μi (h ), μ+i,r = μi (xi,r + 0 ) = +∞, i

where r = 1, mi , denote the boundary values of μi at the above mentioned endpoints. We intentionally write μi (h) not −κc,i in order to stress that the value μi (h) is not prescribed in this proof. Integrating Eq. (16) on every interval Ii,r , one gets

⎧ μ− i,1 ds ⎪ ⎪ = x + ci, 0 , ⎪ ⎪ μi (x ) χi (s; λi ) ⎪ ⎪ ⎪ μi ( x ) ⎨ ds  , = x + ci,r − + χ ( μi,r i s; λi ) ⎪ ⎪ ⎪

μi ( x ) ⎪ ⎪ ds  ⎪ ⎪ ⎩− μ+ χi (s; λi ) = x + ci,mi , i,m

x ∈ Ii,0 ; x ∈ Ii,r , r = 1, mi − 1;

(51)

x ∈ Ii,m . i

i

Substituting x = 0 + 0, x = xi,r+1 − 0, x = h − 0 into the first, second, and third lines, respectively, of (51), one derives ci, 0 ,

 , and c , where r = 1, m − 1. Substituting the found c into (51), then substituting x = x − 0, x = x + 0, x = x ci,r +0 i,r i,1 i,r i i,m i,m i

i

into the first, second, and third lines, respectively, of the formula found at the previous step, and replacing μ± in (51) in i,k accordance with (50), one arrives at the formulas

κs,i ⎧ ds ⎪ 0 < xi,1 = , ⎪ ⎪ χ ( −∞ ⎪ i s; λi ) ⎪

+∞ ⎨

⎪ ⎪ ⎪ ⎪ ⎪ ⎩

0 < xi,r+1 − xi,r = 0 < h − xi,m = i

−∞ +∞

μi ( h )

ds

χi (s; λi ) ds

χi (s; λi )

, r = 1, mi − 1;

(52)

.

In detail this procedure is presented in [12]. Formulas (52) express distances between adjacent zeros of vi . Since the left-hand sides in (52) are finite, then the righthand sides are also finite. Thus all the improper integrals in the right-hand sides converge. We note that the mentioned convergence can also be derived form the fact that the right-hand side of (23) converges. Summing up all the terms in (52), one finally arrives at

xi,1 + xi,2 − xi,1 + xi,3 − xi,2 + . . . + xi,m − xi,m −1 + h − xi,m = i

i

i

κs,i −∞

ds + (mi − 1 ) χi (s; λi )

+∞

−∞

+∞ ds + χi (s; λi ) μi ( h )

ds

χi (s; λi )

.

Formula (21) results from the previous formula. The computations made allows one to claim that solution vi ≡ vi (x; λi , α ii ) to the Cauchy problem (6), (8) globally exists for x ∈ I¯. Uniqueness of this solution and its continuous (and differentiability) on x ∈ I¯, λi ∈ i , and αii ∈ R+ result from smoothness of the right-hand side of Eq. (6) w.r.t. vi , λi , and α ii [34].  v ( h )

Proof of Theorem 1. It follows from (9) that μi (h ) = vi (h ) = −κc,i . The DE (22) results from (21) under condition (20). i It follows from the computation made in the proof of statement 2 that every solution (eigenvalue) to the problem Pi

i . satisfies DE (22) with certain mi = m Let us prove that every solution to DE (22) is an eigenvalue of the problem Pi . Let λi =

λi be a solution to (22) with

. We consider the Cauchy problem (6), (8) with λ =

m = m λ . Bearing in mind statement 2, one finds that solution v ≡ i

i

i

i

i

vi ( x;

λi , αii ) of this Cauchy problem globally exists and is unique for x ∈ I¯. Using solution vi , one can construct functions θi = v2i and μi = vi /vi . Clearly, θi (0 ) = A2i and μi (0 ) = κs,i . At this step we do not assume that μi (h ) = −κc,i . Setting, for certainty, μi (h ) = vi (h )/vi (h ) = a > −κc,i and using the functions θ i and μi , one can construct the expression

κs,i −∞

ds

i − 1 )Ti (

+ (m λi ) + χi (s;

λi )

+∞ a

ds

χi (s;

λi )

= h.

(53)

S.V. Tikhov and D.V. Valovik / Commun Nonlinear Sci Numer Simulat 75 (2019) 76–93

89

Since solution to the Cauchy problem (6), (8) is unique, then the integrand in (53) coincides with the analogous in (22).

i and subtracting (22) from (53), one finds Bearing in mind that λi =

λi satisfies (22) with mi = m

+∞ a

ds − χi (s;

λi )

+∞ −κc,i

ds =− χi (s;

λi )

a

ds

−κc,i

χi (s;

λi )

= 0.

(54)

By virtue of evident estimation

a

ds

−κc,i

χi (s;

λi )

> 0,

one gets that assumption μi (h ) = a > −κc,i is incorrect. Similarly, one can prove that it is not true that μi (h ) = a < −κc,i . Taking into account these two facts, one obtains that (54) is fulfilled if and only if a = −κc,i . This means that

λi is an eigenvalue. The computations made in the proof of statement 2 lead to the result about zeros xi,k of the eigenfunction vi , see formulas (52).  Proof of statement 3. It follows from (16) and (23) that Ti (λi ) is positive and depends continuously on (λi , αii ) ∈ i × R+ . √ −1 −1 −1 −1 The rough inequalities 2λi 2 ln λi + O(λi 2 )  Ti (λi )  2 2λi 2 ln λi + O(λi 2 ) result from the obvious estimations

Ti∗  Ti (λi )  where Ti∗ =

√ 2Ti∗ ,

(55)

 +∞

ds √ −∞ |s2 +ε −λ |+ 2α C i ii i l,i

above can be replaced with 1.

[21]. Using more delicate calculations, it can be shown that

√

2 in the right-hand side



Proof of Theorem 2. The proof of this theorem can be found in [21]. For this reason here we present only a sketch of the proof. It follows from (25) that limλi →+∞ Ti = 0. This means that there exists an integer mi  0 such that for any integer m  mi Eq. (22) has at least one positive solution λ =

λ . Thus there exists an infinite number of positive solutions. i

i,m

Let λi ∈ [εs,i + δ, εl,i − δ ], where δ > 0 is a sufficiently small number. In this case, one can pass to the limit αii → +0 in (22). The resulting equation takes the form

 ( λ ; m ) ≡ i i i

κs,i −∞

ds + (mi − 1 ) i (s; λi ) χ

+∞

−∞

ds + i (s; λi ) χ

+∞ −κc,i

ds i (s; λi ) χ

= h,

i (s; λi ) = s2 + εl,i − λi . Calculating the integrals in the left-hand side and taking tan , one arrives at the linear DE (11), where χ where one can take λi ∈ (ε s,i , ε l,i ). This leads to the link between solutions to the problems Pi and P0,i . Existence of solutions

λi without linear counterparts results from formula (55). Two sided asymptotical inequality results from (22), (24), and (25). If an eigenfunction vi has two or more zeros inside I¯, then it has a (local) extremum vi, ext inside I. Let this extremum be achieved at x = x∗ , that is, vi (x∗ ) = 0. Using (12) at this point, one gets a biquadratic equation with respect to vi, ext . Solving this equation and assuming that λi is sufficiently big, one comes to item 4.  Proof of Theorem 4. The claim of Theorem 4 is a direct consequence of Theorem 3 applied to the Cauchy problem (1), (3) under additional assumption about GUS of the Cauchy problems (6), (8) for (x, λi , αii ) ∈ I¯ × i × A. Validity of this additional assumption is guaranteed by statement 2. ∗ Continuity of the solution to the Cauchy problem (1), (3) w.r.t. (x, λ, α ) ∈ I¯ ×  × A∗ and its differentiability w.r.t. x ∈ I¯ result from smoothness of the right-hand side of system (1) with respect to ui , λi , and α ij [34].  Proof of statement 4. It follows from Theorem 4 that if λ ∈ ∗ and α → α , then

ui (x; λ, α ) → vi (x; λi , αii ),

ui (x; λ, α ) → vi (x; λi , αii )

(56)

uniformly with respect to x for x ∈ I¯, where ui (x; λi , α) is a solution to the Cauchy problem (1), (3) and vi (x; λi , α ii ) is a solution to the Cauchy problem (6), (8). Strictly speaking, it is necessary to prove separately that ui → vi uniformly with respect to x for x ∈ I¯ as α → α . However, this fact results from GUS of the Cauchy problem (1), (3) and smoothness of the right-hand sides of Eq. (1) with respect to ui . This fact can also be extracted from the proof of Theorem 3 [34]. It follows from formulas (56) that if λ ∈ ∗ and α → α , then

90

S.V. Tikhov and D.V. Valovik / Commun Nonlinear Sci Numer Simulat 75 (2019) 76–93

τi (x; λ, α ) → θi (x; λi , αii )

(57)

uniformly with respect to x for x ∈ I¯; here τ i (x; λ, α), θ i (x; λi , α ii ) are defined by formulas (30), (13), respectively. If x lies sufficiently close to any point xi,r , where r = 1, mi , then inequality (32) is obviously satisfied. Indeed, λi is bounded, other quantities in this formula are positive, and quantity ηi2 (x ) is sufficiently big. Thus it only remains to prove validity of inequalities (32) on some closed sets I¯i,r , where every I¯i,r lies strictly inside Ii,r . Well, if x ∈ I¯i,r , then ηi (x) is bounded and continuous. Let us introduce (closed) segments I¯i,r , where every I¯i,r lies strictly inside Ii,r , see formulas (17). Every function μi (x) is bounded on I¯i,r . mi  mi   I¯ I¯ . Thus functions ηi (x) and μi (x) are defined and bounded on the (closed) set xi . In this Let xi = r=0 i,r r=0 i,r case, by virtue of formulas (56), one obtains that

ηi (x; λ, α ) → μi (x; λi , αii )

(58)

uniformly with respect to x for x ∈ xi as soon as

λ ∈ ∗

and α → α .

By virtue of formulas (57) and (58) for α → α it is true that

ηi2 + ai − λi +

n 

αi j τ j → μ2i + ai − λi + αii θi ,

j=1

uniformly with respect to x for (x, λ ) ∈ xi ×  . ∗



μ2i

It follows from formulas (14), (16) that + εl,i − λi + αii θi  2αiiCi for x ∈ xi . In this case, by virtue of Theorem 4, one ∗ can choose positive constants αi0j for i = j such that inequalities (32) are satisfied for (x, λ, α ) ∈ xi ×  × A∗ . At the same

time α ij (i = j) are chosen such that if x ∈ xi \ xi , then inequalities (32) are satisfied. As above we assume that αi∗j  αi0j for i = j. It follows from the considerations made that δ i can be chosen such that |δ i | are less than any fixed δ > 0 as soon as |α − α | is sufficiently small. Uniting the found results, one arrives at the conclusion that inequalities (32) are satisfied for (x, λ, α) ∈ xi × ∗ × A∗ .  Proof of Theorem 5. It follows from Theorem 4 that the Cauchy problem (1), (3) is globally unique solvable for (x, λ, α ) ∈ ∗ I¯ ×  × A∗ . Let its solution be given by formulas (28). Since all Ai = 0, then all ui ≡ 0. It follows from statement 4 that function ηi (x ) =

ui (x ) ui ( x )

monotonically decreases on every interval Ii,r , r = 0, mi , see, for-

mula (29). If mi = 0, then Ii,0 = I. We remind that since ui ≡ 0, then ui (xi,r ) = 0 for all r = 0, mi . At the endpoints of every interval Ii,r , one has

lim

τi (x ) = A2i ,

ηi,0 = lim ηi (x ) = κs,i ,

lim

τi ( x ) = 0 ,

η = lim ηi (x ) = ±∞,

lim

τi ( x ) =

x→+0 x→xi,r

x→h−0

B2i ,

x→+0

± i,r

x→xi,r ±0

ηi,mi +1 = lim ηi (x ) = −κc,i , x→h−0

where r = 1, mi and Bi are unknown real constants. In Section 2.3 there were introduced (continuous) one-to-one mappings gi,r that allow one to define functions τ i (x) with x ≡ gi,r (ηi ) on every interval Ii,r . Then one can define functions

ςi,r (ηi ; λ ) = ηi2 + εl,i − λi +

n 

  αi j τ j gi,r (ηi ) ,

(59)

j=1

which is the right-hand side of the corresponding equation in (31). Index r corresponds to the interval Ii,r on which the expression ς i,r (ηi ; λ) is considered. Now one can integrate equations ηi = ςi,r (ηi ; λ ) on intervals Ii,r . The integration procedure is very similar to that, which is used in the proof of statement 2. Integrating equation ηi = ςi,r (ηi ; λ ) on every interval Ii,r , one gets

⎧ ηi,−1 ds ⎪ = x + ci,0 , ⎪ ⎪ ς ⎪ ηi (x ) i,0 (s; λ ) ⎪ ⎪ ⎪ ⎨ η (x ) ds

i

−

=x+c ,

i,r ςi,r (s; λ ) η ⎪ ⎪ ⎪ ⎪

ηi (x ) ⎪ ds ⎪ ⎪ = x + ci,mi , ⎩− + ςi,mi (s; λ ) ηi,m + i,r

x ∈ Ii,0 ; x ∈ Ii,r , r = 1, mi − 1;

(60)

x ∈ Ii,mi .

i

Substituting x = 0 + 0, x = xi,r+1 − 0, x = h − 0 into the first, second, and third lines, respectively, of (60), one derives ci,0 , ci,r , and ci,m , where r = 1, mi − 1. Substituting the found ci,r into (60), then substituting x = xi,1 − 0, x = xi,r + 0, x = xi,mi + 0 i

S.V. Tikhov and D.V. Valovik / Commun Nonlinear Sci Numer Simulat 75 (2019) 76–93

91

± into the first, second, and third lines, respectively, of the found formulas, and replacing ηi,0 , ηi,r , and ηi,mi +1 with κ s,i , ± ∞, and −κc,i , respectively, one arrives at the formulas

⎧

κs,i ds ⎪ ⎪ 0 < xi,1 = , ⎪ ⎪ ς −∞ i,0 (s; λ ) ⎪ ⎪

+∞ ⎨ ds 0 < xi,r+1 − xi,r = , r = 1, mi − 1; ⎪ −∞ ςi,r (s; λ ) ⎪ ⎪

+∞ ⎪ ds ⎪ ⎪ . ⎩ 0 < h − xi,mi = −κc,i ςi,mi (s; λ )

(61)

Formulas (61) express distances between adjacent zeros of ui . Since the left-hand sides in (61) are finite, then the righthand sides are also finite. Thus all the improper integrals in the right-hand sides converge. Summing up all the terms in (61), one finally arrives at

xi,1 + xi,2 − xi,1 + xi,3 − xi,2 + . . . + xi,mi − xi,mi −1 + h − xi,mi =

κs,i −∞

mi −1 +∞  ds + ςi,0 (s; λ ) −∞ r=1

ds + ςi,r (s; λ )

+∞

−κc,i

ds

ςi,mi (s; λ )

.

Taking into account notation (33)–(35), the found formula results in (36). It follows from the considerations made above that every solution (eigentuple) of the problem P(α) satisfies system of ¯ i. Eq. (36) with certain mi = m ¯ 1, . . . , λ ¯ n ) be a solution Now let us prove that every solution to (36) is an eigentuple of the problem P(α). Let λ¯ = (λ ¯ i . We consider the Cauchy problem (1), (3). Taking into account Theorem 4, one finds that solution to (36) with mi = m ∗ ui ≡ ui (x; λ¯ , α ) of the above mentioned Cauchy problem exists, is unique and defined globally for x ∈ I¯ and (λ¯ , α ) ∈  × A∗ . Using this solution, one can construct functions τi = u2i and ηi = ui /ui . Clearly, τi (0 ) = A2i and ηi (0 ) = κs,i . We do not assume that the condition ηi (h ) = −κc,i is fulfilled. Let, for certainty, ηi (h ) = ui (h )/ui (h ) = a > −κc,i . Using τ i and ηi , one can construct the expression

κs,i

ds

−∞

ςi,0 (s; λ¯ )

¯ i −1 m

+



T (λ¯ ; Ii,r ) +

r=1

+∞

ds

ςi,m¯ i (s; λ¯ )

a

= h,

(62)

which is similar to one from (36). Since the solution to the Cauchy problem (1), (3) is unique, then the integrands in (62) coincide with respective integrands in (36). At the same time λ = λ¯ satisfies (36). Subtracting the respective line of (36) from (62), one obtains

+∞

ds

ςi,m¯ (s; λ¯ )

a

−

+∞

−κc,i

i

ds

=−

ςi,m¯ (s; λ¯ ) i

a

ds

−κc,i

ςi,m¯ i (s; λ¯ )

= 0.

(63)

This implies that (63) is valid if and only if a = −κc,i . This means that λ¯ is an eigentuple of the problem P(α).



Proof of Corollary 2. Let us consider the difference

di =

mi  r=0

T (λ; Ii,r ) − i (λi ; mi ) =

 

mi −1

+

r=1

+∞ −∞

ds − ςi,r (s; λ )



κs,i

ds − ςi,0 (s; λ )

−∞ +∞

−∞

ds χi (s; λi )



+

κs,i

ds

−∞

χi (s; λi )

+∞

−κc,i

ds − ςi,mi (s; λ )

+∞ −κc,i

ds

χi (s; λi )

,

(64)

where the functions χ i (s; λi ) and ς i,r (s; λ) are defined in (16) and (59), respectively. Writing separately one of the terms in the right-hand side of (64), for example, for i = i and r = r  , where 1  r   mi − 1, one finds

+∞ −∞

ds − ςi ,r (s; λ )

+∞ −∞

ds = χi (s; λi )

+∞

−∞

   αi i θi (x ) − τi (x ) −  αi j τ j (x ) ds, ςi ,r (s; λ )χi (s; λi )

(65)

  where  = nj=1, j=i , x ≡ gi ,r (s ), x ≡ gi ,r (s ); here the functions gi,r (s ) are constructed with the help of the functions μi in the same manner as the functions gi,r (s). Since ui → vi and ui → vi uniformly w.r.t. x for x ∈ I¯ when λ ∈ ∗ and α → α , then gi,r (s ) → gi,r (s ) as α → α . Thus the term θi i (x ) − τi (x ) tends to 0 uniformly w.r.t. x for x ∈ I¯ when λ ∈ ∗ and α → α ; functions τ j (x), where j = i , are bounded  ∗ for all (x, λ, α ) ∈ I¯ ×  × A∗ . Since α → α , then α ij → 0 for i = j. This means that the quantities  αi j τ j (x ) can be made arbitrarily small as soon as αi j are sufficiently small, where j = i . Then it follows from (65) that

92

S.V. Tikhov and D.V. Valovik / Commun Nonlinear Sci Numer Simulat 75 (2019) 76–93

  

+∞ −∞

=

 

+∞ αi i θi (x ) − τi (x ) −  αi j τ j (x )  M ds d s  ςi ,r (s; λ )χi (s; λi ) 2αi i Ci − δ −∞ χi (s; λi )





M 2αi i Ci − δ

+∞ −∞



ds

(s2 + εl,i − λi )2 + 2αi i Ci

,

(66)

where quantities M and δ depend on λ, α. It follows from the considerations made that M → 0 as soon as α → α ; it follows from statement 4 that δ can be made arbitrarily small as soon as α → α . The integral in the right-hand side of (66) converges, see proof of Theorem 2. Thus, the difference (65) tends to 0 as α → α . Since the computations made are true for every summand in (64) and for all i, then one arrives at the conclusion that di → 0 as α → α . Similarly one can prove that differences

κs,i −∞

ds − ςi,0 (s; λ )

κs,i −∞

ds χi (s; λi )

also tend to 0 as α → α for all i.

and

+∞

−κc,i

ds − ςi,mi (s; λ )

+∞ −κc,i

ds

χi (s; λi )



Proof of Theorem 6. By virtue of Corollary 2, for any  > 0 there exists   > 0 such that absolute values of the left-hand sides in formulas (38) are less than  as soon as |α − α | <   . Zeros of the right-hand sides in (38) are eigenvalues of the problems Pi . It follows from Theorem 2 that every problems Pi has infinitely many (positive) eigenvalues

λi,m with an accumulation point at infinity. Among eigenvalues

λi,m there are infinitely many eigenvalues

λi,mr satisfying Corollary 1. Let us fix the intervals ∗i = (εs,i , λ∗i ), where λ∗i are chosen sufficiently big in order that the intervals (εl,i , λ∗i ) certainly contain eigenvalues

λ defined in Corollary 1. Due to Corollary 1 for any eigenvalue

λ ∈ ∗ there exists an interval i,mr

i,mr

i

Ui,m  containing

λi,mr such that the right-hand side corresponding to

λi,mr in (38) takes values of different signs on oppor site endpoints of this interval. Since the left-hand sides can be made arbitrarily small, the right-hand sides do not depend on α ij for i = j, are continuous with respect to λi and change signs crossing

λi,mr , then in the above introduced intervals for ¯ i,m satisfying system (36).  sufficiently small α ij with i = j there exist values λ r

Acknowledgements This work was supported by the Russian Science Foundation [grant number 18-71-10015]. References [1] [2] [3] [4] [5] [6] [7] [8] [9] [10] [11] [12] [13] [14] [15] [16]

[17] [18] [19] [20] [21] [22] [23] [24] [25] [26] [27] [28]

Valovik DV. Nonlinear multi-frequency electromagnetic wave propagation phenomena. J Opt 2017;19(11). article no. 115502 (16 pages). Shibata T. Variational methods for nonlinear multiparameter elliptic eigenvalue problems. Nonlinearity 1997;10(5):1319–29. Vainberg MM. Variational methods for the study of nonlinear operators, Holden-Day series in mathematical physics. 1st. ed. Holden-Day; 1964. Ambrosetti A, Rabinowitz PH. Dual variational methods in critical point theory and applications. J Funct Anal 1973;14(4):349–81. Rabinowitz PH. Some global results for nonlinear eigenvalue problems. J Funct Anal 1971;7(3):487–513. Browne PJ. A completeness theorem for a nonlinear multiparameter eigenvalue problem. J Differ Equ 1977;23(2):285–92. Krasnosel’skii MA. Topological methods in the theory of nonlinear integral equations. Pergamon Press, Oxford – London – New York – Paris; 1964. Vainberg M, Trenogin V. Theory of branching of solutions of nonlinear equations. Monographs and textbooks on pure and applied mathematics; Wolters-Noordhoff B.V.; 1974. Turyn L. Perturbation of linked eigenvalue problems. Nonlinear Anal 1983;7(1):35–40. Valovik DV. On the problem of nonlinear coupled electromagnetic te-tm wave propagation. J Math Phys 2013;54(4):042902(14). Smirnov YG, Valovik DV. Problem of nonlinear coupled electromagnetic te-te wave propagation. J Math Phys 2013;54(8):083502(13). Valovik DV. Integral dispersion equation method to solve a nonlinear boundary eigenvalue problem. Nonlinear Anal Real World Appl 2014;20(12):52–8. Malkin IG. Some problems in the theory of nonlinear oscillations. Cleveland, Ohio: U.S. atomic energy commission, technical information service; 1956. Gokhberg IT, Krein MG. Introduction in the theory of linear nonselfadjoint operators in hilbert space. Amer. Math. Soc.; 1969. Kato T. Perturbation theory for linear operators. Heidelberg: Springer; 1966. Boardman AD, Egan P, Lederer F, Langbein U, Mihalache D. Third-Order nonlinear electromagnetic TE and TM guided waves. Elsevier sci. Publ. North-Holland, Amsterdam London New York Tokyo; 1991. Reprinted from Nonlinear Surface Electromagnetic Phenomena, Eds. H.-E. Ponath and G. I. Stegeman. Valovik DV. On the existence of infinitely many nonperturbative solutions in a transmission eigenvalue problem for nonlinear helmholtz equation with polynomial nonlinearity. Appl Math Model 2018;53:296–309. Schürmann HW, Smirnov YG, Shestopalov YV. Propagation of te-waves in cylindrical nonlinear dielectric waveguides. Phys Rev E 2005;71(1):016614(10). Akhmediev NN, Ankevich A. Solitons, nonlinear pulses and beams. London: Chapman and Hall; 1997. Shen YR. The principles of nonlinear optics. New York–Chicester–Brisbane–Toronto–Singapore: John Wiley and Sons; 1984. Smirnov YG, Valovik DV. Guided electromagnetic waves propagating in a plane dielectric waveguide with nonlinear permittivity. Phys. Rev. A 2015;91(1). 013840 (6 pages). Binding PA, Browne PJ. Asymptotics of eigencurves for second order ordinary differential equations, i. J Differ Equ 1990;88(1):30–45. Binding PA, Browne PJ. Asymptotics of eigencurves for second order ordinary differential equations, ii. J Differ Equ 1991;89(2):224–43. Turyn L. Sturm–liouville problems with several parameters. J Differ Equ 1980;38(2):239–59. Shibata T. New asymptotic formula for eigenvalues of nonlinear sturm-liouville problems. J Anal Math 2007;102(1):347–58. Binding PA. On the use of degree theory for nonlinear multiparameter eigenvalue problems. J Math Anal Appl 1980;73(2):381–91. Atkinson EV, Mingarelli AB. Multiparameter eigenvalue problems. sturm-Liouville theory. USA: CRC Press; 2011. Podlevskii BM. Newton’S method as a tool for finding the eigenvalues of certain two-parameter (multiparameter) spectral problems. Comput Math Math Phys 2008;48(12):2140–5.

S.V. Tikhov and D.V. Valovik / Commun Nonlinear Sci Numer Simulat 75 (2019) 76–93 [29] [30] [31] [32] [33] [34] [35] [36] [37] [38] [39] [40] [41] [42] [43] [44] [45] [46] [47] [48]

93

Faierman M. Two-Parameter eigenvalue problems in ordinary differential equations. Longman Scientific and Technical, UK; 1991. Sleeman BD. Multiparameter spectral theory and separation of variables. J Phys A 2008;41(1):015209(20). Amer MA. Constructive solutions for nonlinear multiparameter eigenvalue problems. Comput Math Appl 1998;35(11):83–90. Browne PJ, Sleeman BD. Nonlinear multiparameter sturm-liouville problems. J Differ Equ 1979;34(1):139–46. Adams MJ. An introduction to optical waveguides. Chichester – New York – Brisbane – Toronto: John Wiley & Sons; 1981. Pontryagin LS. Ordinary differential equations. London: Pergamon Press; 1962. Eleonskii PN, Oganes’yants LG, Silin VP. Cylindrical nonlinear waveguides. Sov Phys JETP 1972;35(1):44–7. Landau LD, Lifshitz EM, Pitaevskii LP. Course of theoretical physics (vol. 8). electrodynamics of continuous media. Oxford: Butterworth-Heinemann; 1993. Eleonskii PN, Oganes’yants LG, Silin VP. Structure of three-component vector fields in self-focusing waveguides. Sov Phys JETP 1973;36(2):282–5. Boardman AD, Twardowski T. Theory of nonlinear interaction between te and tm waves. J Opt Soc Am B 1988;5(2):523–8. Boardman AD, Twardowski T. Transverse-electric and transverse-magnetic waves in nonlinear isotropic waveguides. Phys Rev A 1989;39(5):2481–92. Valovik DV. Novel propagation regimes for te waves guided by a waveguide filled with kerr medium. J Nonlinear Opt Phys Mater 2016;25(4). 1650051 (17 pages). Xie P, Zhang Z-Q. Multifrequency gap solitons in nonlinear photonic crystals. Phys Rev Lett 2003;91(21). 213904 (4 pages). Skryabin DV, Biancalana F, Bird DM, Benabid F. Effective kerr nonlinearity and two-color solitons in photonic band-gap fibers filled with a raman active gas. Phys Rev Lett 2004;93(14). 143907 (4 pages). Tricomi FG. Differential equations. New York: Blackie & Son Limited; 1961. Theotokoglou EE, Panayotounakos DE. Nonlinear asymptotic analysis of a system of two free coupled oscillators with cubic nonlinearities. Appl Math Model 2017;43:509–20. Cuevas J, Kevrekidis PG, Saxena A, Khare A. Pt-symmetric dimer of coupled nonlinear oscillators. Phys Rev A 2013;88(3). 032108 (11 pages). Borkowski L, Perlikowski P, Kapitaniak T, Stefanski A. Experimental observation of three-frequency quasiperiodic solution in a ring of unidirectionally coupled oscillators. Phys Rev E 2015;91(6). 062906 (7 pages). Tao M. Simply improved averaging for coupled oscillators and weakly nonlinear waves. Commun Nonlinear Sci Numer Simul 2019;71:1–21. In V, Kho A, Neff JD, Palacios A, Longhini P, Meadows BK. Experimental observation of multifrequency patterns in arrays of coupled nonlinear oscillators. Phys Rev Lett 2003;91(24). 244101 (4 pages).