Spectra of Elliptic Operators on Quantum Graphs with Small Edges

: We consider a general second order self-adjoint elliptic operator on an arbitrary metric graph, to which a small graph is glued. This small graph is obtained via rescaling a given ﬁxed graph γ by a small positive parameter ε . The coefﬁcients in the differential expression are varying, and they, as well as the matrices in the boundary conditions, can also depend on ε and we assume that this dependence is analytic. We introduce a special operator on a certain extension of the graph γ and assume that this operator has no embedded eigenvalues at the threshold of its essential spectrum. It is known that under such assumption the perturbed operator converges to a certain limiting operator. Our main results establish the convergence of the spectrum of the perturbed operator to that of the limiting operator. The convergence of the spectral projectors is proved as well. We show that the eigenvalues of the perturbed operator converging to limiting discrete eigenvalues are analytic in ε and the same is true for the associated perturbed eigenfunctions. We provide an effective recurrent algorithm for determining all coefﬁcients in the Taylor series for the perturbed eigenvalues and eigenfunctions.


Introduction
Spectral properties of elliptic operators on graphs with small edges are of a special interest and attract quite a lot of attention. Small edges is a specific singular geometric perturbation, which can be introduced due to the nature of graphs. One of early results on graphs with small edges was about a model of a woven membrane, see [1,2], which was shown to be approximated by a two-dimensional operator on an Euclidean domain covered by such membrane.
The graphs with finitely many small edges attracted more attention and one of early results stated that that a general vertex condition in a graph can be approximated in the norm resolvent sense via ornamenting by small edges supporting a magnetic field and with delta-coupling at the vertices [3]. An essential progress was made in very recent works [4,5] just few years ago. Here, Schrödinger operators on general graphs with arbitrarily placed small edges were considered. The main obtained result stated that under a certain nonresonance condition, see Condition 3.2 in [4], the perturbed operator converged to a certain limiting operator in the norm resolvent sense and the estimate for the convergence rate was established. The limiting operators was defined on a graph, in which the small edges were replaced by the vertices to which they shrink and at such vertices certain limiting boundary conditions were determined. These results were further developed in [6,7], where a general elliptic operator with varying coefficients was considered on a graph with small edges rescaled by means of a single small parameter. The coefficients in the differential expression and in the boundary conditions were allowed to depend analytically on the same small parameter. It was shown that certain parts of the resolvent of such operator depended analytically in ε and were represented by converging Taylor series. This allowed to represent the perturbed resolvent by a converging Taylor-like series with effective estimates for the remainders.
A next natural question is the behavior of the spectrum. This question was addressed in [4] for the aforementioned Schrödinger operator on a general graph with small edges. The convergence of the spectrum and corresponding spectral projectors was established. In particular, it was shown that the non-resonance condition was important and without it, the convergence of the spectrum could fail. We also mention few recent papers [8][9][10][11], where the resolvents and spectra were studied for some toy models represented by very simple graphs with small edges.
In the present work, we continue studying the model proposed in [6,7]. Namely, we consider a general self-adjoint second order elliptic operator on an arbitrary graph, to which a small graph is glued. The small graph is obtained via rescaling a given fixed graph by a small positive parameter ε, see Figures 1 and 2. The coefficients of the differential expression are varying and can additionally depend on ε. For the coefficients on fixed edges, this dependence on ε is analytic, while the coefficients on small edges can be even meromorphic. The boundary conditions are of a general matrix form with the matrices analytically depending on ε. The limiting operator is known thanks to the results in [6,7]. We show that the spectrum of the perturbed operator converges to that of the limiting operator and we provide an estimate for the distance between these two spectra. We also establish the convergence of the corresponding spectral projectors. Our second main result states that the eigenvalues of the perturbed operator converging to limiting discrete eigenvalues are analytic in ε. A similar analyticity property is established for the associated eigenfunctions. We provide an effective recurrent algorithm for determining the coefficients in the Taylor series for both eigenvalues and eigenfunctions. Once the coefficients are found, the sums of the Taylor series are exactly the eigenvalues and the eigenfunctions of the perturbed operator and in this sense, we can say that these eigenvalues and eigenfunctions are found explicitly. The paper is organized as follows. In the next section, we describe the problem, introduce auxiliary notations, formulate our main results and discuss their main features. The third section is devoted to proving the convergence of the spectrum and of the associated spectral projectors. In the fourth section, we prove the analyticity of the perturbed eigenvalues and the associated eigenfunctions, while in the fifth section, we describe an algorithm for determining the coefficients of their Taylor series.

Problem
Let Γ be a finite graph, that is, it contains finitely many vertices and edges. We assume that this graph contains no isolated vertices and edges and is metric, that is, each edge is equipped with a length, direction and a variable on it and on each edge, the usual Lebesgue measure is introduced. The graph Γ can contain edges of both finite and infinite lengths.
By M 0 we denote an arbitrary vertex in the graph Γ, while e i , i = 1, . . . , d 0 , are the edges incident to M 0 ; each possible loop in the family {e i } i=1,...,d 0 is counted twice. We introduce one more arbitrary finite metric graph γ containing only edges of finite lengths and no isolated vertices or edges, see Figure 2.
We contract the graph γ in ε −1 times, where ε is a small positive parameter, that is, each edge e in the graph γ of a length |e| is replaced by an edge of the length ε|e|, while all vertices remain unchanged. The graph obtained under such contraction is denoted by γ ε .
The graph we study is obtained by gluing the graph γ ε to Γ at the vertex M 0 . This gluing is defined as follows. First we choose n arbitrary vertices M j , j = 1, . . . , n, in the graph γ ε . Using the same number n, we partition the edges e i , i = 1, . . . , d 0 , of the graph Γ into n disjoint non-empty groups {e i } i∈J j , j = 1, . . . , n, and replace the vertex M 0 in the graph Γ by its n copies, one copy for each group {e i } i∈J j . Here, J j are some disjoint non-empty sets of indices and n j=1 J j = {1, . . . , d 0 }, 1 n d 0 . The gluing then is made by considering the union of the graphs Γ and γ ε and identifying the vertices M j , j = 1, . . . , n, with the aforementioned copies of the vertex M 0 in the graph Γ, see Figure 1. This new graph is denoted by Γ ε . It consists of a fixed part Γ and of a small graph γ ε . In the graph Γ ε , the edges e i , i = 1, . . . , d 0 , are not incident to the same vertex M 0 . Instead of this, each group {e i } i∈J j is incident to the vertex M j in the graph γ ε .
Hereafter, we identify the graphs Γ and γ ε with the corresponding subgraphs in Γ ε . Each function defined on Γ ε is also supposed to be defined on Γ and γ ε and vice versa. On each edge in Γ and γ we introduce arbitrarily a direction and an associated variable. The chosen directions then is transferred to the graph Γ ε ; the same is done for the variables on the edges in the graph Γ.
The main object of our study is an operator H ε in L 2 (Γ ε ) with the differential expres-sionĤ This expression is defined on each edge of the graph Γ ε and the symbol x stands for the variable on a considered edge. The symbol i denotes the imaginary unit, V are some real bounded measurable functions defined on the graphs Γ and V γ . By S ε we denote a linear operator mapping L 2 (γ) onto L 2 (γ ε ) and acting as on each edge e ε in the graph γ ε .
We assume that V γ ∈ L 2 (γ), and that these functions are analytic in ε in the sense of the corresponding norms. The operator H ε is supposed to be uniformly elliptic, namely, we assume that there exists a fixed c H > 0 independent of x ∈ Γ, ξ ∈ γ and ε such that Thanks to the assumed analyticity in ε, it is sufficient to impose the above condition for ε = 0 and then it holds for sufficiently small ε with a possible less constant c.
where x i is the variable on the edge e i . For each vertex M ∈ Γ ε of a degree d(M) > 0, the following general vertex condition is imposed: where A M (ε) and B M (ε) are some matrices of size d(M) × d(M) analytic in ε. The matrix is supposed to be self-adjoint, where which is obviously equivalent to for all sufficiently small ε.
For an arbitrary graph we denoteẆ 2 2 ( · ) := e∈ · W 2 2 (e). The domain of the operator H ε is a subspace of a Sobolev spaceẆ 2 2 (γ ε ) and this subspace is formed by the functions obeying boundary conditions (4). The action of the operator H ε on such functions is defined by differential expression (1). According [6] (Theorem 2.1), the operator H ε is self-adjoint.

Assumptions and Notations
Here, we introduce additional notations and assumptions, which will be employed then to formulate our main results. We attach edges of infinite lengths (leads) e ∞ i , i ∈ J j , j = 1, . . . , n, to the vertices M j , j = 1, . . . , n, in the graph γ. The obtained graph is denoted by γ ∞ , see Figure 3. In this graph, the vertices M j serve as the origins for the attached edges e i , i ∈ J j . On the graph γ ∞ , namely, in the space L 2 (γ), we introduce an operator H γ with the differential expression where, we recall, e i are the edges in the graph Γ incident to the vertex M 0 , and ξ is a variable on the graph γ ∞ . The vertex conditions read as where the vectors U M (u) and U M (u) are introduced in the same way as in (2)  2 (γ) satisfying the imposed vertex conditions. It was shown in [6] (Section 4.2) that the operator H γ is self-adjoint and its essential spectrum coincides with [0, +∞).
The main condition we assume in the work is as follows.
(A) The operator H γ has no embedded eigenvalue at the bottom of its essential spectrum. The above assumption means that zero, which is the bottom of the essential spectrum, is not an eigenvalue of the operatorĤ γ . In other words, the boundary value problem has no non-trivial solutions square integrable on γ ∞ . In view of the definition of the differential expressionĤ γ in (6), each non-trivial solutions to problem (8) is to be a linear function on the leads e ∞ i . Therefore, Condition (A) excludes the existence of non-trivial solutions to problem (8) vanishing identically on e ∞ i . However, this problem can have non- . . , n, with some constants depending on the choice of the edge e ∞ i . In this case one usually says that the operator H γ has a virtual level at the bottom of its essential spectrum.
On each function u defined and continuous on the edges e ∞ i in the vicinity of the points M j we introduce the following vector where u e ∞ i stands for the restriction of u to the edge e ∞ i . In terms of the introduced notations, Condition (A) can be equivalently reformulated as follows: each non-trivial solution ϕ to problem (8) obeys the inequality We proceed to auxiliary notations. Let P Γ : L 2 (Γ ε ) → L 2 (Γ) and P γ ε : L 2 (Γ ε ) → L 2 (γ ε ) be the operators restricting a given function to the subgraphs Γ and γ ε of the graph Γ ε : In the sense of the decomposition the identity P Γ ⊕ P γ ε = I Γ ε holds, where I Γ ε stands for the identity mapping in L 2 (Γ ε ).
As we have said above, under Condition (A), the operator H γ can have a virtual level at the bottom of its essential spectrum. By ϕ (j) , j = 1, . . . , k, we denote linearly independent bounded non-trivial solutions to problem (8) obeying condition (9). We clearly have k d 0 .
If the operator H γ has no virtual level at the bottom of its essential spectrum, we let k := 0. We denote: Ψ (j) := U γ (ϕ (j) ), j = 1, . . . , k. Hereafter we suppose that the functions ϕ (j) are chosen so that the vectors Ψ (j) are orthonormalized in where M ∈ γ ∞ , e i (M) are the edges incident to the vertex M, the numbers ν i (M) are defined as in (5). Employing Formula (12) for M = M j , we need to continue the functions . . , n, and we do this as V . . . . . .
Then we let where the symbol 0 in the first row of the matrix A M 0 denotes the zero matrix of the size k × (d 0 − k), while in the second row the same symbol stands for the zero matrix of the size (d 0 − k) × k. The first matrix in the definition of the matrix B M 0 is of the size d 0 × d 0 and the symbols 0 stand for the zero matrices of respectively the sizes The numbers ν i (M 0 ) are defined as in (5) and, we recall, e i are the edges in the graph Γ incident to the vertex M 0 .
Let H 0 be an operator in L 2 (Γ) with the differential expression and with the vertex conditions where the matrices A  (14). According [6] (Theorem 2.1), the operator H 0 is self-adjoint. Then its spectrum consists of its essential part and a discrete spectrum.
We shall make use of the following spaces of continuous functions on graphs: In the case when the functions V γ ∈Ċ(γ) and are analytic in ε in the norms of these spaces, we define the spaceĊ 2 . The spectrum of an operator is denoted by σ ( · ). Given a self-adjoint operator A and two real numbers a, b such that a < b and a, b / ∈ σ (A), by L [a,b] (A) we denote a spectral projector of A associated with the segment [a, b].

Main Results
Now we are in position to formulate our main results. Our first result states the convergence of the spectrum. Theorem 1. Let Condition (A) be satisfied. As ε → +0, the spectrum of the operator H ε converges to the spectrum of the operator H 0 . Namely, for each bounded interval I on the real axis we have where C I is some constant independent of ε but depending on the choice of the segment I. For arbitrary real numbers a, b such that a < b, a, b / ∈ σ (H 0 ), the convergence of spectral projectors holds: Let λ 0 be a discrete eigenvalue of the operator H 0 , that is, an isolated eigenvalue of a finite multiplicity . By ψ (j) 0 , j = 1, . . . , , we denote the associated eigenfunctions orthonormalized in L 2 (Γ 0 ). Our second result states the analyticity of the perturbed eigenvalues converging to λ 0 and of the corresponding eigenfunctions. ε , j = 1, . . . , , can be chosen so that they are analytic in ε in the following sense: are analytic in the norm ofĊ 2 (γ).
In this work, we also describe a way of determining the coefficients of the Taylor series for the perturbed eigenvalues and eigenfunctions described in the previous theorem. This is our third main result. However, in order to formulate it, we need to introduce many additional notations and this is why it is more convenient to provide this theorem in Section 5, see Theorem 3.

Discussion and Possible Generalizations of Main Results
Here, we briefly discuss our main results and their possible generalizations. Our first result, Theorem 1, establishes the convergence of the spectrum of the operator H ε to that of the operator H 0 . Inequality (15) states the convergence of the spectrum in each compact subset in the complex plane and moreover, the distance between the perturbed and limiting spectra is estimated. Relation (16) means that we also have a convergence of the corresponding spectral projectors. The limiting spectral projector in fact corresponds to the operator H 0 only and there is no contribution from the subgraph γ ε . These convergence results are of same nature as in the case of the uniform resolvent convergence of the selfadjoint operators. However, the results in [6] on the approximation of the resolvent of H ε do not imply immediately the statement of Theorem 1 just by applying classical theorems on norm resolvent convergence from the spectral theory. The reason is that our operator H ε acts in the Hilbert space L 2 (Γ ε ) depending on ε and this dependence is singular in the sense that the subgraph γ ε shrinks to a single vertex. Moreover, it is very essential to have Condition (A) here since without it the convergence of the spectrum can fail, see examples in [4] and paper [12]. Since the graph Γ can contain also the edges of infinite lengths, the spectra of both operators H 0 and H ε can also involve a non-empty essential component. Theorem 1 then states also the convergence of such components as well as of the corresponding spectral projectors.
Our second theorem states the analytic dependence on ε of the perturbed eigenvalues converging to limiting discrete operators. The same result holds for the associated perturbed eigenfunctions. The restriction of these eigenfunctions to the subgraph Γ is analytic in the sense of the norms inẆ 2 2 (Γ) andĊ 2 (Γ). The restriction to γ ε is not analytic but it becomes analytic after applying the operator S −1 ε . This means that we consider the restriction of the eigenfunctions to γ ε and then we pass to the rescaled variable ξ = xε −1 , which ranges over γ. Then the considered restriction becomes analytic in the norm ofĊ 2 (γ).
The analytic dependence on ε of the perturbed eigenvalues and eigenfunctions imply that they can be represented by converging Taylor series (56)-(58). In Section 5 we describe an algorithm for determining all coefficients in these three series. Although in Section 5 we consider only the case of a simple limiting eigenvalue λ 0 , our technique also work perfectly for the case of a multiple eigenvalue λ 0 . However, in this case the formulae become too bulky and this is the reason why we restrict ourselves by considering only the case a simple limiting eigenvalue. We stress that since series (56)-(58) converge, we can regard the result of Section 5 as a way of finding explicitly the perturbed eigenvalues and the eigenfunctions. This is true in the sense that they are represented by converging series and there is a direct and effective way of determining their coefficients.
To demonstrate our results, we mention a simple example of a star graph with a small edge considered in [9,10], see Figure 4a. Here, the edges e ± are finite and of the lengths a ± , while the edge e ε is small and its length is ε. The differential expression for the operator H ε on the edges e ± reads aŝ where V 0 = V 0 (ξ) is some real valued continuous function. At the vertex M 0 we impose the Kirchhoff condition, that is, while the boundary vertices are subject to the Dirichlet condition. Then the limiting operator H 0 is defined one the graph formed by the vertices e − and e + , see Figure 4b, with the differential expressionĤ 0 := − d 2 dx 2 subject to the Dirichlet condition at all vertices. The eigenvalues of the latter can be found explicitly: Let λ 0 = π 2 m 2 a −2 − be a simple eigenvalue corresponding to some m ∈ N, then according to our result, there exists a unique eigenvalue λ ε of the operator H ε converging to λ 0 as ε → +0. The eigenvalue λ ε is simple and analytic in ε. Then our approach described in Section 5 gives the following formula for the leading terms in the Taylor series for λ ε : If all these small graphs are rescaled by the same small parameter ε, all our results remain true. In such situation, Condition (A) is to be assumed separately for each small graph. If each of the small graphs γ (j) ε is rescaled by means of a small parameter ε (j) and the parameters ε (j) are assumed to be independent, then under certain additional conditions for the coefficients in the differential expression and the boundary conditions it is still possible to obtain the results analyticity of the eigenvalues and the eigenfunctions with respect to all small parameters. These conditions are the same as those in [6], which ensured a similar analyticity of the resolvent.
One more way of generalization is to assume that the coefficients in the differential expression and the boundary conditions are infinitely differentiable in ε or just have power in ε asymptotic expansions instead of the analyticity property. In this case the convergence result in Theorem 1 remains true. The perturbed eigenvalues λ (j) ε in Theorem 2 now becomes infinite differentiable in ε or having power asymptotic expansions in ε, while a similar statement of the associated eigenfunctions can fail. In a particular case, when the limiting eigenvalue λ 0 is simple, the corresponding perturbed eigenfunction is nevertheless infinite differentiable in ε or has a power asymptotic expansion in ε. Such statements can be proved by a technique different from that used in the present work. The results of Section 5 also remains true with the only modification that now series for the eigenvalues and eigenfunctions are to be treated as asymptotic in ε. If λ 0 is a multiple eigenvalue, then the technique used results of Section 5 still allows one to construct power asymptotic series for all perturbed eigenvalues converging λ 0 . However, the corresponding eigenfunctions can have a more complicated dependence on ε and in general, power asymptotic series are no longer true even in the asymptotic sense.
In conclusion, let us also mention the following application of our main results. Quantum graphs often arise and are employed, for instance, as models of thin quantum waveguide and nanotubes and many other applications, see [13] (Chapter 7). Then graphs with small edges correspond to the case, when the lengths of some tubes or waveguides are much less than the lengths of the others. Usually, in such regime, it is extremely troublesome to employ numerical methods for solving boundary value problems on such graphs since one has to employ the numerical schemes with a very small step. At the same time, our results show how to find effectively the eigenvalues and eigenfunctions of the considered model up to an arbitrary small error. Namely, it is sufficient to employ the partial sums of series (56)-(58), the coefficients of which can be found by the algorithm presented in Theorem 3.
We fix an arbitrary bounded segment I ⊂ R and we define its subset I δ := {λ ∈ I : dist λ, σ (H 0 ) δ} for some fixed δ > 0, which will be chosen later. We are going to show that as ε is small enough, for all λ ∈ I δ with δ = Cε 1 2 the resolvent (H ε − λ) −1 is welldefined; here C is a constant independent of ε. This fact will imply that σ (H ε ) ∩ I ⊂ I \ I δ and this will yield estimate (15).
For λ ∈ I δ we consider the equation for the resolvent (H ε − λ)u = f with an arbitrary f ∈ L 2 (Γ ε ) and we rewrite it as According Theorem 2.1 in [6] and Theorem 2 in [7], the resolvent (H ε − i) −1 is welldefined and it can be represented as where the direct sum is understood in the sense of decomposition (10) and A 1 (ε) is a bounded operator in L 2 (Γ ε ) obeying the estimate where C is some constant independent of ε. We apply the resolvent (H ε − i) −1 to Equation (17) and then we substitute representation (18) into the obtained identity: We have: where I is the identity mapping in L 2 (Γ ε ) and I Γ is the identity mapping in L 2 (Γ) and for λ ∈ I δ the resolvent (H 0 − λ) −1 is well-defined. We substitute then (21) into (20) and invert the operator (H 0 − λ)(H 0 − i) −1 ⊕ I γ ε , where I γ ε is the identity mapping in L 2 (γ ε ). Then we get: As λ ∈ I δ , we have a standard estimate for the resolvent of the self-adjoint operator H 0 : where the resolvent is regarded as an operator in L 2 (Γ). Then, by (19) we immediately get where C is some fixed constant independent of ε, δ and λ. The right hand side in the above estimate is less than 1 2 provided δ = 2Cε 1 2 . Then estimate (23) yields the existence of the a bounded inverse operator I + A 3 (λ, ε) −1 in L 2 (Γ ε ). This existence allows us to solve Equation (22) as Hence, the resolvent (H ε − λ) −1 is well-defined and its action on arbitrary function f ∈ L 2 (Γ ε ) is given by the right hand side of the above identity.
We proceed to studying the convergence of the spectral projectors. The operator H 0 is obviously lower-semibounded and in view of the established convergence of the spectrum we can choose a fixed real negative constant µ such that σ (H 0 ) ⊂ [µ + 1, +∞) and σ (H ε ) ⊂ [µ + 1, +∞) for all sufficiently small ε. Then both the resolvents (H 0 − µ) −1 and (H ε − µ) −1 are well-defined. Moreover, it is possible to reproduce the proof of Theorem 2.1 in [6] for the operator (H ε − µ) −1 and to establish an identity similar to (10). Namely, the only point we need to confirm in the proof in [6] is the invertibility of the matrix X 0 − µG 0 used in [6] (Lemma 5.1), where X 0 , G 0 are some fixed self-adjoint matrices and G 0 is positive definite. This is obviously true provided µ is negative and its absolute value is large enough. Then the aforementioned identity similar to (10) reads as where are some bounded operator analytic in ε and where c i are some linear functionals on L 2 (Γ) and we recall that ϕ (i) are non-trivial solutions of problem (8). Having decomposition (10) in mind, we introduce an unitary operator U ε : Making the unitary transformation of the operator H ε by means of U ε , we rewrite (24) as The above definitions and identities (25) imply where A 8 : L 2 (Γ) → L 2 (Γ), A 9 : L 2 (γ) → L 2 (Γ), A 10 : L 2 (Γ) → L 2 (γ) and A 11 : L 2 (γ) → L 2 (γ) are bounded operators analytic in z 2 . Hence, the operators A 6 , A 7 are also analytic in z.
Since the operator H ε is self-adjoint, the unitary equivalent operator U ε H ε U −1 ε is self-adjoint in L 2 (Γ) ⊕ L 2 (γ). The established analyticity of the operators A 6 , A 7 and representation (27) mean that the operator where C is some constant independent of ε. We also observe that the σ (H ε ) = σ (U ε H ε U −1 ε ) and the spectra of the operator U ε H ε U −1 ε and of its resolvent A similar correspondence holds for the spectra of the operator H 0 and of its resolvent (H 0 − µ) −1 . In particular, this means that the spectrum of (U ε H ε U −1 ε − µ) −1 converges to that of (H 0 − µ) −1 .

Analyticity of Eigenvalues and Eigenfunctions
In this section, we prove Theorem 2. The first part of the theorem on existence and convergence of the perturbed eigenvalues λ (j) ε is implied immediately by Theorem 1. In the considered case the spectral projector L [λ 0 −δ,λ 0 +δ] (H ε ) is the projector on the space spanned over the eigenfunctions associated with the eigenvalues of H ε converging to λ 0 , while L [λ 0 −δ,λ 0 +δ] (H 0 ) is the projector on the eigenspace associated with the eigenvalue λ 0 . The convergence of these projectors implied by (16) means that the total multiplicity of the perturbed eigenvalues λ (j) ε converging to λ 0 coincides with the multiplicity of λ 0 . We proceed to proving the analyticity of the eigenvalues λ and it is analytic in ε 1 2 . We consider finitely many eigenvaluesλ (j) ε converging to a limiting eigenvalue (λ 0 − µ) −1 of a finite multiplicity. Then according to the results in [15] (Chapter VII, Section 3.1), we can apply the statements established in [15] (Chapter II, Section 6) for analytic self-adjoint operators in finite-dimensional spaces. This implies that the eigenvaluesλ ε , we then conclude that they are analytic in ε 1 2 . We also observe that the eigenfunctions ψ (j) ε are recovered from the eigenfunctionsψ ε,γ ∈ L 2 (Γ) ⊕ L 2 (γ) by the following formulae valid due to (31), (26): As a next step, we are going to prove that the eigenvalues λ (j) ε and the functions P Γ ψ ε,Γ are analytic in ε in the norms of the spacesẆ 2 2 (Γ) ∩Ċ 2 (Γ) andĊ 2 (γ). First of all we observe that owing to formulae (27), (29) and eigenvalue Equation (31) we have: ε,γ is analytic in ε 1 2 in the spaceĊ 2 (γ). The stated analyticity property means that where Λ γ,i are some analytic in ε functions. We consider the eigenvalue equation for λ (j) ε and ψ (j) ε as a boundary value problem for the differential equation subject to boundary conditions (3). The above equation is considered separately on subgraphs Γ and γ ε ; on the latter subgraph we also pass to the rescaled variable ξ = xε −1 . This leads us to the following equations: whereĤ γ (ε) is the following differential expression: The coefficients of both differential expressionsĤ(ε) andĤ γ (ε) are analytic in ε. Having this fact in mind, we substitute representations (32) into (33) and in view of the A similar procedure can be done with boundary conditions (3), in which we can also separate an analytic in ε part and a similar one multiplied by ε 1 2 . Considering then the obtained boundary value problems for Ψ (j) γ,i , we see immediately that they can be written as two equations where Ψ (j) 0 is an eigenfunction of the operator H 0 . We multiply the first equation in (35) by Ψ (j) 1 in L 2 (Γ ε ) and integrate by parts employing the second equation in (35): L 2 (Γ) = 0 as ε → +0, the obtained identity imply immediately that Λ

Taylor Series
In this section, we describe an algorithm for determining all coefficients in the Taylor series for the perturbed eigenvalues λ ε and the associated eigenfunctions ψ ε .
We substitute series (56), (57) into the eigenvalue equation H ε ψ ε = λ ε ψ ε , rewrite it as a boundary value problem on Γ ε and consider this problem separately on the subgraph Γ and γ ε . This leads us to equations similar to those on Γ in (34) with Λ (j) 1,Γ = 0. In these equations we expand the coefficients of the differential expressions into the power series in ε and equate the coefficients at the like powers of ε and we arrive at the following recurrent system of equations: where the differential expressions L Γ,p are defined as In the same way we substitute series (57) into boundary conditions (3), expand the matrices A M and B M into the power series in ε and equate the coefficients at the like powers of ε. This gives: where the matrices A The way of finding boundary value problems for the coefficients of series (58) follows the same lines but with an additional trick. Namely, first we introduce an auxiliary graph γ ex by attaching additional unit edges e ex i , i ∈ J j , j = 1, . . . , n, to the each vertex M j , j = 1, . . . , n, in the graph γ. The boundary vertices, being the end-points of the edges e ex i and not coinciding with M j , are denoted by M ex i , i ∈ J j , j = 1, . . . , n. We observe that the set of the vertices in the graph γ ex not coinciding with M ex i , i = 1, . . . , n, is in fact the set of the vertices of the graph γ ∞ . We introduce the notations We stress that since the eigenvalues λ ε and the eigenfunctions ψ ε are well-defined, the obtained boundary value problems for the functions ψ Γ,i and ψ γ,i are solvable. Using these boundary value problems, we are going to study the structure of the functions ψ Γ,i and ψ γ,i and to find the constants λ i . In order to do this, we first introduce auxiliary notations and mention some statements proved in [6].
First, we consider Γ as a separate graph and for each vertex M ∈ Γ we introduce matrices U M by formulae (11), wherẽ The matrices U M turn out to be unitary, see [6] (Section 4). Then, as above, for such vertices M we let P M to be the projector in C d(M) onto the eigenspace of the matrix U M associated with the eigenvalue −1, and P ⊥ M := E d(M) − P M . Finally, for M ∈ Γ, and for M ∈ γ ∞ we define: We also denote: The following lemma was proved in [6], see Lemma 5.2 in that work.
is solvable inẆ 2 2 (γ ex ) if and only if for each j = 1, . . . , k. Under these conditions, there exists the unique solution u * obeying the identities U γ (u * ), Ψ (j) The general solution of problem (42) reads as where c j are arbitrary constants.
A similar statement for the operator H 0 on the graph Γ can be proved exactly in the same way as the above lemma was proved in [6]. This statement is as follows.

Lemma 2.
Let λ 0 be a simple eigenvalue of the operator H 0 and ψ 0 be the associated eigenfunction normalized in L 2 (Γ). Given an arbitrary family of vectors g M ∈ P The general solution of problem (44) reads as where c is an arbitrary constant.
Equation (40) with i = 0 for ψ γ,0 and the first boundary condition in (39) are homogeneous. Hence, the function ψ γ,0 is a linear combination of the functions ϕ (i) ; we recall that the latter functions correspond to the virtual level at the bottom of the essential spectrum of the operator H γ . We have where a i,0 are some constants to be determined. The function ψ Γ,0 = ψ 0 is apriori known and in view of the boundary condition for this function at the vertex M 0 we see that U M 0 (ψ 0 ), Ψ (p) C d 0 = 0, p = k + 1, . . . , d 0 .
The functions ψ γ,i are solutions of boundary value problems (39)-(41) represented as where φ γ,i are solutions to the same problems as for ψ γ,i but satisfying the orthogonality condition The constants a i,p are determined by formula The constants λ i are given by Formula (55).

Conclusions
Here, we stress two important aspects of our results. The first of them concerns the analyticity of the eigenvalues. This is a rather surprising fact since the small edges are an example of a singular perturbation, which can not be treated by the methods of the regular perturbation theory presented, for instance, in [15]. Under singular perturbations, one can usually expect only the possibility to construct some asymptotic series but not to prove their convergence [16]. However, the specific feature of our problem is that in our case, series (56)-(58) do converge, and even uniformly in small ε. The second important point is that we can effectively find all coefficients in Taylor series (56)-(58) by the recurrent procedure described in Section 5. This algorithm is an adaption of the method of matching asymptotic expansions widely used for partial differential equations on Euclidean domains and manifolds [16]. The present work is the first where this method is employed and adapted for general graphs with small edges for finding eigenvalues and eigenfunctions.