Effective Behaviour of Critical-Contrast PDEs: Micro-resonances, Frequency Conversion, and Time Dispersive Properties. I

A novel approach to critical-contrast homogenisation for periodic PDEs is proposed, via an explicit asymptotic analysis of Dirichlet-to-Neumann operators. Norm-resolvent asymptotics for non-uniformly elliptic problems with highly oscillating coefficients are explicitly constructed. An essential feature of the new technique is that it relates homogenisation limits to a class of time-dispersive media.


Introduction
The research aimed at modelling and engineering metamaterials has been recently brought to the forefront of materials science (see, e.g., [19] and references therein).It is widely acknowledged that these novel materials acquire non-classical properties as a result of a careful design of the microstructure, which can be assumed periodic with a small enough periodicity cell.The mathematical machinery involved in their modelling must therefore include as its backbone the theory of homogenisation (see e.g.[7,4,43]), which aims at characterising limiting, or "effective", properties of small-period composites.A typical problem here is to study the asymptotic behaviour of solutions to equations of the type where for all ε > 0 the matrix A ε is Q-periodic, Q := [0, 1) d , non-negative symmetric, and may additionally be required to satisfy the condition of uniform ellipticity.
On the other hand, the result sought (i.e., the "metamaterial" behaviour in the limit of vanishing ε) belongs to the domain of the so-called time-dispersive media (see, e.g., [67,34,68,35]).For such media, in the frequency domain one faces a setup of the type where A is a constant matrix and B(ω 2 ) is a frequency-dependent operator in L 2 (R d ) taking the place of −ω 2 in (1), if, for the sake of argument, in the time domain we started with an equation of second order in time.If, in addition, the matrix function B is scalar, i.e., B(ω 2 ) = β(ω 2 )I with a scalar function β, the problem of the type appears in place of the spectral problem after a formal division by −β(ω 2 )/ω 2 , with frequency-dependent (but independent of the spatial variable) matrix A(ω 2 ).Thus, the matrix elements of A(ω 2 ), interpreted as material parameters of the medium, acquire a nontrivial dependence on the frequency, which may lead to their taking negative values in certain frequency intervals.The latter property is, in turn, characteristic of metamaterials [71].It is therefore of paramount interest to understand how inhomogeneity in the spatial variable (see (1)) can lead, in the limit ε → 0, to frequency dispersion, and, in particular, to uncover the conditions on A ε sufficient for this.A result, which from the above perspective can be seen as negative, is provided by the homogenisation theory in the uniformly strongly elliptic setting (i.e., both A ε and (A ε ) −1 are bounded uniformly in ε).Here one proves (see [73,10] and references therein) the existence of a constant matrix A hom such that solutions u ε to (1) converge to u hom satisfying −div A hom ∇u hom − ω 2 u hom = f, which leaves no room for time dispersion.This negative result also carries over to vector models, including the Maxwell system.If the uniform ellipticity assumption is dropped, the asymptotic analysis of the problems (1) becomes more complicated.By employing the technique of two-scale convergence, first Zhikov [74,76] then Bouchitté,Bourel,and Felbacq [12,13,14] obtained a related effective problem of the form (2). The former works treat the critical-contrast model of the type (1), while the latter are devoted to an associated scattering problem.Here, by "critical contrast" one means that the components of a composite dielectric medium have material properties in a proper contrast to each other, governed by the size of periodicity cell (see Section 2 for further details).More recently, Kamotski and Smyshlyaev [44] developed a general two-scale compactness argument for the analysis of "degenerate" homogenisation problems, of which critical-contrast problem (1) is a particular case.In a parallel development, prompted by [57], the two-scale convergence approach has been used [46,15,49,50] to mathematically justify the emergence of artificial magnetism in the case of "resonant" metallic inclusions, i.e. in the case when the conductivity is high and scaled appropriately with the microstructure size, with the "split-ring" geometry of the inclusions.In yet another development, operator-theoretic techniques were used by Hempel and Lienau [41] to prove that the limit spectrum in the case of critical-contrast periodic media has a band-gap structure, see also [37] for a result concerning the asymptotics of the integrated density of states as well as [30] for the asymptotic spectral analysis of periodic Maxwell problems with the wavenumber-frequency pair situated in the vicinity of the lowest of the light-lines of the material components.
Although well received, the above results fall short of establishing a rigorous one-to-one correspondence between homogenisation limits for critical-contrast media and time dispersion in the effective medium.This is due to the following: (i) the additional assumptions imposed only permit to treat a limited set of models (curiously, excluding even the one-dimensional version of the problem, cf.[21,20]); (ii) the control of the convergence error, for the solution sequences as well as the spectra of the problems, is lacking, due to the rather weak convergence of solutions claimed.A more general theory, akin to that of Birman and Suslina [10,11] in the moderate-contrast case, is therefore required.The present paper attempts to suggest such a theory.
The benefits of the novel unified approach as developed henceforth are these, in a nutshell: (1) Being free from additional assumptions on the geometry and PDE type, it can be successfully applied in a consistent way to diverse problems motivated by applications; (2) It can be viewed as a natural (albeit non-trivial) generalisation of the approach of Birman and Suslina in the uniformly elliptic case; (3) The analysis is shown to be reducible by purely analytical means to an auxiliary uniformly elliptic problem; the latter, unlike the original problem, is within the reach of robust numerical techniques; (4) The error bounds are controlled uniformly via norm-resolvent estimates (yielding the spectral convergence as a by-product); (5) Not only is the relation of the composite to the corresponding effective time-dispersive medium made transparent (showing the artificial introduction of second ("fast") variable via the two-scale asymptotics to be unnecessary from the technical point of view), but the approach can be also seen to offer a recipe for the construction of such media with prescribed dispersive properties from periodic composites whose individual components are non-dispersive.The analytical toolbox we propose also allows us1 : (i) to explicitly construct spectral representations and functional models for both homogenisation limits of critical-contrast composites and the related timedispersive models; (ii) on this basis, to solve direct and inverse scattering problems in both setups.It thus paves the way to treating the inverse problem of constructing a metamaterial "on demand", based on its desired properties.We shall therefore reiterate that the present paper can be seen as an example of how surprisingly far one can reach by a consistent application of the existing vast toolbox of abstract spectral theory.
In [24,26] (see also an earlier paper [27] dealing with critical-contrast homogenisation on R) we considered a rather simple model of a high-contrast graph periodic along one axis.A unified treatment of critical-contrast homogenisation was proposed and carried out in three distinct cases: (i) where neither the soft nor the stiff component of the medium is connected; (ii) where the stiff component of the medium is connected; (iii) where the soft component of the medium is connected.
In the present paper we turn our attention to the PDE setup, and we focus on the scalar case, leaving the treatment of vector problems, in particular the full Maxwell system of electromagnetism as well as the systems of two-and three-dimensional elasticity, to a future publication.In view of keeping technicalities to a bare minimum, and at the same time making the substance of the argument as transparent as possible, we consider two classical models, see Section 2 for details.The main ingredients of the theory, which is formulated in abstract, yet easily applicable, terms, remain virtually unchanged in more general models.One such generalisation is briefly discussed in Section 6.
The analytical backbone of our approach is the so-called generalised resolvent, or in other words the resolvent of the original operator family sandwiched by orthogonal projections to one of components of the medium ("soft" one, see Section 2 for details).In its analysis, we draw our motivation from the celebrated general theory due to Neumark [54,55] and the follow-up work by Strauss [62].An explicit analysis of Dirichlet-to-Neumann (DN) maps (separately for the components comprising the medium) becomes necessary to facilitate the use of the well-known Kreȋn resolvent formula.The corresponding analysis is based on a version of Birman-Kreȋn-Višik theory [8,47,48,72] as proposed by Ryzhov [59].In adopting the mentioned approach, we deviate from the boundary triples theory utilised previously in our ODE analysis [24,26].Still, the theory developed in [59] can be seen as a "version" of the latter, in that it attempts to write the second Green's identity in operator-theoretic terms as where A is an operator in a Hilbert space H (e.g., a Laplacian would lead to the classical Green's identity with Γ 0 and Γ 1 defined as traces of the function and of its normal derivative on the boundary, respectively) and (H, Γ 0 , Γ 1 ) is the boundary triple, consisting of an auxiliary Hilbert space (e.g., the L 2 space over the boundary) and a pair of operators onto H. Unfortunately, it is well known that the direct approach via (3), although possible, encounters problems in the PDE setting (we refer the reader to a review [33] of the state of the art in the theory, see also Section 2 below), thus necessitating an alternative.For this we have selected the boundary triple theory of [59], as a natural fit for the analysis we have in mind.
In order to put the results of our work into the correct context, we mention that already in the extensively studied setup of double-porosity models (see Model I, Section 2) we are not only able to develop the operatortheoretical approach, but also to extend the existing results on spectral convergence in at least two ways: firstly, by ascertaining the rate convergence as O(ε 2/3 ), and secondly, by disposing of the assumption on the eigenvectors of the Dirichlet Laplacian on the soft component, which was previously considered essential.Moreover, the end product of our analysis in this case, Theorem 5.9, puts the frequency-dispersive (and thus time-dispersive) properties of the homogenised medium under the spotlight, allowing to write the corresponding "Zhikov function" B(z) explicitly and thus putting on a firm ground the possibility to view media with double porosity as frequency-converting devices.Furthermore, an asymptotic development of our approach allows us to upgrade the order of convergence estimates from O(ε 2/3 ) to O(ε α ), α > 1, which we shall carry out in the forthcoming second part [25] of our work.
Finally, we compare our results to earlier approaches to the spectral analysis of scalar PDE problems with critical contrast, which thus far have been applied in the particular geometric setup of "disperse" soft inclusions in a connected stiff medium, see [74], [22].Similarly to [22], our machinery goes well beyond the two-scale convergence techniques of [74], in that the convergence of spectra is now a simple consequence of the convergence of solutions to the original PDE in the operator-norm topology.In contrast to the norm-resolvent asymptotics of [22], where order O(ε) estimates were obtained for the difference between the resolvents of the original problems and the approximating operators, the estimate of the present paper, see Theorem 4.6, provides an explicit expression for a modified approximating homogenised operator, so that the new convergence estimates are of the higher order O(ε 2 ).In addition, in contrast to [74] and [22], we consider arbitrary geometries, where the soft inclusions need not be disconnected and the stiff component need not be connected.Yet another novel feature of the present work is an explicit description of the related effective macroscopic time-dispersive media, see the discussion above Section 5.This establishes an explicit link between critical-contrast composites and metamaterials, which the earlier works were somewhat lacking.
In conclusion, we would like to mention that although our techniques necessarily rely upon what might be seen as a purely abstract framework, the exposition of the paper is effectively complemented for the benefit of more applied-minded audience by that of our earlier paper [24], which attempts to make our ideas more transparent, by resorting to the asymptotic analysis of the eigenvectors, as opposed to the resolvents, of the corresponding self-adjoint operators.

Problem setup and some preliminaries
In the present work, we consider the problem (1) under the following assumptions: where Two distinct models will be of a particular interest to us, see Fig. 1.For simplicity but without loss of generality, in both of them the coefficient a 2 will be assumed identically equal to unity.In Model I, which is unitary equivalent to the model of [41], [37], we will assume further that Q soft is a simply connected domain with a smooth boundary Γ, such that the distance from its boundary to the boundary of the unit cube Q is positive, cf.[74], [22].In Model II, in contrast, we will assume that Q stiff is a simply connected domain with a smooth boundary Γ, separated from the boundary of the unit cube Q.The assumption that Γ is smooth in both cases cannot in general be dropped, as we will require DN maps pertaining to this boundary to be well-defined as pseudo-differential operators of order one [42], [36]; also see e.g.[3] and references therein for the treatment of DN maps in the case of Lipschitz boundary, which can be used to generalise our approach to this case.It has to be noted, however, that recent developments in the application of pseudo-differential calculus to the extension theory on non-smooth domains, see [1] and references therein, indicate that nonsmooth boundaries may also be treated using our strategy.To summarise, the setup of Model I is therefore that of soft inclusions in a stiff connected "matrix", Model II represents the "dual" case of stiff inclusions in a soft connected matrix.
In what follows, we present a systematic treatment of Models I and II, followed by the list of amendments (see Section 6) allowing us to obtain the same results in the general case (4), with non-piecewise constant coefficient a 2 in (4).This choice of exposition is justified, since the treatment of the general case (which does not seem to be well-motivated from the point of view of real-world applications) is in no way different to the treatment of Models I and II.
In both cases we shall deal with the resolvent (A ε −z) −1 of a self-adjoint operator in L 2 (R d ) corresponding to the problem (1), so that solutions of the latter are expressed as u ε = (A ε − z) −1 f with z = ω 2 .Here the operator family A ε can be thought of as defined by the coercive, but not uniformly coercive, forms In what follows, we will always assume that z ∈ C is separated from the spectrum of the original operator family, more precisely, we assume that z ∈ K σ , where After we have established the operator-norm asymptotics of (A ε − z) −1 for z ∈ K σ , the result is extended by analyticity to a compact set the distance of which to the spectrum of the leading order of the asymptotics is bounded below by σ.
In view of dealing with operators having compact resolvents, it is customary to apply either a Floquet (see [37]) or Gelfand (see [38]) transform, of which we use the latter, in line with [10].
2.1.Gelfand transform and direct integral.Consider the Gelfand transform of a function The Gelfand transform thus defined is a unitary operator between L 2 (R d ) and L 2 (Q × Q ), and the inverse transform is expressed by the formula where U is extended to R d × Q by Q-periodicity in the spatial variable.For the scaled version of the above transform, for u ∈ L 2 (R d ) we set and the inverse of the transform ( 7) is given by As in [10], an application of the Gelfand transform G to the operator family A ε corresponding to the problem (1) yields the two-parametric operator family subject to periodic boundary conditions on the boundary of εQ and defined by the corresponding closed coercive sesquilinear form.For each ε > 0, the operator A ε is then unitary equivalent to the von Neumann integral (see e.g.[9, Chapter 7]) Next, we introduce a unitary rescaling for the operator family where the unitary Φ ε acts on u ∈ L 2 (εQ) by the following rule: 2 The formula ( 5) is first applied to continuous functions U with compact support, and then extended to the whole of L 2 (R d ) by continuity.
To simplify notation, we choose to keep the same symbol ε , where t = τ /ε ∈ ε −1 Q , hoping that this does not lead to a confusion.The operator family subject to periodic boundary conditions on ∂Q, and it is uniquely defined by its form.
Clearly, one can first apply a unitary rescaling to the operator family A ε , followed by an application of the Gelfand transform (5), leading to the same operator family A (τ ) ε , so that the fibre representation holds, where ∼ = stands for unitary equivalence.Although the main body of the paper makes no explicit reference to the operators ε , the order in which we apply the (scaled) Gelfand transform and the unitary rescaling, and thus the very appearance of the operator A (t) ε , is required for the transparency of the analysis of Section 5.4.1.Therein, we find it instructive to utilise the scaled inverse Gelfand transform, in view of the natural appearance of an integral over an expanding dual cell ε −1 Q , which is, in its turn, required for the main result of the section, Theorem 5.5.
Similar to [37] and facilitated by the abstract framework of [59], instead of the form-based definition of the operators A (τ ) ε we will consider them as operators of transmission problems, see [60] and references therein.Henceforth we treat the cube Q as a torus by identifying its opposite sides, and Q soft and Q stiff are viewed as subsets of this torus so that, in particular, ∂Q soft = ∂Q stiff = Γ.For each ε, τ, the mentioned transmission problem is formulated as finding a function that solves (in the variational, or weak, sense) the boundary-value problem (BVP) with f ∈ H : where u + , ∂u/∂n + , u − , ∂u/∂n − are the limiting values of the function u and its normal derivatives on Γ from inside Q soft and Q stiff , respectively; the vectors n + and n − = −n + are outward normals to Γ with respect to Q soft and Q stiff .By a classical argument the weak solution of the above problem is shown to be equal to (A 2.2.Boundary operator framework.Following [59] (cf.[6], [16] and references therein for alternative approaches), which is based on the ideas of the classical Birman-Kreȋn-Višik theory (see [8,47,48,72]), the linear operator of the transmission boundary-value problem is introduced as follows.Let H := L 2 (Γ), and consider the "τ -harmonic" lift operators Π soft and Π stiff defined on φ ∈ H via3 with periodic conditions on ∂Q added for each component as required.These operators are first defined on φ ∈ C 2 (Γ), in which case the corresponding solutions u φ can be seen as classical.The standard ellipticity estimate allows one to extend both τ -harmonic lifts to bounded operators from H to L 2 (Q stiff(soft) ), so that the functions u φ are treated as distributional solutions4 of the respective BVPs, see e.g., [53, Theorem 3.2, Chapter 5], [52,Theorem 4.25].The solution operator Π : is defined as follows: Πφ := Π soft φ ⊕ Π stiff φ, φ ∈ H. Consider the self-adjoint operator family A 0 (where we drop the superscript (τ ) and the subscript ε for brevity) to be the Dirichlet decoupling of the operator family A (τ ) ε , i.e. the operator of the Dirichlet BVP on both Q stiff and Q soft , with periodic boundary conditions on ∂Q where appropriate, defined by the same differential expression as A ran Π stiff(soft) .In the same way we introduce the operator Γ 0 and its null extension Γ 0 to the domain dom A 0 ran Π.
The solution operators S stiff z , S soft z of the BVPs Γ stiff 0 u φ = φ, and Γ soft 0 u φ = φ, are defined as linear mappings from φ to u φ , respectively.These operators are bounded from H to L 2 (Q stiff ) and L 2 (Q soft ), respectively, and admit the representations Clearly, the operator S soft z is ε-independent (which is a result of applying the rescaling Φ ε ), and the operator S stiff z admits the obvious estimate uniformly in τ ∈ Q , z ∈ K σ , in the uniform operator-norm topology, owing to the fact that Finally, the solution operator and is bounded.Introducing the orthogonal projections P soft and P stiff from H onto L 2 (Q soft ) and L 2 (Q stiff ), respectively, one has the obvious identities Next, fix self-adjoint in L 2 (Γ) (and, in general, unbounded) operators Λ soft , Λ stiff defined on domains dom Λ stiff(soft) , (for problems considered in the present article one has dom Λ stiff = dom Λ soft = H 1 (Γ), where H 1 (Γ) is the standard Sobolev space pertaining to the smooth boundary Γ).Still following [59], we define the "second boundary operators" Γ stiff and the action of Γ for all f ∈ L 2 (Q stiff(soft) ), φ ∈ dom Λ stiff (soft) .
Under the additional assumption that dom Λ soft = dom Λ stiff =: D, we also introduce the operator Γ 1 as the formal 6 sum of Γ soft 1 and Γ stiff 1 .The rigorous definition of Γ 1 is as follows: dom Γ 1 := dom A 0 ΠD, ( (The PDE interpretation of this definition, which is, in fact the jump of the co-normal derivative on the boundary, is given below on p. 9.) We remark that the operators Γ 1 , Γ soft(stiff) 1 thus defined are not assumed to be either closed or closable.
For purposes of our analysis, we set Λ soft(stiff) to be the DN maps pertaining to the components Q soft and Q stiff , respectively.More precisely, for the problem 7 we define Λ soft as the operator mapping the boundary values φ of u φ to the trace of its co-normal derivative 8 where, as before, n + is the outward unit normal to Γ relative to Q soft .This operator is first defined on sufficiently smooth φ ∈ H.It is a classical result [36,37,65,66] that it is then extended to a self-adjoint operator with domain D = H 1 (Γ) (moreover, it is a pseudo-differential operator of order one).Similarly, on the stiff component Q stiff we consider the problem and define Λ stiff as the operator mapping the boundary values φ of u φ to the ε −2 -weighted trace of its co-normal derivative where, as before, n − is the outward unit normal to Γ relative to Q stiff .In the same way as Λ soft , this is extended to a self-adjoint operator with domain D = H 1 (Γ).Note that the operator Λ stiff thus defined explicitly depends on ε.The only dependence on ε is in the definition of the weighted co-normal derivative and is multiplicative.
In view of what follows, it is important to remark that, in general, the operator Λ := Λ soft + Λ stiff cannot be guaranteed to be self-adjoint on H 1 (Γ), since it is a sum of unbounded self-adjoint operators, albeit defined on the same domain.In our setup, however, the situation is better, as described in the following statement.
Lemma 2.1.For ε < ε 0 , where ε 0 is an independent constant, the operator Proof.Using the argument of [37,Lemma 2], one has that for each τ the operators Λ stiff(soft) are perturbations of the corresponding Dirichlet-to-Neumann maps for τ = 0 by uniformly bounded operators.Therefore, it is sufficient to establish the claim of the lemma at an arbitrary single value of τ.From the definition of the operator Λ stiff one has Λ stiff = ε −2 Λ stiff , where Λ stiff is the "unweighted" DN map of the problem (21), i.e., the operator mapping the boundary values φ of u φ to the traces of its co-normal derivative ∂ (τ ) n u φ , which 6 Note that Γ stiff 1 and Γ soft 1 have different domains. 7Recall that Q soft and Q stiff are viewed as subsets of the torus obtained by identifying opposite sides of the cube Q. 8 This definition is inspired by [59].Note that the operator thus defined is the negative of the classical Dirichlet-to-Neumann map of e.g.[37].
For our choice of Λ soft(stiff) we have D = H 1 (Γ), which allows us to consider Γ 1 on dom A 0 ΠH 1 (Γ), see (19).One then writes [59] the second Green identity where the A is the null extension (see [58]) of the operator A 0 to dom A 0 ran Π.Thus, the triple (H, Γ 0 , Γ 1 ) is closely related to a boundary quasi-triple [6] (see also [5]) for the transmission problem considered; cf.[16,17] and references therein for the way to introduce a "proper" boundary triple of [40,32].Unlike the analysis of [24], [26] and [27], this latter version of the theory does not suit our needs, owing to the PDE setup we are dealing with here.
The calculation of Π * in [59] shows that Π * = Γ 1 A −1 0 and also that Γ 1 introduced above acts as follows: on the right-hand side are defined above as (weighted) co-normal derivatives; we reiterate that the normal vector is always assumed to be external, both for Q stiff and Q soft .The transmission problem at hand therefore (at least, formally) corresponds to the "matching" (or "interface", "transmission") condition Γ 1 u = 0. 2.3.M -function and generalised resolvents.On the basis of the triple (H, Γ 0 , Γ 1 ), we next introduce the corresponding version of the Dirichlet-to-Neumann map, namely the so-called M -function, which maps the boundary data Γ 0 u to Γ 1 u for all u ∈ ker(A − z), cf.(11), (12), and note the its additivity property (Lemma 2.4) with respect to the decomposition of Q into the stiff (Q stiff ) and soft (Q soft ) components.We then recall a suitable version [59] of the Kreȋn formula for the resolvents of operators describing BVPs for the equation Au = zu subject to boundary conditions of the form β 0 Γ 0 u + β 1 Γ 1 u = 0, where β 0 , β 1 are linear operators on the boundary.In the case of the operators A (τ ) ε , see (9), such resolvent expressions are an essential ingredient in the derivation (Section 3) of the corresponding norm-resolvent asymptotics as ε → 0.
The next result of [59] summarises the properties of the M -function that we will need in what follows.
2. M is an analytic function with values in the set of closed operators in H, densely defined on the z-independent domain dom Λ.
4. For u z ∈ ker(A − zI) ∩ {dom A 0 Π dom Λ} the following formula holds: Alongside M (z), we define M stiff (z) and M soft (z) pertaining to stiff (Q stiff ) and soft (Q soft ) components of the composite by the formulae A straightforward application of Proposition 2.3, together with ( 13), ( 16), yields the following statement.
Lemma 2.4.The identity The value of the above lemma is clear: in contrast to A (τ ) ε , which cannot be additively decomposed into "independent" parts pertaining to the soft and stiff components, owing to the transmission interface conditions on the common boundary Γ of the two, the M -function turns out to be additive (cf., e.g., [28], where this additivity was observed and exploited in an independent, but closely related, setting of scattering).
In what follows we will observe that the resolvent (A (τ ) ε − z) −1 can be expressed in terms of M (z) via a version of the celebrated Kreȋn formula, thus reducing the asymptotic analysis of the resolvent to that of the corresponding M -function.
Alongside the transmission problem ( 9), whose boundary conditions can now be (so far, formally) represented as Γ 1 u = 0, u ∈ dom A 0 ran Π, in what follows we consider a wider class of problems, formally given by transmission conditions of the type where β 1 is a bounded operator on H and β 0 is a linear operator defined on dom β 0 ⊃ dom Λ.In general, the operator β 0 Γ 0 + β 1 Γ 1 is not defined on the domain dom A. This problem is taken care of by the following assumption, which we apply throughout: We remark that by Proposition 2.3 the operators β 0 +β 1 M (z) are closable for all z ∈ ρ(A 0 ), and the domains of their closures coincide with dom β 0 + β 1 Λ.
For any f ∈ H and φ ∈ dom Λ, the equality to a bounded operator from H ß to H, for which we keep the notation β 0 Γ 0 + β 1 Γ 1 for the sake of convenience.
In our analysis of the problems (9) we shall make use of the following version of the Kreȋn formula, cf.[32, Prop.2, Sec.2] and references therein.
is the resolvent at the point z of a closed densely defined operator A β0,β1 , with domain and defined by the formula In particular, the resolvent of the (self-adjoint) operator of the transmission problem (9), which corresponds to the choice β 0 = 0, β 1 = I, admits the following characterisation ("Kreȋn formula"): In this case, one clearly has follows from [59, Remark 3.9], which shows that (A 0,I − z) −1 f is the (unique) weak solution of the corresponding BVP.
We remark that the operators β 0 and β 1 above can be assumed z-dependent, as this does not affect the related proofs of [59].In this case, however, the operator function R β0,β1 (z) is shown to be the resolvent of a z-dependent operator family.In what follows, we will face the situation where R β0,β1 (z) is a generalised resolvent in the sense of [54,55]; the operator A β0,β1 then acquires a dependence on z in the description of its domain, cf.[62].

Norm-resolvent asymptotics
In the present section we make use of the Kreȋn formula (30) to obtain a norm-resolvent asymptotics of the family A (τ ) ε .In doing so we compute the asymptotics of M (z) −1 based on a Schur-Frobenius type inversion formula, having first written M (z) as a 2 × 2 operator matrix relative to a carefully chosen orthogonal decomposition of the Hilbert space H.In our analysis of operator matrices, we rely on [69].
3.1.Diagonalisation of the M -operator on the stiff component.Our overall strategy is based on the ideas of [26].The proof of the main result is, however, obtained in a different way, although the rationale (see also [24]) is essentially the same.Towards the end of the section we outline a sketch of an alternative line of reasoning, which mimics the machinery of [26].First, we consider the M -function M stiff (z), and observe that where S stiff z is the solution operator of the problem and Γ stiff 1 is defined similarly to (17), (18), with Λ stiff replaced by Λ stiff = ε 2 Λ stiff , cf. proof of Lemma 2.1, and A stiff 0 replaced by A stiff 0 , the self-adjoint operator of the problem (32) subject to Dirichlet boundary condition φ = 0.It is easily seen that Γ stiff 1 computes traces of the co-normal derivative where n is the external unit normal to Q stiff .(Henceforth we often drop the subscripts + and − in the notation n + and n − whenever the choice is clear from the context.We reiterate that the normal is always chosen in the outward direction to the domain considered.)The M -function M stiff (z) is therefore the classical z-dependent DN map pertaining to the magnetic Laplacian on Q stiff (with periodic boundary conditions on ∂Q, if appropriate).Notice that it does not depend on the parameter ε.The representation (26) applied to M stiff (z) together with (31) implies that with a uniform estimate on the remainder term.The asymptotic expansion above follows from the obvious (due to the Rayleigh principle) estimate (cf.( 15)) The second term of the asymptotic expansion (33) cannot be dropped since Λ stiff is singular at τ = 0 and, by continuity of DN maps (see [36,37]), near-singular in the vicinity of τ = 0.In fact, in the setup of Model II it is singular at all values of quasimomentum τ .It is expected that the term O(ε 2 ) in ( 33) is asymptotically irrelevant from the point of view of the analysis of the resolvent (A The above argument suggests that we need to separate the singular and non-singular parts of Λ stiff .In order to do so, we recall that Λ stiff is a self-adjoint operator in H, with domain H 1 (Γ).Clearly, its spectrum ("Steklov spectrum") consists of the so-called Steklov eigenvalues, defined as values µ such that the problem or equivalently the problem (∇ + iτ ) has a non-trivial solution.(Solutions to (34) are the boundary values of the corresponding solutions to (35).)If Γ is Lipschitz, the Steklov spectrum is discrete, accumulates 10 to −∞, and obeys a Weyl-type asymptotics, see e.g.[39] and references therein.
The least (by absolute value) Steklov eigenvalue µ τ clearly vanishes at τ = 0 (the corresponding solution of ( 35) is a constant).By the analysis of [37, Lemma 7], it is quadratic in τ in the vicinity of τ = 0 in the case of Model I and is identically zero for all τ in the case of Model II.The eigenvalue µ τ is simple, and we henceforth denote the associated normalised eigenfunction of Λ stiff ("Steklov eigenvector") by ψ τ ∈ H. Introduce the corresponding orthogonal H-projection P (τ ) := •, ψ τ H ψ τ and notice that it is a spectral projection relative to Λ stiff .Next, we decompose the boundary space H: ⊥ := I − P (τ ) is the mutually orthogonal projection to P (τ ) , which yields the following matrix representation for Λ stiff (and therefore a similar representation for Λ stiff = ε −2 Λ stiff ): where Λ stiff ⊥ := P (τ ) ⊥ is treated as a self-adjoint operator in P (τ ) ⊥ H.Likewise, a matrix representation can be written for M stiff and M stiff , since they are obtained as bounded additive perturbations of Λ stiff and Λ stiff , respectively, see (33) or Proposition 2.3.

First asymptotic result for
ε .Our strategy is to write the operator M (z) as a block-operator matrix relative to the decomposition (36), followed by an application of the Schur-Frobenius inversion formula, see [69] for the current state of the art in the area.The analysis is governed by the properties of Λ, since by Proposition 2.3 the operator M (z) is its bounded perturbation.Although P (τ ) is not in general a spectral projection of Λ (unless τ = 0), the fact that P (τ ) H is one-dimensional simplifies the analysis.
For evaluating M (z) −1 we use the Schur-Frobenius inversion formula [69, Theorem 2.3.3](cf.[61]) In applying this formula to (38), we use the facts that M (z) is closed and A −1 is bounded uniformly in τ for z ∈ K σ by Proposition 2.3 as the inverse of a (scalar) R-function [56, Section 2.10] (also known as Herglotz, or Nevanlinna, function). 11As for the Schur-Frobenius complement S of the block A, it is closed since D is closed (see below) and EA −1 B is bounded.We next show that S is boundedly invertible, with a uniformly small norm of the inverse.In doing so, we use the fact that by construction there is a bounded operator B on P ⊥ H on the same domain.In addition, B = ε −2 B 1 with an ε-independent boundedly invertible B 1 , so that the uniform bound Since B 1 is boundedly invertible, for each v one can pick u := B −1 It follows that AB −1 1 is bounded, and hence AB −1 ≤ Cε 2 .Therefore, Returning to (39), we obtain with a τ -uniform estimate for the remainder term, where, we recall, A = P (τ ) M (z)P (τ ) .Comparing the result of substituting ( 40) into (30) with (29), where β 0 = P (τ ) ⊥ , β 1 = P (τ ) , we arrive at the following statement.
Theorem 3.1.There exists C > 0 such that for all z ∈ K σ , τ ∈ Q , ε < ε 0 one has Proof.For the resolvent (A 30) is applicable, in which for M (z) −1 we use (40).As for the resolvent is clearly applicable.Moreover, for this choice of β 0 , β 1 , the operator 11 By a suggestion of a referee, we provide an alternative rough estimate, which yields a uniform bound on A −1 , independently of complex analysis and, in particular, analytic properties of Herglotz functions.Note that by Proposition 2.3, one has , and therefore, for some τ -independent constants c 1 , c 2 > 0, where we have used the fact that the operator A 0 is bounded below by a positive constant independent of τ as well as the continuity of Π stiff ψτ H with respect to τ, which in turn is a direct consequence of [36,Proposition 2.2].The required bound for A −1 = (P (τ ) M (z)P (τ ) ) −1 now follows immediately.
It can be shown that the operator A P (τ ) ⊥ ,P (τ ) in the above theorem is self-adjoint.We do not require this fact in the analysis to follow and therefore skip the proof of this statement.We note further that the statement of Theorem 3.1 gives an answer to the problem of homogenisation we set out to analyse.Yet, the result in this form is rather inconvenient, as (P (τ ) M (z)P (τ ) ) −1 is not written in an explicit form; the operator A P (τ ) ⊥ ,P (τ ) , which is the norm-resolvent asymptotics of the family ε , is therefore not easily analysed.We next address this issue.

3.3.
Generalised resolvent on the soft component.First, we rewrite the norm-resolvent asymptotics obtained in Theorem 3.1 in a block-matrix form, this time, however, relative to the decomposition This allows us to express the asymptotics of (A which is derived from Theorem 3.1 as follows.
Remark 3.3.The right-hand side of ( 43) is well defined for all z ∈ K σ by the proof of Theorem 3.1, or, in other words, due to the analytic properties of the operator functions M soft and B (τ ) .
The proof of Lemma 3.2 follows immediately from (30) and Lemma 2.4.
Comparing the statements of Lemma 3.2 and Proposition 2.5, one realises that the generalised resolvent R (τ ) ε is the solution operator of a spectral BVP, with the parameter z appearing not only in the differential equation but also in the boundary condition.Namely, we have the following statement.
The boundary condition (44) can be written in the more conventional form where n is the outward normal to the boundary Γ of Q soft .Equivalently, where A soft −B (τ ) (z),I for any fixed z is the operator in L 2 (Q soft ) defined by Proposition 2.5 relative to the triple (H, Π soft , Λ soft ).This operator is maximal anti-dissipative for z ∈ C + and maximal dissipative for z ∈ C − , see [63].
Proof.Follows immediately from Proposition 2.5 applied to the triple (in the sense of [59]) (H, Π soft , Λ soft ) with the following choice of the operators parameterising boundary conditions: 12 We remark that P (τ ) ⊥ , E = 0, D = I in (39)) with respect to the decomposition H = P (τ ) H ⊕ P (τ ) ⊥ H.
Remark 3.5.The above lemma is a realisation of the corresponding general result of [62].Theorem 3.1 now implies a uniform asymptotics for the generalised resolvents R (τ ) ε as ε → 0.
Theorem 3.6.The operator family R (τ ) ε (z), see (42), admits the following asymptotics in the operator-norm topology for z ∈ K σ and τ ∈ Q : eff (z) is the solution operator for the following spectral BVP on the soft component Q soft : ⊥ − P (τ ) B (τ ) (z)P (τ ) and β 1 = P (τ ) .The boundary condition in (45) can be written in the more conventional form where ⊥ −P (τ ) B (τ ) (z)P (τ ) ,P (τ ) , for any fixed z, is the operator in L 2 (Q soft ) defined by Proposition 2.5 relative to the triple (H, Π soft , Λ soft ), where the term "triple" is understood in the sense of [59].This operator is maximal anti-dissipative for z ∈ C + and maximal dissipative for z ∈ C − , see [63].
Proof.On the one hand, by Theorem 3.1, the resolvent (A and therefore where we have used (16).
On the other hand, by Proposition 2.5 the generalised resolvent of the BVP (45) admits the form by an application of the Schur-Frobenius inversion formula (39).Comparing the right-hand sides of ( 46) and ( 47) completes the proof.
Theorem 3.6 can be further clarified by the following construction.Consider the "truncated" boundary space 13 H := P (τ ) H (note, that in our setup H is one-dimensional).Introduce the truncated τ -harmonic lift on H by Πsoft := Π soft | H and the truncated operator Λsoft := P (τ ) Λ soft | H.We prove the following result. 13In what follows we consistently supply the (finite-dimensional) "truncated" spaces and operators pertaining to them by the breve overscript.
holds, where Ssoft z is the solution operator of the problem and M soft is the M -operator defined in accordance with (25), (27) relative to the triple ( H, Πsoft , Λsoft ).

The "effective" resolvent
eff (z) is represented as the generalised resolvent of the problem 3. The triple ( H, Γsoft 0 , Γsoft 1 ) is the classical boundary triple [40,32] for the operator A max defined by the differential expression −(∇ + iτ ) 2 on the domain dom A max = dom A soft 0 ran Πsoft .Here Γsoft Next, we check that Msoft (z) = P (τ ) M soft (z)P (τ ) | H. Indeed, Proposition 2.3 implies and therefore from which the first claim follows.
We now turn our attention to the operator Γsoft Πsoft H, where we have taken into account that Λsoft is bounded.Its action is set by Using the definitions of Πsoft and Λsoft , one has Γsoft where we have taken into account that φ ∈ H. Comparing the latter expression with the definition of Γ soft 1 , we obtain Γsoft , as required.The triple property now follows from (24) and dimensional considerations. 14Given that H is one-dimensional, one is tempted to refer to the M -matrix here as the Weyl-Titchmarsh m-coefficient.

Simplified asymptotics for (A (τ )
ε − z) −1 .Equipped with Theorems 3.6 and 3.7, we are now set to provide a convenient representation for the asymptotics of (A Theorem 3.9.The resolvent (A (τ ) ε − z) −1 admits the following asymptotics in the uniform operator-norm topology: , where the operator R (τ ) eff (z), which in general depends on ε, has the following representation relative to the decomposition Here K z := Γ soft 0 | Nz , where N z := ran S soft z P (τ ) , z ∈ C ± .Proof.First, we note that since ran(S soft z P (τ ) ) is finite-dimensional (one-dimensional, to be precise), the operator K (τ ) z is well defined as a bounded linear operator from N z to H, where N z is equipped with the standard norm of L 2 (Q soft ).We proceed by representing the operator (A P (τ ) ⊥ ,P (τ ) − z) −1 , see Theorem 3.1, in a block-operator matrix form relative to the orthogonal decomposition ⊥ ,P (τ ) − z) −1 P soft , as established by Theorem 3.6, and write Here in the first equality we use Proposition 2.5 along with (41), in second equality we use the fact that Γ soft 0 S soft z = I, and in the third equality we use Proposition 2.5 again, see also (47).Invoking ( 14) completes the proof in relation to the bottom-left element of (50).
Passing over to the top-right entry in (50), we write , and the claim pertaining to the named entry follows by a virtually unchanged argument.
Finally, for the bottom-right entry in (50) we have which completes the proof.
The representation for R (τ ) eff (z) given by Theorem 3.7 allows us to further simplify the leading-order term of the asymptotic expansion for (A (τ ) ε − z) −1 .Indeed, consider the operator P (τ ) B (τ ) (z)P (τ ) in (48); since B (τ ) = −M stiff by definition, we invoke (33) to obtain ), with a uniform estimate for the remainder term.Here the truncated τ -harmonic lift Πstiff is introduced as Πstiff := Π stiff | H relative to the same truncated boundary space as above, H = P (τ ) H. The truncated DN map Λstiff on H is defined by Λstiff := P (τ ) Λ stiff | H, and Theorem 3.7 allows us to use [26,Theorem 6.3] directly.As a result, we obtain By a classical result of [62] (see also [54,55]), the operator R (τ ) hom (z) is a generalised resolvent, so it defines a z-dependent family of closed densely defined operators in L 2 (Q soft ), which are maximal anti-dissipative for z ∈ C + and maximal dissipative for z ∈ C − .Thus, we have proven the following result.
Theorem 3.10.The resolvent (A (τ ) ε − z) −1 admits the following asymptotics in the uniform operator-norm topology: , where the operator R (τ ) hom (z), which is allowed to depend on ε, has the following representation relative to the decomposition (52) Here K (τ ) z is as defined in Theorem 3.9 and the generalised resolvent R (τ ) hom (z) is defined by (51).The above theorem provides us with the simplest possible leading-order term of the asymptotic expansion for (A (τ ) ε − z) −1 .However, it is not obvious whether this leading-order term R (τ ) hom (z) is a resolvent of some self-adjoint operator in the space L 2 (Q soft ) ⊕ Πstiff H ⊂ H.In the next section we prove that this is the case, but prior to doing so we discuss the relationship between our approach to proving Theorem 3.10 and the approach of [26].
3.5.Discussion of strategies for proving Theorem 3.10.We recall that in [26] we adopted an approach to an ODE problem closely related to the one studied here, whereby we derived a result of the same type as in Theorems 3.9 and 3.10, under the only prerequisite that a statement akin to Theorem 3.6 holds.The question thus naturally arises, whether this strategy could also be applied in the PDE setting, or it is limited to the ODE setup.Here we present a sketch of the argument giving the former answer to this alternative. 15ndeed, assume that one has obtained the asymptotic behaviour of the generalised resolvent R (τ ) ε (z), so that there exists a generalised resolvent R (τ ) eff (z) (which can still be ε-dependent) such that the difference between the two admits a uniform estimate by O(ε r ), r > 0 (cf.Theorem 3.6).Then one can still proceed as in the proof of Theorem 3.9, see also [26,Theorem 7.1], to establish the block-matrix representation for (A ε (z).In order to prove the asymptotic expansion for the former, one then needs to establish the boundedness of a suitably redefined operator K (τ ) z .For this one would need to prove the boundedness of the trace operator on the space S z H, which is obviously impossible, since the mentioned space necessarily contains non-smooth functions from ran Π soft .However, if one redefines K (τ ) z to be the trace operator restricted to , the situation can be partially salvaged.Indeed, either by the main a priori estimate of [60] or by a direct analysis based on the compactness properties of Neumann-to-Dirichlet maps, one then proves that Γ soft 0 | Rz is a bounded operator.This is so since, on the one hand, R z ⊂ H 2 (Q soft ) and, on the other hand, R z ⊂ ker(A soft − z), where, as for the operator A in (24), A soft is the null extension of A soft 0 onto dom Γ soft 1 , whence, by a standard ellipticity argument, the L 2 norm of functions in R z is equivalent to their H 2 norm (uniformly in z ∈ K σ ).The boundedness of the trace operator on H 2 is then invoked to complete the proof.
In pursuing the strategy of [26] in the PDE context, one would also have to establish that K (τ ) z is a bounded operator on ran(R (τ ) eff (z) − (A soft 0 − z) −1 ) as well, in order to permit the estimate sought.For the problem we study here this range is finite-dimensional, but only so for the concrete choice of R (τ ) eff given by Theorem 3.1, which constitutes precisely the part of the argument we are trying to avoid here.Therefore, the condition that ran R (τ ) eff (z) − (A soft 0 − z) −1 is finite-dimensional must appear as an additional assumption.As a useful alternative, it suffices to assume that ran(R Either way, we see that the strategy of [26] is indeed applicable in the PDE setup, albeit it requires an additional assumption 16 on the effective operator R (τ ) eff .This is the reason why in the present paper we have chosen a different approach, despite the obvious appeal of a strategy whereby one only needs to establish an asymptotic expansion for the family of generalised resolvents pertaining to the soft component.
Having both above approaches in mind, one notices that the overall strategy in [26] and the present paper can be seen as rooted in the analysis of the mentioned generalised resolvents.In a nutshell, one can understand this by the following abstract discourse.
3.6.The nature of our strategy in abstract terms.Assuming, for the sake of argument, that R (τ ) ε (z) has a limit, as ε → 0, in the uniform operator topology for z in a domain D ⊂ C, and, further, that the resolvent (A (τ ) ε − z) −1 also admits such limit, one clearly has where R (τ ) 0 and eff are the limits introduced above.The powerful idea of simplifying the required analysis by passing to the resolvent "sandwiched" by orthogonal projections onto a carefully chosen subspace stems from the pioneering work of Lax and Phillips [51], where the resulting sandwiched operator is shown to be the resolvent of a dissipative operator.This idea was later successfully extended to the case of generalised resolvents in [32], as well as in [2] with the scattering theory in mind.
The function R (τ ) 0 defined by ( 53) is a generalised resolvent, whereas A eff is its out-of-space self-adjoint extension (or Strauss dilation [62]).By a theorem of Neumark [54] (cf.[55]) this dilation is defined uniquely up to a unitary transformation of a special form, which leaves the subspace H soft := L 2 (Q soft ) = P soft H intact, provided that the minimality condition holds.This can be reformulated along the following lines: one has minimality, provided that there are no eigenmodes in the effective media modelled by the operator A (τ ) ε , and therefore in the medium modelled by the operator A (τ ) eff as well, such that they "never enter" the soft component of the medium.A quick glance at the setup of critical-contrast composites helps one immediately convince oneself that this must be true.It then follows that the effective medium is completely determined, up to a unitary transformation, by R (τ ) 0 (z).We must admit that the Neumark-Strauss general theory is not directly applicable in our setting.Part of the reason for this is that R (τ ) ε (z) in general does not converge (we have shown it to admit an asymptotic expansion instead).Even in Model II, where one can obtain a limit proper, one still needs to prove that the resolvents (A (τ ) ε − z) −1 converge as well.Therefore, in our analysis we only use the general theory presented above as a guide.In fact, we manage to compute the required asymptotics of the resolvents independently, thus eliminating the non-uniqueness due to the mentioned unitary transformation.4. Effective model for the homogenisation limit 4.1.Self-adjointness of the asymptotics.This part of the paper draws upon our earlier results, see [26] and references therein.Definition 4.1.Consider the Hilbert space H = H soft ⊕ H, where H is an auxiliary Hilbert space.Let Π be a bounded and boundedly invertible operator17 from H to H. We point out that the operator Π may depend on τ, which we conceal in the notation, in line with the convention we adopt for Π stiff(soft) .It is assumed throughout that this dependence of Π on τ is continuous.For each τ ∈ Q , we set where A max is a closed operator in H soft and Γsoft 0 : dom A max → H is a surjective mapping defined in Theorem 3.7.Clearly, dom A (τ ) thus defined is dense in H.We introduce a linear operator A (τ ) on this domain by setting where B (τ ) is assumed to be a bounded operator in H and Γsoft We have the following statement (cf.[26]) for each τ ∈ Q .

Lemma 4.3. The operator
Proof.On the one hand, for all (u, u) ∈ dom A (τ ) and for (v, v) ∈ dom A (τ ) one has where the second Green's identity has been used.
On the other hand, u u , and the claim follows by comparing the right-hand sides of ( 55) and (56).
In fact, A (τ ) is self-adjoint if and only if B (τ ) is self-adjoint.This follows from the next theorem via an explicit construction of the resolvent A (τ ) − z −1 .
Theorem 4.4.Fix τ ∈ Q , and assume that B (τ ) = B (τ ) * .Then A (τ ) is self-adjoint, and its resolvent is defined at all z ∈ C ± by the following expression, relative to the space decomposition H = Here R(z) is a generalised resolvent in H soft defined as the solution operator of the BVP and the bounded rank-one operator K z is defined by K z := Γsoft 0 | Nz , where, as before (cf.Theorem 3.9), N z = ran S soft z P (τ ) .

General homogenisation result.
Here we formulate the ultimate version of our homogenisation result, having collected all the necessary notation and definitions in one place for easy reference.
Compared to the general setup of Section 4.1, here we identify, by means of an obvious unitary transform, the auxiliary Hilbert space H with C 1 .Then the result of Corollary 4.5 together with that of Theorem 3.10 leads to the following explicit description of the homogenised resolvent for (A Let µ τ be the least (by absolute value) Steklov eigenvalue of the restriction of the BVP to the stiff component Q stiff of the medium, see (35).The multiplicity of µ τ is one, and the associated normalised Steklov eigenvector is denoted by ψ τ ∈ H 1 (Γ).The corresponding orthogonal projection is P (τ ) := •, ψ τ H ψ τ and its mutually orthogonal one is P (τ ) Recall that Λ stiff is the DN map pertaining to the stiff component as defined by ( 21)- (22).The τ -harmonic lift operators Π stiff(soft) of the stiff (respectively, soft) component are defined in (10).The "truncated" DN map and τ -harmonic lift operators are Λstiff = P (τ ) Λ stiff | H and Πstiff(soft) = Π stiff(soft) | H. Furthermore, we note, see (37), that Λstiff is the operator of multiplication by ε −2 µ τ .
The operator A max , see Theorem 3.7, is defined on the domain as the null extension, by a one-dimensional subspace, of the Dirichlet magnetic Laplacian A soft 0 given by the differential expression where u| Γ ∈ P (τ ) H is the trace of the function u and j Q is the unitary operator from the space L{Π stiff ψ τ } (supplied with the standard norm of hom the action of the operator is set by Here Π τ is a boundedly invertible operator from P (τ ) H to C 1 , and (20).We remark that, due to the identity Λstiff = Γstiff 1 Πstiff , the second component of the action of the operator A (τ ) hom can be written as where, in the second term in brackets on the right-hand side, ∂ n is the co-normal derivative on Q stiff , see (22).
The main result of the present work, which follows from Corollary 4.5, is as follows.
Theorem 4.6.The resolvent (A (τ ) ε − z) −1 admits the following estimate in the uniform operator norm topology: where Θ is a partial isometry from H onto H hom : on the subspace H soft it coincides with the identity, and each function from The estimate ( 65) is uniform in τ ∈ Q and z ∈ K σ .n Ψ τ , where Ψ τ is the τ -harmonic lift to Q stiff of the function ψ τ , or, equivalently yet again, ε −2 µ τ where µ τ is the least Steklov eigenvalue of the stiff component, see (35).
Clearly, ψ τ = Ψ τ | Γ , so one is only looking to determine µ τ and Ψ τ .The complexity of determining these objects in Models I and II is incomparable.Indeed, in Model II one immediately has µ τ =0 and Thus, Λstiff ≡ 0, cancelling out the second term in the second component of the definition (64).Moreover, the rank-one operator Π τ sends exp(−iτ so that its action is given by where e 1 is the unit vector in C 1 .Thus, and every object in ( 63) and ( 64) has been computed explicitly.
In contrast, for Model I one cannot explicitly compute either µ τ or Ψ τ (in fact, even in the ODE setup of [26] this computation yields highly non-trivial functions of τ ).The best one can hope for is getting hold of asymptotic expansions for both quantities.Indeed, [37] provides a result of this kind, proving that µ τ is quadratic in τ to the leading order and that Ψ τ admits an expansion ) .Therefore, one would ideally wish to have µ τ , ψ τ and Ψ τ in ( 63) and ( 64) replaced by either the leading-order terms of their asymptotic expansions or finite sums of these.It is clear that in this situation one cannot hope for the error estimate of the same order O(ε 2 ) as in the case of Model II, but the first question is whether the result is stable under this asymptotic procedure.The answer to the posed question is not obvious and will be addressed in the next section.

Stability considerations for Model I.
In the present paper we shall only treat the zero-order stability problem: we seek to replace ψ τ and Ψ τ in ( 63) and ( 64) by their zero-order terms in τ , namely, ψ 0 = |Γ| −1/2 1| Γ and Ψ 0 = |Γ| −1/2 1| Q stiff , i.e. the 0-harmonic lift of ψ 0 .Only the setup of Model I will be considered.For the analysis in this setting it suffices to use the results of [37], whereas in the more general cases (non-constant symbol a 2 (•/ε) of the operator A (τ ) ε on the stiff component, or higher-order stability problems), one would require an advanced asymptotic argument, to which the second part [25] of this paper will be devoted.
The main idea behind our construction below is to prove that the resolvent , where the definitions of β 0 and β 1 are based on the projections P (0) and P ⊥ rather than P (τ ) and P ⊥ .An obvious guess would appear to be β 0 = P ⊥ and β 1 = P (0) , however we instead arrive at the choice β 0 = P (0) ⊥ + Λ ∆ and β 1 = P (0) , where Λ ∆ is a Hermitian operator in P (0) H.We shall not prove this asymptotics for all values of quasimomentum τ at once but first for |τ | ≤ Cε 2/3 , and then the result for all τ ∈ Q will follow from Theorem 3.10 by using the triangle inequality, with an error estimate of order O(ε 2/3 ), see Theorem 4.13.Remark 4.7.Asymptotics leading to better estimates of the remainder (all the way up to O(ε 2 )) are available by resorting to higher order stability problems, see [25] for details.
The starting point of our argument is the following statement.
Lemma 4.11.Consider the operator P (0) Λ stiff P (0) + Λ ∆ on P (0) H.For sufficiently small τ, it is a negative Hermitian operator, and the leading order in τ of its eigenvalue is a negative-definite quadratic form Proof.By Lemma 4.10, we have which yields the required asymptotics.The negative-definiteness now follows from [37], where it is proved that the quadratic in τ leading-order term of the Taylor series for µ τ is negative definite 19 .
We have thus constructed an operator in P (0) H that has the same (to the leading order) eigenvalue as the operator P (τ ) Λ stiff P (τ ) in the "shifted" subspace P (τ ) H and is therefore a good approximation of the latter.Theorem 4.12.For z ∈ K σ , one has where β 0 = P (τ ) ⊥ , β 1 = P (τ ) and β 0 = P (0) Proof.For both resolvents we use the representation (29) provided by Proposition 2.5.As the first terms are identical, the only difference is between the functions which we compare next.The expression Q (z) − Q(z), within the choice for β j , β j , j = 0, 1, as in the formulation of the theorem, is given by where in the second equality we use the second Hilbert identity, and in the third one the estimate P (τ ) −P (0) = O(|τ |), see the proof of Lemma 4.8, as well as the fact that by the Schur-Frobenius inversion formula (39), utilising the triangular form of P ⊥ + P (0) M (z) + Λ ∆ , one has = P (0) P (0) We estimate the three terms of the product in (72) one by one.Note first that P (τ ) bounded uniformly in ε and τ, which is seen by writing, see (39), and using the fact that P (τ ) M (z)P (τ ) −1 is uniformly bounded, see the argument following (39).Furthermore, one has where we use ( 26), (33), and the analysis is immediately reduced to Lemma 4.8.The lemma is applied directly to (P (τ ) −P (0) )Λ stiff P (0) , which results in the O(|τ | 3 /ε 2 ) term in (75).As for the term (P (τ ) −P (0) )Λ soft P (0) , we proceed as the proof of Lemma 4.8, obtaining in place of ( 69) where Π soft,0 is the 0-harmonic lift and ψ 0 is the same as above.Applying the operator P (τ ) − P (0) to the expression (76) yields a O(|τ | 2 ) estimate in (75).Finally, we estimate the operator in (73).On the one hand, for a constant c > 0, which we choose below, it is bounded uniformly when |τ | ≤ cε.Indeed, notice first that due to (33) one has P (0) M (z)P (0) + Λ ∆ = P (0) M soft (z)P (0) + P (0) Λ stiff P (0) + Λ ∆ + zP (0) Π * stiff Π stiff P (0) + O(ε 2 ).Furthermore, using the fact that M soft (z) = ( z)(S soft z ) * S soft z , we obtain (cf.( 60)) The above estimate and the continuity of the operator Π stiff with respect to τ (see e.g.[36, Proposition 2.2]), and hence of the of the rank-one operator Π stiff | P (0) H , yield the claim for |τ | ≤ cε.On the other hand, in the case |τ | ≥ cε we argue, similarly to (75), that P (0) M stiff (z)P (0) + Λ ∆ is approximated with P (0) Λ stiff P (0) + Λ ∆ , which in turn is estimated by Lemma 4.11, while P (0) M soft (z)P (0) is uniformly bounded.As a result, for |τ | > cε, we have where O(1)-terms are uniformly bounded in the operator norm independently of ε, τ.Choosing the constant c to be sufficiently large yields an O(ε 2 /|τ | 2 ) estimate for such ε, τ.
Remark 4.14.The result of Theorem 4.13 does not seem to be order-sharp, owing to the fact that in proving it we have resorted to the triangle inequality.Yet the bottom-right element of the matrix in (74) seems to prevent the possibility of obtaining a better estimate than that claimed by Theorem 4.12.The question of whether such an estimate can be obtained remains open, to the best of our knowledge.We hope to return to this matter in [25].
4.5.Relationship with the Birman-Suslina spectral germ.Here we discuss the relationship between our results and the spectral germ, which lies at the centre of the spectral approach to homogenisation [10], [11] in the moderate-contrast case.In particular, we refer the reader to the recent paper [64] (cf.[23] for a close development using a different approach), where the setup of moderate-contrast homogenisation in perforated media is treated in full detail.We shall consider a special case of [64], which is our Model I (possibly with non-constant weight a 2 (•/ε)) with the soft component replaced by voids, supplied with Neumann boundary conditions on Γ.
The following result of [64] then holds: there exists an operator S in L 2 (Q stiff ) given by called the spectral germ, where q τ is a positive definite quadratic form on Q × Q (which we henceforth also refer to as the spectral germ) that governs the asymptotic behaviour of the problem under consideration in the norm-resolvent sense: where is the Gelfand transform of the original operator.We compare q τ with µ * τ • τ of Lemma 4.11, which is a negative definite quadratic form in τ .It can be shown that these forms are related as follows: The proof of this fact in the case of a general smooth weight a 2 (•/ε) will appear in the second part [25] of the present paper.The value of this result is that if the Birman-Suslina spectral germ for the perforated medium (with the soft component deleted) is known via the analysis of [64] (see also references therein), then the matrix µ * is also known, thus making all the objects appearing in the formulation of Theorem 4.6 explicit.Since the expressions for P (0) and Π 0 in this case are straightforward, the norm-resolvent asymptotics of the operator family A (τ ) ε of Model I is essentially reduced to the same question of controlling the spectral germ as in the case of moderate-contrast perforated medium.From this point of view, our result can be seen as a natural, albeit non-trivial, generalisation of the spectral approach of Birman and Suslina.One significant difference, however, is that in our setup only the O(ε 2/3 ) error bound is obtained, compared to O(ε) in the perforated case.Including additional corrector terms, however, allows one to improve the order of convergence, to which the second part [25] of our work is devoted.

Time-dispersive formulations
In this section we continue the study of Models I and II, for which we have constructed the resolvent asymptotics, in view to obtain equivalent time-dispersive formulations on the stiff component of the medium.In order to achieve this, we first introduce the orthogonal projection P of H hom onto H hom H soft , the latter space being C 1 in both models under consideration.Following this, we determine the corresponding Schur-Frobenius complement P(A (τ ) hom − z) −1 P, see [61, p. 416] and the state of the art in [69].
n u = −(∂u/∂n + iτ • nu)| Γ is the co-normal boundary derivative for Q soft and µ * τ • τ is the negative definite quadratic form in Q × Q provided by Lemma 4.11.
The operator hom is thus the leading-order asymptotics as ε → 0, in the norm-resolvent sense, for the family A (τ ) ε ; the corresponding error term is estimated as O(ε 2/3 ) in the H → H operator norm (where H = L 2 (Q)).We remark that in this case the homogenised operator A (τ ) hom still depends on ε.This is to be expected, as so does (under the chosen scaling) the homogenised operator (ε −2 S + 1) −1 of [10], see (78) above.Nevertheless, the setup of moderate-contrast homogenisation of [10] does allow to get rid of the factor ε −2, by employing a non-unitary transform (which is precisely what is done in [10] and follow-up papers).In the case of critical contrast, it proves impossible to "hide" the parameter ε in this way.hom for all τ ∈ Q (cf. the result of [26,Example (0)] and also of [27]; note the form that the Datta-Das Sarma [31] boundary conditions of the mentioned works admit in the PDE setup).Let where u| Γ , κ are as in (80) and The action of the operator is set by n u = −(∂u/∂n + iτ • nu)| Γ is the co-normal derivative for Q soft .We remark that in this case the homogenised operator A (τ ) hom does not depend on ε, and therefore the asymptotics claimed by Theorem 4.6 is in fact the limit proper; the error is estimated as O(ε 2 ).5.3.Asymptotic dispersion relations.Due to similarities in the formulae, we shall consider Models I and II simultaneously.To this end, for τ ∈ Q we set where in Models I and II we have ψ = ψ 0 and ψ = ψ τ , µ * = 0, respectively.The problem of calculating P(A where z is such that the resolvent (A (τ ) hom − z) −1 exists (which is the case, in particular, for z ∈ K σ ).In order to express u from (85) we represent it as a sum of two functions: one of them is a solution to the related inhomogeneous Dirichlet problem, while the other takes care of the boundary condition.More precisely, consider the solution v to the problem or, equivalently, one has u = zβκ −1 (A soft 0 − z) −1 v, where A soft 0 is, as above, the Dirichlet operator in L 2 (Q soft ) associated with the differential expression −(∇ + iτ ) 2 .We now write the "boundary" part of the system (85) as where Thus P(A (τ ) hom − z) −1 P is the operator of multiplication in C 1 by (K(τ, z) − z) −1 .The formula (88) shows, in particular, that the dispersion function K is singular only at eigenvalues of the Dirichlet Laplacian on the soft component.It allows to compute K in terms of the spectral decomposition of A soft 0 , cf. [74].In order to see this, we represent the action of the resolvent (A soft 0 − z) −1 as a series in terms of the normalised eigenfunctions {ϕ where λ (τ ) j , j = 1, 2, 3, . . ., are the corresponding eigenvalues and v is defined in (86).In each case we consider the problem (85), where operator Γ τ depends on the specific example at hand.

Model I.
Here the eigenvalues of the (τ -dependent) Dirichlet operators A soft 0 are independent of τ, and we write λ (τ ) j = λ j , τ ∈ Q , where the corresponding eigenfunctions satisfy ϕ (τ ) By the second Green's identity we calculate where .
Thus, ultimately one has which is treated, for each τ ∈ Q , as a meromorphic function in z.
Remark 5.1.Formula (91) reveals the meaning of the additional assumption,20 imposed both in [41] and [37], that 1| Q soft has to be non-orthogonal to all eigenfunctions of the Dirichlet Laplacian pertaining to the soft component of the medium.We also remark that the limiting spectrum of our Model I, following from Theorem 4.6, coincides with the set of z such that K I (τ, z) = z for some τ ∈ Q , thus allowing to use (92) to derive the results of [41,37].and the second Green's identity, and hence, proceeding as in the proof of Lemma 4.8, cf. in particular the formulae ( 69), (70), see also Remark 4.9, we obtain Thus, one has where Remark 5.3.In the proof of Lemma 5.7 we also use the representation which follows from (89), ( 93), (94).
Model II.Here we have and we set Ψ τ := Π soft ψ τ ≡ v, cf.(86).One has, by the second Green's identity, Therefore, one has Rearranging the terms and using Green's formula, we obtain where, in accordance with our convention, n is the normal to Γ pointing into Q stiff .Expressing V τ in terms of Ψ τ from (98) and using the definition of the M -matrix on the soft component M soft , we obtain The formula (99) yields, in particular, an asymptotic description for small τ (which will be discussed in [25]), showing that in the case of Model II the integrated density of states admits a limit proper as ε → 0, as opposed to the case of Model I, see [37] for details.We thus and arrive at the following statement.
Lemma 5.4.In both Model I and Model II, the action of the Schur-Frobenius complement P(A (τ ) hom − z) −1 P is represented as the operator of multiplication by (K(τ, z) − z) −1 , where the dispersion function K is given by the following formulae: • For Model I, K = K I (τ, z) defined by (92) or equivalently by (95); • For Model II, K = K II (τ, z) defined by (99).

5.4.
Effective macroscopic problems in R d with frequency dispersion.Here we shall interpret the Schur-Frobenius complements constructed in the previous section as a result of applying the Gelfand transform, see Section 2.1, to a homogenised medium in R d .5.4.1.Preparatory material.In order to make our result more useful for applications, we seek to obtain the whole-space setting as the homogenised medium, see [23] and references therein, contrary to the approach of [64] where the homogenised medium is represented as a perforated one, albeit described by an operator with constant (homogenised) symbol away from the Neumann perforations.
To this end, we put t = ε −1 τ (which corresponds to the inverse unitary rescaling to the one of Section 2), then unitarily immerse the L 2 space of functions of t into the L 2 space of functions of t and x corresponding to the stiff component of the original medium, by the formula and write the effective problem (87) in the form The solution operator for (101), namely is the composition of a projection operator in L 2 εQ stiff × ε −1 Q onto constants in x and multiplication by the function K(εt, z) − z −1 , as follows: for all z such that K(εt, z) − z is invertible, in particular, for z ∈ K σ .We map the operator of (102) unitarily to the rank-one operator in L 2 εQ × ε −1 Q ) defined as Note that we have thus used a sequence of two unitary embeddings, where the first one, given by (100), essentially maps the one-dimensional space C back to L 2 (Q stiff ), whereas the second one "eliminates" the inclusions and results in a homogeneous medium.For this reason we present these embeddings separately, instead of combining them into a single unitary operator.
The sought representation in R d of the Schur-Frobenius complement of Lemma 5.4 is obtained by sandwiching the operator (103) with the Gelfand transform G defined in (7) and its inverse G * given by (8), so that the overall operator is given by

This results in the mapping
which yields the required effective problem on R d .In deriving (104) we use the fact that the inner product of the Gelfand transform of a function and the identity coincides with its Fourier transform.Here is the Fourier transform of F. Combining Theorems 4.6 and 4.13, we arrive at the following statement.
Theorem 5.5.The direct integral of Schur-Frobenius complements , (4), and by , and the exponent r is as in Theorem 5.5.Note that if one applies a unitary rescaling to the original medium, so that Ω ε stiff is transformed to Ω 1 stiff := ∪ n∈Z d (Q stiff + n), cf. the setups of [41,37], then the family Θ ε hom is mapped to an ε-independent unitary operator.
We will now discuss the implications of the above result for the effective time-dispersive medium in each of the two models.5.4.2.Model I. Our first result in this direction is the following lemma.Define K I (τ, z) by the formula (cf.(92)) where (cf.( The following estimate allows us to replace K I by K I in the expression for the operator asymptotics we are seeking for Model I. Lemma 5.7.For all F ∈ L 2 (R d ), set (cf. ( 104)) 21 Here the unitary operator corresponding to the equivalence does not depend on ε, see Corollary 5.6. 22An explicit description of the operator Θ ε hom is somewhat ugly and depends on whether one deals with Model I or Model II.In both cases it is based on a unitary operator Vτ which is block-diagonal with respect to the decomposition of L 2 (Q stiff ) into span{ητ } and its orthogonal complement.In the Model I ητ = |Q stiff | −1 1 Q stiff , and in Model II ητ = Π stiff ψτ −1 Π stiff ψτ , cf.Theorem 4.6.The operator Vτ sends ητ to |Q stiff | −1 1 Q stiff and acts as an arbitrary unitary operator on its orthogonal complement.The operator Θ ε hom is then given by Θ ε hom = G −1 Φ * ε Vτ ΦεG stiff , where G is the Gelfand transform of Section 2.1, Φε is the unitary rescaling defined in Section 2.1, Φε is a version of this unitary rescaling obtained by a restriction to L 2 (εQ stiff ), and G stiff is an appropriate version of the Gelfand transform on Ω ε stiff .
Then following estimate holds: Proof.For each z ∈ C \ R, consider the function Using the expressions (95), ( 96), (106), we write where We shall now bound the right-hand side of (109) in terms of ε and express Υ ε z − Υ ε z in terms of D, which will yield the required estimate (108).
Proceeding with this plan, notice that m(•, τ ) is an R-function for each τ ∈ Q , and | m(z, τ )|/| z| is estimated away from zero uniformly in z ∈ C \ R, τ ∈ Q .Indeed, for each τ ∈ Q denote by E τ the spectral measure of the operator A soft 0 and by ν τ the scalar measure on R defined on Borel sets B by the formula Assume first that z > 0. Then by the spectral theorem (see e.g.Combining this with (113) yields (Υ , where the constant C > 0 is independent of F, as claimed. In the next lemma we extend the integration in τ in the formula (107) to the whole space R d , thus eliminating the dependence on ε in the asymptotics.Lemma 5.8.One has the estimate Proof.We follow the strategy of the proof of Lemma 5.7.Notice that for ε > 0, z ∈ K σ , and x ∈ R d one has where Proceeding as in the proof of Lemma 5.7 via the Parceval identity, we conclude that Therefore, the following bound holds: and the claim follows.
Therefore, for all F ∈ L 2 (R d ), one has Taking into account Corollary 5.6, the following result has thus been proved.
Theorem 5.9.In the case of Model I, one has the following estimate: Here Ω ε stiff is the stiff component of the original composite, see ( 1)-( 4), P ε stiff is the orthogonal projection of L 2 (R d ) onto L 2 (Ω ε stiff ), A hom is an elliptic operator in L 2 (R d ) with the constant symbol (120), the function B(z) is defined by (118), and Θ ε hom is a unitary operator from L 2 (R d ) to L 2 (Ω ε stiff ).Remark 5.10.Note that an estimate of this type can be obtained by Theorem 4.6 with a better error than O(ε 2/3 ).In this case, however, the operator A hom may no longer be shown to be elliptic, see, e.g., [24] where it is shown that spatial dispersion (i.e.spatial non-locality) may appear in the description of A hom .5.4.3.Model II.In the case of Model II, there appears to be no obvious way to proceed beyond the statement of Theorem 5.5 in the analysis of (A (τ ) ε − z) −1 sandwiched onto the stiff component of the medium.In particular cases this can be accommodated, however, see, e.g., Examples (0) and (2) of [26], where a related ODE setup leads to a finite-difference frequency-dispersive homogenised operator on the stiff component.This owes to the explicit dependence of K II (τ, z) on τ in models stemming from ODEs.(In a nutshell, the function K II is shown to depend explicitly on cos τ only.) We would like to also point out a relation between the approach and the results of the present paper pertaining to Model II and the setup of [18] (see also references therein), which is, however, beyond the scope of the present paper and shall be scrutinised elsewhere.

Concluding remarks
Here we indicate the applicability of our argument to the more general case, where the symbol a 2 (•/ε) of A ε is no longer constant on the stiff component Q stiff of the medium.It is of course questionable whether this general setup is appealing enough from the point of view of applications in materials science to vindicate its in-depth study, since it requires to select a rather exotic medium with material parameters depending on the spatial variable on the stiff component, while at the same time the major influence on the results is surely exerted by the presence of the soft component, as evident from in particular the results of Section 5 above.
The material of Sections 3 and 4, up to and including Theorem 4.6, is formulated in an abstract way and requires no amendments, provided that the DN maps introduced in Section 2 are redefined accordingly.In particular, the co-normal derivative ∂ n , and provided the asymptotic analysis of the DN map Λ stiff (and in particular of its least eigenvalue) is carried out.We shall address this in the second part of the paper [25].In particular, Theorem 4.13 continues to hold, pending the necessary modification of the one-dimensional operator Λ ∆ .
Finally, the analysis of Section 5 continues to hold with no modifications.
, one has A 0 = A stiff 0 ⊕ A soft 0 relative to the orthogonal decomposition H = L 2 (Q stiff ) ⊕ L 2 (Q soft ), where A stiff 0 and A soft 0 are the operators of the Dirichlet BVP on Q stiff and Q soft , for the differential expressions −ε −2 (∇ + iτ ) 2 and −(∇ + iτ ) 2 , respectively.All three operators A 0 , A soft 0 and A stiff 0 are self-adjoint and positive-definite.Moreover, there exists a bounded inverse A −1 0 for all τ ∈ Q , and dom A 0 ∩ ran Π = {0}.Denoting by Γ stiff(soft) 0the left inverse 5 of Π stiff(soft) , we introduce the trace operator Γ

Lemma 3 . 4 .
The generalised resolvent R (τ ) ε (z) is the solution operator of the following BVP on L 2 (Q soft ):

4. 3 .
Models I and II.The definition of the homogenised operator in the previous section essentially relies on the knowledge of the following objects: (i) the normalised first Steklov eigenvector ψ τ ; (ii) the result of the calculation Π τ ψ τ ; and (iii) the truncated DN map of the stiff component Λstiff or, equivalently, ε −2 ∂ (τ )

5. 2 .
Model II.Asymptotically equivalent setup.Applying Theorem 4.6 together with the results of Section 4.3, we obtain the following explicit description of the homogenised operator A (τ )

Corollary 5 . 6 .
105) is O(ε r )-close in the uniform operator-norm topology, uniformly in z ∈ K σ , to an operator unitary equivalent 21 to the (pseudo-differential) operator Υ ε z .Here r = 2/3, K = K I for Model I and r = 2, K = K II for Model II.The direct integral (105) is the composition of the original resolvent family (A ε −z) −1 applied to functions supported by the stiff component in R d and the orthogonal projection onto the same stiff component.Applying a version of the scaled transform to the restriction of functions in L 2 (R d ) to the stiff component of the original ε medium, we obtain the following result for the operator family A ε .Denote by Ω ε stiff

n
pertaining to the stiff component of the medium Q stiff is defined as∂ (τ ) n u = −a 2 (∇u + iτ u) • n| Γ .The analysis of Section 4.4 stays intact after an obvious modification to the definition of the operator Λ stiff based on the above re-definition of the operator ∂ (τ ) max → H is a surjective mapping defined in Theorem 3.7.Recall that ( H, Γsoft 0 , Γsoft 1 ) is a classical boundary triple for A max , see Theorem 3.7.Remark 4.2.The formula (54) is the PDE version of a formula in Definition 7.3 of [26].