A variational approach to solitary gravity-capillary interfacial waves with infinite depth

We present an existence and stability theory for gravity-capillary solitary waves on the top surface of and interface between two perfect fluids of different densities, the lower one being of infinite depth. Exploiting a classical variational principle, we prove the existence of a minimiser of the wave energy $\mathcal{E}$ subject to the constraint $\mathcal{I}=2\mu$, where $\mathcal{I}$ is the wave momentum and $0<\mu<\mu_0$, where $\mu_0$ is chosen small enough for the validity of our calculations. Since $\mathcal{E}$ and $\mathcal{I}$ are both conserved quantities a standard argument asserts the stability of the set $D_\mu$ of minimisers: solutions starting near $D_\mu$ remain close to $D_\mu$ in a suitably defined energy space over their interval of existence. The solitary waves which we construct are of small amplitude and are to leading order described by the cubic nonlinear Schr\"odinger equation. They exist in a parameter region in which the `slow' branch of the dispersion relation has a strict non-degenerate global minimum and the corresponding nonlinear Schr\"odinger equation is of focussing type. We show that the waves detected by our variational method converge (after an appropriate rescaling) to solutions of the model equation as $\mu \downarrow 0$.


The model
We consider a two-layer perfect fluid with irrotational flow subject to the forces of gravity, surface tension and interfacial tension. The lower layer is assumed to be of infinite depth, while the upper layer has finite asymptotic depth h. We assume that density ρ of the lower fluid is strictly greater than the density ρ of the upper fluid. The layers are separated by a free interface {y = η(x, t)} and the upper one is bounded from above by a free surface {y = h + η(x, t)}. The fluid motion in each layer is described by the incompressible Euler equations. The fluid occupies the domain Σ(η) ∪ Σ(η), where Σ(η) := (x, y) ∈ R 2 : − ∞ < y < η(x, t) , Σ(η) := (x, y) ∈ R 2 : η(x, t) < y < h + η(x, t) , and η = (η, η). Since the flow is assumed to be irrotational in each layer, there exist velocity potentials φ and φ satisfying ∆φ = 0 in Σ(η), ∆φ = 0 in Σ(η).
Since β − (1 + ρ) β < β + (1 + ρ)β, we have that λ ± (k) → ∞ as |k| → ∞. In view of the behaviour at 0, we conclude that λ − (k) is minimised at some k = k 0 > 0. In order to find solitary waves we will assume the following non-degeneracy conditions. The first part of the assumption is introduced in order to avoid resonances. The second part is introduced in order to obtain the inequality (1.14) below. This in turn dictates the choice of model equation (the cubic nonlinear Schrödinger equation). We note that these conditions are satisfied for generic parameter values, but that there are exceptions; see Figures 3 and 4. Set ν 0 = λ − (k 0 ) and note that v 0 = (1, −a) is an eigenvector to the eigenvalue For future use we also introduce the matrix-valued function which satisfies g(k 0 )v 0 = 0 and (due to the second part of Assumption 1.1 and evenness) for ||k| − k 0 | 1, where c is a positive constant. Bifurcations of nonlinear solitary waves are expected whenever the linear group and phase speeds are equal, so that ν (k) = 0 (see Dias & Kharif [12,Section 3]). We therefore expect the existence of small-amplitude solitary waves with speed near ν 0 , bifurcating from a linear periodic wave train with frequency k 0 ν 0 . Making the Ansatz where 'c.c.' denotes the complex conjugate of the preceding quantity, and expanding in powers of µ one obtains the cubic nonlinear Schrödinger equation for the complex amplitude A, in which and A 3 and A 4 are functions of ρ, β and β which are given in Proposition 3.27 and Corollary 3.24. At this level of approximation a standing wave solution to (1.15) of the form A(X, T ) = e iνNLST φ(X) with φ(X) → 0 as X → ±∞ corresponds to a solitary water wave with speed Proof. Let v(k) be a smooth curve of eigenvectors of F (k) −1 P (k) corresponding to the eigenvalue Evaluating the first equation at k = k 0 and using that λ − (k 0 ) = 0, we find that Taking the scalar product of the second equation with v(k), evaluating at k = k 0 and using the previous equality, we therefore find that where we have also used that g(k 0 )v 0 = 0. This concludes the proof since λ − (k 0 ) > 0 and F (k 0 ) and g(k 0 ) are positive definite.
It follows that a necessary and sufficient condition for (1.15) to possess solitary standing waves is that the coefficient in front of the cubic term is negative.
The following lemma gives a variational description of the set of such solutions (see Cazenave [8,Section 8] The set of complex-valued solutions to the ordinary differential equation These functions are precisely the minimisers of the functional E NLS : H 1 (R) → R given by the constant 2ν NLS is the Lagrange multiplier in this constrained variational principle and

Main results
The main result of this paper is an existence theory for small-amplitude solitary-wave solutions to equations (1.1)-(1.8) under Assumptions 1.1 and 1.3. The waves are constructed by minimising the energy functional E subject to the constraint of fixed horizontal momentum I; see Theorem 3.4 for a precise statement. As a consequence of the existence result we also obtain a stability result for the set of minimisers; see Theorem 3.5. Before describing our approach in further detail, we note that the above formulation of the hydrodynamic problem has the disadvantage of being posed in a priori unknown domains. It is therefore convenient to reformulate the problem in terms of the traces of the velocity potentials on the free surface and interface. We denote the boundary values of the velocity potentials by [19] and Benjamin & Bridges [3] (see also [9,10]) we set the natural choice of canonical variables is (η, ξ), where η = (η, η), ξ = (ξ, ξ). We formally define Dirichlet-Neumann operators G(η) and G(η) which map (for a given η) Dirichlet boundary-data of solutions of the Laplace-equation to the Neumann boundary-data, i.e.
; see Section 2 for the rigorous definition. Note that G only depends on η, whereas G depends on η and η. The boundary conditions (1.3)-(1.4) imply that we can recover Φ and Φ from ξ using the formulas under assumption (1.18). Moreover, the total energy and horizontal momentum can be reexpressed as and respectively, where we have abbreviated . (1.23) Note that We now give a brief outline of the variational existence method. We tackle the problem of finding minimisers of E(η, ξ) under the constraint I(η, ξ) = 2µ in two steps.

Minimise
. Because ξ η minimises E(η, ·) over T µ there exists a Lagrange multiplier γ η such that Hence Furthermore we get , where 25) with N (η) := G(η) −1 and see Proposition 2.19. For J µ (η) we obtain the representation We address the problem of minimising J µ using the concentration-compactness method. The main difficulties are that the functional is quasilinear, nonlocal and nonconvex. These difficulties are partly solved by minimising over a bounded set in the function space, but we then have to prevent minimising sequences from converging to the boundary of this set. This is achieved by constructing a suitable test function and a special minimising sequence with good properties using the intuition from the nonlinear Schrödinger equation above.
Our approach is similar to that originally used by Buffoni [4] to study solitary waves with strong surface tension on a single layer of fluid of finite depth, and later extended to deal with weak surface tension [5,6,13], infinite depth [4,14], fully localised three-dimensional waves [7] and constant vorticity [15]. Our main interest is in investigating the nontrivial modifications needed to deal with multi-layer flows. We give detailed explanations when needed (see in particular the discussion of the vector-valued Dirichlet-Neumann operators in the next section) and refer to the above papers for the details of the proofs when possible.
Note that we could also have considered a bottom layer with finite depth. This introduces an additional dimensionless parameter in the problem (the ratio between the depths of the two layers), which allows for other phenomena (for example, the slow speed can have a minimum at the origin). We refer to [25] for a discussion of the dispersion relation and numerical computations of solitary waves in the finite depth case. One of the reasons why we chose to look at the infinite depth problem is that it entails some technical challenges which invalidates the use of certain methods which are widely used to find solitary waves in hydrodynamics. In particular, the idea originally due to Kirchgässner [18] of formulating the steady water wave problem as an ill-posed evolution equation and applying a centre-manifold reduction cannot be used. The variational method that we use is less sensitive to these issues. Note however that Kirchgässner's method has been extended to deal with the issues due to infinite depth by several authors (see [2] and references therein) and this could have been used in order to construct solitary waves also in our setting. These methods give no information about stability, however.
As far as we are aware, there are no previous existence results for solitary waves in our setting. However, Iooss [16] constructed small-amplitude periodic travelling-wave solutions of problem (1.1)-(1.8) in two situations. The first situation is when the parameters are chosen so that ν 2 = λ + (k) or ν 2 = λ − (k) for some wavenumber k = 0 which is not in resonance with any other wavenumber (i.e. λ ± (nk) = ν 2 for all n ∈ Z) and λ ± (k) = 0 (where the sign is chosen such that λ ± (k) = ν 2 ). The second situation is the 1 : 1 resonance, that is when k is a non-degenerate critical point of λ ± . In both situations he proved the existence of small amplitude waves with period close to 2π/k using dynamical systems techniques. The second situation includes our setting, but is somewhat more general (the critical point is e.g. not assumed to be a minimum). There are also a number of papers dealing with solitary or generalised solitary waves (asymptotic to periodic solutions at spatial infinity) in the related settings where either one or both of the surface and interfacial tension vanishes (see [1,2,11,17,23,24] and references therein). The variational method presented in this paper does not work in those settings since it requires both surface tension and interfacial tension. Finally, let us conclude this section by mentioning that our assumptions exclude two possibilities which could be interesting for further study (by variational or other methods), that is when λ − has a degenerate global minimum at k 0 (see Figure 4) or when the minimum value is attained at two distinct wave numbers ( Figure 3). Also, when Assumption 1.1 is satisfied, but the corresponding nonlinear Schrödinger equation is of defocussing type (so that Assumption 1.3 is violated), one would expect the existence of dark solitary waves.

The functional-analytic setting
The goal of this section is to introduce rigorous definitions of the Dirichlet-Neumann operators G(η) and G(η) and their inverses N (η) and N (η), as well as the operators G(η) and K(η).

Lower fluid
In order to define G(η) and N (η), we first introduce suitable function spaces on which these operators are well-defined. We begin by recalling the definition of the Schwartz class S(Ω) for an open set Ω ⊂ R n : Definition 2.1.

(i) LetḢ
The following result is classical and the proof is therefore omitted (we do however present a proof of a similar result for the upper domain later; see Proposition 2.7).

Proposition 2.2.
(i) The trace map u → u| y=η defines a continuous mapḢ 1 (Σ(η)) →Ḣ 1 2 (R) and has a continuous right inverseḢ Using Proposition 2.2 and the definition of G(η), we find that for some constant c > 0 which depends on η W 1,∞ (R) . From this we immediately obtain the following result.
is defined as the inverse of G(η).

Upper fluid
We next discuss the same questions for the upper fluid. Here we have the additional difficulty that both boundaries are free. Choose h 0 ∈ (0, 1). In order to prevent the boundaries from intersecting, we consider the class (iii) Let X be the Hilbert space .
(iv) Let Y be the Hilbert space .

Note that we have the inclusions
The reason for introducing the space X is that it is the natural trace space associated withḢ 1 (Σ(η)). Since this is not completely standard, we include a proof.
Conversely, given This means that u is the element ofḢ 1 (Σ 0 ) whose partial derivatives have Fourier transforms It is clear from these formulas that the map Proposition 2.8. The space Y can be identified with the dual of X using the duality pairing . Definition 2.9. For η ∈ W , the bounded linear operator G(η) : X → Y is defined by where ·, · denotes the Y × X pairing and φ j , j = 1, 2, is the unique function inḢ 1 for all ψ ∈Ḣ 1 (Σ(η)) with ψ| y=η = 0 and ψ| y=1+η = 0 As in the case of the lower fluid, we obtain that for some constant c > 0 which depends on h 0 and η W 1,∞ (R) , and the following consequence.
Definition 2.11. For η ∈ W , the Neumann-Dirichlet operator N (η) : Y → X is defined as the inverse of G(η).

Further operators
We now proceed with the rigorous definition of the operators G(η), N (η) and K(η). Recall that the definition of G(η) involves various combinations of the components of G(η) (cf. (1.23)). We can formally write but since the definition of the function space X involves the condition Φ s − Φ i ∈ H 1 2 (R) which couples the components Φ s and Φ i , the definition of the components G ij requires some care. Note however that (H , by Definition 2.6 and (2.2) with Φ s = 0. It follows that Recall that we formally defined the operator G(η) by .
However, we need to extend it to a larger space in order to define K(η). We record some lemmas which enable us to do this.
Proof. The first part follows from the facts that ). The second part now follows from the fact that Recall that ξ is defined in terms of Φ and Φ through (1.17). Conversely, we can formally recover Φ and Φ from ξ under the assumption (1.18) through (1.20). We now investigate these relations in more detail. We begin defining appropriate function spaces for ξ and G(η)ξ.
(i) LetX be the Hilbert space (ii) LetỸ be the Hilbert space equipped with the inner product .
An argument similar to Proposition 2.8 shows thatỸ is dual toX. and Proof. By definition we have that This defines an element ofḢ 1 2 (R) by Corollary 2.15 and the continuity of Similarly, It is easily seen that all of the involved operators are bounded. The final formula follows by straightforward algebraic manipulations.
Proof. Assume that ξ ∈X. A direct computation then shows that where we have used Lemma 2.17. Similarly, We have to show that the last expression is actually an element ofḢ − 1 2 (R). To see this, we note that by the definition of Y and Definition 2.9. On the other hand The boundedness of G(η) follows from the above formulas and Lemma 2.17.
Proof. We begin by showing that N (η) defines an operatorỸ →X.
We are now finally ready to discuss the operator K(η).
(i) LetX be the Hilbert space equipped with the inner product .
(ii) LetY be the Hilbert space equipped with the inner product .
Note,X andY are each other's duals and that (H Proof. The fact that K(η) is a bounded operator fromY toX follows by noting that ∂ x is an isomorphism fromX toỸ and fromX toY . The lower bound follows by settingξ = (ξ, −ξ) and noting that ).
This also shows that K(η) is an isomorphism.
It will be useful to write K(η) in the form

Analyticity and higher regularity
In this section we discuss the analyticity of the operators K(η) and K(η) as functions of η and η respectively. We also discuss how they act on higher order Sobolev spaces, assuming that η is sufficiently regular. We begin by considering the second operator using the method explained by Groves & Wahlén [15]. First note that N (η) is given by where φ ∈Ḣ 1 (Σ(η)) is a weak solution of the boundary-value problem The upper domain can be treated in a similar way. Set and let F (x, y ) = (x, y + f (x, y )). The function u(x, y) = φ(F (x, y)) then solves the boundary value problem ∇ · ((I + Q)∇u) = 0 0 < y < 1, (I + Q)∇u · (0, 1) = Φ s , y = 1, .
Proceeding as before, we obtain the following result. The next theorem follows from the above lemmas and the definitions of the involved operators. It is also possible to study these operators in spaces with more regularity. A straightforward modification of the techniques in [15,20] results in the following theorem.
) are analytic for each s > 0.

Variational functionals
In this section we study the functionals and As a direct consequence of the above formulas and Theorem 2.25 we obtain the following result.
In particular, this lemma implies that We turn now to the construction of the gradients L (η) and L (η) in L 2 (R) and (L 2 (R)) 2 respectively. The following results are proved using the methods explained in [15, Section 2.2.1].
We also find that where L k (η), k = 2, 3, . . ., are the terms in the power series expansion of L(η) at the origin.
The first few terms in the power series expansion of L are given by and (2.15) The first terms in the power series expansion of K(η) will also be needed later (the corresponding gradients are readily obtained from these expressions): (2.16) Note in particular that where P (k) and F (k) are given by equation (1.9). We end this section by recording some useful inequalities.
Proposition 2.29. The estimates Proof. The first estimate is immediate from the form of K(η). The estimates for L(η) follow directly from Proposition 2.21, while those for L 2 (η) follow from (2.17).

Existence and stability
This section contains the main results of the paper. We begin by proving that the functional J µ has a minimiser in U \{0}. This is done by using concentration-compactness and penalisation methods as in [4,5,6,7,14,15] and we refer to those papers for the details of some of the proofs. The outcome is the following result. (i) The set C µ of minimisers of J µ over U \{0} is non-empty.
(ii) Suppose that {η n } is a minimising sequence for J µ on U \{0} which satisfies There exists a sequence {x n } ⊂ R with the property that a subsequence of {η n (x n + ·)} converges in (H r (R)) 2 , 0 ≤ r < 2 to a function η ∈ C µ .
The first statement of the theorem is a consequence of the second statement, once the existence of a minimising sequence satisfying (3.1) has been established. The existence of such a sequence can be proved using a penalisation method [4,7,14,15]. A key part of the proof is the existence of a suitable 'test function' η which satisfies the inequality This implies in particular that any minimising sequence {η n } satisfies this property for n sufficiently large. We construct such a test function in the appendix. Once the existence of the test function has been proved, the remaining steps in the construction of the special minimising sequence satisfying (3.1) are similar to [4,7,14,15], to which we refer for further details. In fact, this special minimising sequence satisfies further properties which will be used below (note that a general minimising sequence satisfies the weaker estimate η n 2 1 ≤ cµ by Proposition 2.29). THEOREM 3.2. Suppose that Assumptions 1.1 and 1.3 hold. There exists a minimising sequence {η n } for J µ over U \{0} with the properties that η n 2 2 ≤ cµ and J µ (η n ) < 2ν 0 µ − cµ 3 for each n ∈ N, and lim n→∞ J µ (η n ) 0 = 0.
The second statement of Theorem 3.1 is proved by applying the concentration-compactness principle (Lions [21,22]) (a form suitable for the present situation can be found in [15,Theorem 3.7]) to a minimising sequence satisfying (3.1). The key step is to show that the function is strictly sub-additive.  (i) The set D µ of minimisers of E over the set is non-empty.
(ii) Suppose that {(η n , ξ n )} ⊂ S µ is a minimising sequence for E with the property that sup n∈N η n 2 < M.
There exists a sequence {x n } ⊂ R with the property that a subsequence of {η n (x n + ·), ξ n (x n + ·)} converges in (H r (R)) 2 ×X, 0 ≤ r < 2 to a function in D µ .
We obtain a stability result as a corollary of Theorem 3.4 using the argument given by Buffoni [4,Theorem 19]. Recall that the usual informal interpretation of the statement that a set V of solutions to an initial-value problem is 'stable' is that a solution which begins close to a solution in V remains close to a solution in V at all subsequent times. The precise meaning of a solution in the theorem below is irrelevant, as long as it conserves the functionals E and I over some time interval [0, T ] with T > 0.
and sup Choose r ∈ [0, 2), and let 'dist' denote the distance in (H r (R)) 2 ×X. For each ε > 0 there exists This result is a statement of the conditional, energetic stability of the set D µ . Here energetic refers to the fact that the distance in the statement of stability is measured in the 'energy space' (H r (R)) 2 ×X, while conditional alludes to the well-posedness issue. At present there is no global well-posedness theory for interfacial water waves (although there is a large and growing body of literature concerning well-posedness issues for water-wave problems in general). The solution t → (η(t), ξ(t)) may exist in a smaller space over the interval [0, T ], at each instant of which it remains close (in energy space) to a solution in D µ . Furthermore, Theorem 3.5 is a statement of the stability of the set of constrained minimisers D µ ; establishing the uniqueness of the constrained minimiser would imply that D µ consists of translations of a single solution, so that the statement that D µ is stable is equivalent to classical orbital stability of this unique solution.
Finally, we can also confirm the heuristic argument given in Section 1.2.
. Furthermore, the speed ν µ of the corresponding solitary wave satisfies Note in particular that since v 0 = (1, −a) with a > 0 (cf. eq. (1.12)) the surface profile η is to leading order a scaled and inverted copy of the interface profile η (cf. Figure 1). The fact that we don't know if the minimiser is unique up to translations is reflected by the lack of control over ω; for the model equation, the minimiser is in fact not unique up to translations (see Lemma 1.4). Using dynamical systems methods (see e.g. [2]), we expect that one can prove the existence of two solutions corresponding to ω = 0 and ω = π above, but without any knowledge of stability. Since the proof of Theorem 3.6 follows [15, Section 5.2] closely, we shall omit it.
The goal of the rest of this section is to prove Theorem 3.3, which follows directly from the strict sub-homogeneity of I µ (see Corollary 3.32). This property is established by considering a 'near minimiser' of J µ over U \{0}, that is a function in U \{0} with for some N ≥ 3. Hence we have L(η), L 2 (η) > cµ (by Proposition 2.29 and the inequality µ 2 /L(η) < 2ν 0 µ) and can identify the dominant term in the 'nonlinear' part of J µ (η). The existence of near minimisers is a consequence of Theorem 3.2. Note that we will work under Assumptions 1.1 and 1.3 throughout the rest of the section, without explicitly mentioning when they are needed. One of the main tools that we will use is the weighted norm and a splitting of η in view of the expected frequency distribution. In fact we split each η ∈ U into the sum of a function η 1 with spectrum near k = ±k 0 and a function η 2 whose spectrum is bounded away from these points. To this end we write the equation where g(k) is given by (1.13). We decompose it into two coupled equations by defining η 2 ∈ (H 2 (R)) 2 by the formula and η 1 ∈ (H 2 (R)) 2 by η 1 = η − η 2 , so thatη 1 has support in S : Here we have used the fact that is a bounded linear operator (L 2 (R)) 2 → (H 2 (R)) 2 .
It will also be useful to express vectors w = (w, w) in the basis {v 0 , v 0 }, where v 0 is the zero eigenvector of the matrix g(k 0 ) (see Section 1.2) and v 0 v. The exact choice of the complementary vector v 0 is unimportant, but in order to simplify the notation later on we choose v 0 = (0, 1). This implies that where c 1 = w and c 2 = w + aw.
The following proposition is an immediate consequence of the definition of η 1 .
As a consequence, η 1 satisfies the equation where In keeping with equation (3.2) we write the equation for η 2 in the form where the decomposition η = η 1 − H(η) + η 3 forms the basis of the calculations presented below. An estimate on the size of H(η) is obtained from (3.4) and Proposition 3.11.
Proposition 3.14. The estimate holds for each η ∈ U .
The above results may be used to derive estimates for the gradients of the cubic parts of the functionals which are used in the analysis below.
Proof. Observe that and estimate the right-hand side of this equation using Propositions 3.11 and 3.14.
An estimate for L 3 (η) is obtained in a similar fashion using Propositions 3.11, 3.13, and 3.14.
Lemma 3.17. Any near minimiserη satisfies the estimates We now have all the ingredients necessary to estimate the wave speed and the quantity |||η 1 ||| α .
Proposition 3.18. Any near minimiserη satisfies the estimates Proof. Combining Lemma 3.12, inequality (3.5) and Lemma 3.17, one finds that 3 2 , from which the given estimates follow by Proposition 3.8.
Proposition 3.22. Any near minimiserη satisfies the estimates Proposition 3.23. Any near minimiserη satisfies the estimates Corollary 3.24. Any near minimiserη satisfies the estimate where . We now turn to the corresponding result for L 3 (η). The following result is obtained by writing expanding the right hand side and estimating the terms using Propositions 3.7 and 3.11, Lemma 3.19 and the identity n(η 1 ,η 1 ,η 1 ) = 0.
Proposition 3.25. Any near minimiserη satisfies the estimate Proposition 3.26. Any near minimiserη satisfies the estimate Proof. Noting that (see Propositions 3.7 and 3.10, Corollary 3.18 and Lemma 3.19) one finds that recalling the definition of H in (3.4). The proof is concluded by estimating (cf. Propositions 3.7 and 3.11, and Lemma 3.19).
Combining Propositions 3.25 and 3.26, one finds that Expanding the right hand side using Lemma 3.20 we then obtain the following result.
Proof. Lemma 3.12 asserts that uniformly over s ∈ [1,2]. The first result follows by estimating (by Propositions 3.22, 3.23 and 3.27) and The second result is derived in a similar fashion.
Lemma 3.29. Any near minimiser satifies the inequality Proof. Note first that an arbitrary function η ∈ U \ {0} satisfies the inequality where we have used that 1.14)). The result now follows from the calculation Proof. The estimates follow by combining Corollary 3.24, Proposition 3.27 and Lemma 3.28, while the inequality forη 1 is a consequence of the first estimate (with s = 1) and Lemma 3.29.
Proposition 3.31. There exists s 0 ∈ (1, 2] and q > 2 with the property that the function is decreasing and strictly negative.
Proof. This result follows from the calculation where we have used Corollary 3.30 and chosen s 0 > 1 and q 0 > 2 so that (3 − q)A 3 + s(4 − q)A 4 , which is negative for s = 1 and q = 2 (by Assumption 1.3), is also negative for s ∈ (1, s 0 ] and q ∈ (2, q 0 ].
The sub-homogeneity of I µ now follows (see Groves

A Test function
In order to show that C µ is non-empty we have to construct a special test-function. Here the eigenvector v 0 = (1, −a) to the eigenvalue λ − (k 0 ) of the matrix F (k 0 ) −1 P (k 0 ) plays an important role (see Section 1.2).