Valentin Samoyeau
email: valentin.samoyeau@univ-nantes.fr

STRICHARTZ ESTIMATES WITH LOSS OF DERIVATIVES UNDER A WEAK DISPERSION PROPERTY FOR THE WAVE OPERATOR

Keywords: April 12, 2016. 2010 Mathematics Subject Classification. Primary, 35B30, 42B37, 47D03, 47D06 dispersive inequalities, Strichartz inequalities, heat semigroup, Schrödinger group, wave operator

Strichartz estimates with loss of derivatives under a weak dispersion property for the wave operator

Valentin Samoyeau

Introduction

The family of so-called Strichartz estimates is a powerful tool to study nonlinear Schrödinger equations. Those estimates give a control of the size of the solution to a linear problem in term of the size of the initial data. The "size" notion is usually given by a suitable functional space L p t L q x . Such inequalities were first introduced by Strichartz in [START_REF] Strichartz | Restrictions of Fourier transforms to quadratic surfaces and decay of solutions of wave equations[END_REF] for Schrödinger waves on the Euclidean space. They were then extended by Ginibre and Velo in [START_REF] Ginibre | Smoothing properties and retarded estimates for some dispersive evolution equations[END_REF] (and the endpoint is due to Keel and Tao in [START_REF] Keel | Endpoint Strichartz estimates[END_REF]) for the propagator operator associated with the linear Schrödinger equation in R d . So for an initial data u 0 , we are interested in controlling u(t, . ) = e it∆ u 0 which is the solution of the linear Schrödinger equation:

i∂ t u + ∆u = 0 u |t=0 = u 0 .
It is well-known that the unitary group e it∆ satisfies the following inequality:

e it∆ u 0 L p L q ([-T,T ]×R d ) ≤ C T u 0 L 2 (R d )
for every pair (p, q) of admissible exponents which means : 2 ≤ p, q ≤ ∞, (p, q, d) = (2, ∞, 2), and

(1.1)

2 p + d q = d 2 .
The Strichartz estimates can be deduced via a T T * argument from the dispersive estimates

(1.2) e it∆ u 0 L ∞ (R d ) |t| -d 2 u 0 L 1 (R d ) .
If sup T >0 C T < +∞, we will say that a global-in-time Strichartz estimate holds. Such a globalin-time estimate has been proved by Strichartz for the flat Laplacian on R d while the local-in-time estimate is known in several geometric situations where the manifold is non-trapping (asymptotically Euclidean, conic, or hyperbolic, Heisenberg group); see [BT07, Bou11, HTW06, ST02, BGX00] or with variable coefficients [START_REF] Robbiano | Strichartz estimates for Schrödinger equations with variable coefficients[END_REF][START_REF] Tataru | Outgoing parametrices and global Strichartz estimates for Schrödinger equations with variable coefficients[END_REF].

The situation for compact manifolds presents a new difficulty, since considering the constant initial data on the torus u 0 = 1 ∈ L 2 (T) yields a contradiction in (1.2) for large time.

Burq, Gérard, and Tzvetkov [START_REF] Burq | Strichartz inequalities and the nonlinear Schrödinger equation on compact manifolds[END_REF] and Staffilani and Tataru [START_REF] Staffilani | Strichartz estimates for a Schrödinger operator with nonsmooth coefficients[END_REF] proved that Strichartz estimates hold on a compact manifold M for finite time if one considers regular data u 0 ∈ W 1/p,2 (M). Those are called "with a loss of derivatives": e it∆ u 0 L p L q ([-T,T ]×M) ≤ C T u 0 W 1/p,2 (M) .

An interesting problem is to determine for specific situations, which loss of derivatives is optimal (for example the work of Bourgain [START_REF] Bourgain | Fourier transform restriction phenomena for certain lattice subsets and applications to nonlinear evolution equations[END_REF] on the flat torus and [START_REF] Takaoka | On 2D nonlinear Schrödinger equations with data on R × T[END_REF] of Takaoka and Tzvetkov). For instance, the loss of 1 p derivatives in [START_REF] Burq | Strichartz inequalities and the nonlinear Schrödinger equation on compact manifolds[END_REF] is shown to be optimal in the case of the sphere. An important remark is that, by Sobolev embedding, the loss of 2/p derivatives is straightforward. Indeed, by Sobolev embedding, we have W 2 p ,2 ֒→ L q since 2 p -d 2 = 0 -d q so that (1.3)

e it∆ u 0 L q e it∆ u 0 W 2 p ,2 ≤ u 0 W 2 p ,2
and taking the L p ([-T, T ]) norm yields e it∆ u 0 L p L q ≤ C T u 0 W 2 p ,2 . Therefore Strichartz estimates with loss of derivatives are interesting for a loss smaller than 2/p.

Let us now set the general framework of our study. We consider (X, d, µ) a metric measured space equipped with a nonnegative σ-finite Borel measure µ. We assume moreover that µ is Alfhor regular, that is there exist a dimension d, and two absolute positive constants c and C such that for all x ∈ X and r > 0

(1.4) cr d ≤ µ(B(x, r)) ≤ Cr d ,
where B(x, r) denote the open ball with center x ∈ X and radius 0 < r < diam(X). Thus we aim our results to apply in numerous cases of metric spaces such as open subsets of R d , smooth d-manifolds, some fractal sets, Lie groups, Heisenberg group, . . . Keeping in mind the canonical example of the Laplacian operator in R d : ∆ = 1≤j≤d ∂ 2 j , we will be more general in the following sense: we consider a nonnegative, self-adjoint operator H on L 2 = L 2 (X, µ) densely defined, which means that its domain

D(H) := {f ∈ L 2 , Hf ∈ L 2 } is supposed to be dense in L 2 .
One of the motivation of our paper is to study the connection between the wave equation and the Schrödinger equation. We define the wave propagator cos(t √ H) as follows: for any f ∈ L 2 , u(t, . ) := t → cos(t √ H)f is the unique solution of the linear wave problem

(1.5)      ∂ 2 t u + Hu = 0 u |t=0 = f ∂ t u |t=0 = 0.
One can find the explicit solutions of this problem in [START_REF] Folland | Introduction to partial differential equations[END_REF] for the Euclidean case and in [START_REF] Bérard | On the wave equation on a compact Riemannian manifold without conjugate points[END_REF] for the Riemannian manifold case through precise formula for the kernel of the wave propagator. Up to our knowledge, those explicit solutions are not available in our abstract setting. It would be of great interest to be able to compute exact expression of the solution of the wave equation in such a general case. The remarkable property of this operator comes from its finite speed propagation: for any two disjoint open subsets U 1 , U 2 ⊂ X, and any functions

f i ∈ L 2 (U i ), i = 1, 2, then (1.6) cos(t √ H)f 1 , f 2 = 0 for all 0 < t < d(U 1 , U 2 ). If cos(t √ H
) is an integral operator with kernel K t , then (1.6) simply means that K t is supported in the "light cone" D t := {(x, y) ∈ X 2 , d(x, y) ≤ t}. We assume that H satisfies (1.6). In [START_REF] Coulhon | Gaussian heat kernel upper bounds via the Phragmén-Lindelöf theorem[END_REF], Coulhon and Sikora proved that this property is equivalent to the Davies-Gaffney estimates

(1.7) e -tH L 2 (E)→L 2 (F ) e -d(E,F ) 2 4t
for any two subsets E and F of X, and t > 0.

It is known that -H is the generator of a L 2 -holomorphic semigroup (e -tH ) t≥0 (see [START_REF] Davies | Non-Gaussian aspects of heat kernel behaviour[END_REF]). We will also assume that the heat semigroup (e -tH ) t≥0 satisfies the typical upper estimates (for a second order operator): for every t > 0 the operator e -tH admits a kernel p t with (DUE)

|p t (x, x)| 1 µ(B(x, √ t)) , ∀ t > 0, a.e. x ∈ X.
It is well-known that such on-diagonal pointwise estimates self-improve into the full pointwise Gaussian estimates (see [Gri97, Theorem 1.1] or [CS08, Section 4.2] e.g.)

(UE) |p t (x, y)| 1 µ(B(x, √ t)) exp -c d(x, y) 2 t , ∀ t > 0, a.e. x, y ∈ X.
One can find in [START_REF] Bernicot | Dispersive estimates with loss of derivatives via the heat semigroup and the wave operator[END_REF] and the references therein some examples where the previous estimates hold.

When dealing with Schrödinger equation on a manifold or a more general metric space, the L 1 -L ∞ estimate (1.2) seems out of reach. In [START_REF] Bernicot | Dispersive estimates with loss of derivatives via the heat semigroup and the wave operator[END_REF], the authors show how to replace it by a H 1 -BMO estimate, with the Hardy space H 1 and the Bounded Mean Oscillations space BMO both adapted to the semigroup. We do not recall the definition of those spaces here, but refer to [START_REF] Bernicot | Dispersive estimates with loss of derivatives via the heat semigroup and the wave operator[END_REF] for more details.

For any integer m ≥ 0 and x ∈ R + we set ψ m (x) = x m e -x . It forms a family of smooth functions that vanish at 0 (except when m = 0) and infinity, which allows us to consider a smooth partition of unity, using holomorphic functionnal calculus (and requiring C ∞ 0 -calculus). The main assumption of our work is the following Assumption 1.1. There exist κ ∈ (0, ∞] and an integer m such that for every s ∈ (0, κ) the wave propagator cos(s √ H) at time s satisfies the following dispersion property

(1.8) cos(s √ H)ψ m (r 2 H) L 2 (B)→L 2 ( B) r r + s d-1 2 ,
for any two balls B, B of radius r > 0.

This estimate is microlocalized in the physical space due to the balls B and B at scale r and in frequency at scale 1 r through ψ m (r 2 H), thus respecting the Heisenberg uncertainty principle. The parameter κ is linked to the geometry of the space X (its injectivity radius for instance).

In the Euclidean space R d , the L 2 (B) -L 2 ( B) dispersion phenomenon seems only to depend on the distance d(B, B). Indeed, the intuition is that, in an isotropic medium a wave propagates the same way in all the directions. That is what leads us to think that Assumption 1.1 could be proved without using a pointwise explicit formula of its kernel, but with a more general approach, using functional tools only, that could be extend to other settings. To our knowledge the study of such behavior is not known and could be a good direction to investigate.

We mentioned that the finite speed propagation property (1.6) gives us the idea that after time s the solution to the wave problem (1.5) with initial data supported in a ball of radius r is supported in a ball of radius r + s. Given that r ≤ s (otherwise L 2 functional calculus yields Assumption 1.1) and the fact that waves propagate the same way in all directions in an isotropic medium, if we cover the sphere of radius r + s by N ≃ r+s is the natural dispersion one can hope for such waves, if we look for a uniform estimate (depending only on r, s).

Indeed we also emphasize that Assumption 1.1 is weaker than the one in [START_REF] Bernicot | Dispersive estimates with loss of derivatives via the heat semigroup and the wave operator[END_REF], namely: There exist κ ∈ (0, ∞] and an integer m such that for every s ∈ (0, κ) we have

(1.9) cos(s √ H)ψ m (r 2 H) L 2 (B)→L 2 ( B) r r + s d-1 2 r r + |s -L| d+1 2
, where L = d(B, B), which describes more precisely the dispersion inside the light cone. However (1.9) can be difficult to prove in an abstract setting. That is why we are interested in proving what Assumption 1.1 could imply as far as Strichartz estimates are concerned. Estimate (1.8) should indeed be much easier to prove in concrete examples. For more on Assumption 1.1, see Subsection 2.4. Consequently the results we obtain will be weaker too. We recall Theorem 1.3 from [START_REF] Bernicot | Dispersive estimates with loss of derivatives via the heat semigroup and the wave operator[END_REF] in order to compare it with our Theorem 1.3.

Theorem 1.2 ([BS14]). Suppose (1.4) with d > 1, (DUE) and Assumption (1.9) with κ ∈ (0, ∞].
Two cases occur:

• if κ = ∞ then we have Strichartz estimates without loss of derivatives;

• if κ < ∞ then for every ε > 0, every 0 < h ≤ 1 with h 2 ≤ |t| < h 1+ε we have Strichartz estimates with loss of 1+ε p derivatives.

To prove this, the authors first reduced the H 1 -BMO estimation to a microlocalized L 2 -L 2 estimate, and then showed how dispersion for the wave propagator implies dispersion for the Schrödinger group. Theorem 3.1 is playing that role in the present paper.

Our main theorem follows the routine of [START_REF] Bernicot | Dispersive estimates with loss of derivatives via the heat semigroup and the wave operator[END_REF] to deduce Strichartz inequalities from L 2 -L 2 estimates.

Theorem 1.3. Assume (1.4) with d > 2, (DUE), and Assumption 1.1. Then for every 2 ≤ p ≤ +∞ and 2 ≤ q < +∞ satisfying

2 p + d -2 q = d -2 2 ,
and every solution u(t, . ) = e itH u 0 of the problem

i∂ t u + Hu = 0 u |t=0 = u 0 , we have • if κ = ∞, then u satisfies local-in-time Strichartz estimates with loss of derivatives (1.10) u L p ([-1,1],L q ) u 0 W 2( 1 2 -1 q ),2 ; • if κ < ∞, then for every ε > 0, u satisfies local-in-time Strichartz estimates with loss of derivatives (1.11) u L p ([-1,1],L q ) u 0 W 1+ε p +2( 1 2 -1 q ),2
. We would like to point out that the straightforward loss of derivatives given by Sobolev embeddings when

2 p + d -2 q = d -2 2 is 2 p + 1 - 2 q .
Thus the loss is nontrivial here. For more on the loss of derivatives, see Remarks 3.7 and 3.8. It is interesting to see how a weak dispersion property on the wave propagator implies dispersion for the Schrödinger operator.

The idea of the proof here is similar to the one in [START_REF] Bernicot | Dispersive estimates with loss of derivatives via the heat semigroup and the wave operator[END_REF]. More particularly it is due to a precise tracking of the constants in some key estimations (from [START_REF] Keel | Endpoint Strichartz estimates[END_REF] for instance). The aim of this paper is to give a better understanding of how dispersion for the wave propagator implies dispersion for the Schrödinger equation, and what Strichartz inequalities ensue in some contexts, where we do not have precise dispersive estimates on the wave propagators. In other words if one can compute, even inaccurate, information about the wave propagator in general settings, it would allow to have some knowledge of the Schrödinger equation in that framework.

The organization of the paper is as follow: In Section 2 we set the notations used throughout the paper and recall some preliminary facts concerning the semigroup, Hardy and BMO spaces, as well as some motivations of our hypothesis. Then Section 3 is dedicated to the proofs of the Theorems.

Preliminaries

2.1. Notations. We denote diam(X) := sup x,y∈X d(x, y) the diameter of a metric space X. For B(x, r) a ball (x ∈ X and r > 0) and any parameter λ > 0, we denote λB(x, r) := B(x, λr) the dilated and concentric ball. As a consequence of (1.4), a ball B(x, λr) can be covered by Cλ d balls of radius r, uniformly in x ∈ X, > 0 and λ > 1 (C is a constant only depending on the ambient space). If no confusion arises, we will note L p instead of L p (X, µ) for p ∈ [1, +∞]. For s > 0 and p ∈ [1, +∞], we denote by W s,p the Sobolev space of order s based on L p , equipped with the norm

f W s,p := (1 + H) s 2 f L p .
We will use u v to say that there exists a constant C (independent of the important parameters) such that u ≤ Cv and u ≃ v to say that both u v and v u. If Ω is a set, 1 Ω is the characteristic function of Ω, defined by

1 Ω (x) = 1 if x ∈ Ω 0 if x / ∈ Ω.
Throughout the paper, unless something else is explicitly mentioned, we assume that d > 2 and that (1.4), (DUE), (1.7), and Assumption 1.1 are satisfied.

2.2. The heat semigroup and associated functional calculus. We consider a nonnegative, self-adjoint operator H on L 2 = L 2 (X, µ) densely defined. We recall the bounded functional calculus theorem from [RS72]:

Theorem 2.1. H admits a L ∞ -functional calculus on L 2 : if f ∈ L ∞ (R + )
, then we may consider the operator f (H) as a L 2 -bounded operator and

f (H) L 2 →L 2 ≤ f L ∞ .
From the Gaussian estimates of the heat kernel (UE) and the analyticity of the semigroup (see [START_REF] Carron | Gaussian estimates and L p -boundedness of Riesz means[END_REF]) it comes that for every integer m ∈ N and every t > 0, the operator ψ m (tH) has a kernel p m,t also satisfying upper Gaussian estimates:

(2.1) |p m,t (x, y)| 1 µ(B(x, √ t)) exp -c d(x, y) 2 t , ∀ t > 0, a.e. x, y ∈ X.
We now give some basic results about the heat semigroup thanks to our assumptions. The detailed proofs can be found in Section 2 of [START_REF] Bernicot | Dispersive estimates with loss of derivatives via the heat semigroup and the wave operator[END_REF].

Proposition 2.2. Under (1.4) and (DUE), the heat semigroup is uniformly bounded in every L p -spaces for p ∈ [1, ∞]; more precisely for every f ∈ L p , we have

sup t>0 e -tH f L p f L p .
Moreover, for m ∈ N and t > 0, since ψ m (tH) also satisfies (DUE) we have

sup t>0 ψ m (tH) L p →L p 1.
Let us now define some tools for the Littlewood-Paley theory we need in the sequel. For all λ > 0 we set

ϕ(λ) := +∞ λ ψ m (u) du u , φ(λ) := λ 0 ψ m (v) dv v = 1 0 ψ m (λu) du u .
Remark 2.3. Notice that ϕ is, by integration by parts, a finite linear combination of functions ψ k for k ∈ {0, .., m}. Moreover for every λ > 0,

φ(λ) + ϕ(λ) = +∞ 0 u m-1 e -u du = Γ(m) = constant.
The following theorem will be useful to estimate the L p -norm through the heat semigroup:

Theorem 2.4. Assume (1.4) and (DUE). For every integer m ≥ 1 and all p ∈ (1, ∞), we have

f L p ≃ ϕ(H)f L p + 1 0 |ψ m (uH)f | 2 du u 1 2 L p . So if q ≥ 2 f L q ϕ(H)f L q + 1 0 ψ m (uH)f 2 L q du u 1 2
.

Such a result can be seen as a semigroup version of the Littlewood-Paley characterization of Lebesgue spaces. A proof of this theorem can be found in [START_REF] Bernicot | Dispersive estimates with loss of derivatives via the heat semigroup and the wave operator[END_REF] (look for Theorem 2.8 in [START_REF] Bernicot | Dispersive estimates with loss of derivatives via the heat semigroup and the wave operator[END_REF]).

2.3. Hardy and BMO spaces. We now define atomic Hardy spaces adapted to our situation (dictated by a semigroup) using the construction introduced in [START_REF] Bernicot | New abstract Hardy spaces[END_REF]. Again we sum up the definitions and properties we need without proofs. A more detailed explanation with proofs is provided in [START_REF] Bernicot | Dispersive estimates with loss of derivatives via the heat semigroup and the wave operator[END_REF].

Let M be a large enough integer.

Definition 2.5. A function a ∈ L 1 loc is an atom associated with the ball Q of radius r if there exists a function f

Q whose support is included in Q such that a = (1 -e -r 2 H ) M (f Q ), with f Q L 2 (Q) ≤ (µ(Q)) -1 2 .
That last condition allows us to normalize f Q in L 1 . Indeed by the Cauchy-Schwarz inequality

f Q L 1 ≤ f Q L 2 (Q) µ(Q) 1 2 ≤ 1.
Moreover, (1 -e -r 2 H ) M is bounded on L 1 so every atom is in L 1 and they are also normalized in

L 1 : (2.2) sup a a L 1 1,
where we take the supremum over all the atoms. We may now define the Hardy space by atomic decomposition Definition 2.6. A measurable function h belongs to the atomic Hardy space H 1 ato , which will be denoted H 1 , if there exists a decomposition

h = i∈N λ i a i µ -a.e.
where a i are atoms and λ i real numbers satisfying

i∈N |λ i | < +∞.
We equip the space H 1 with the norm

h H 1 := inf h= i λ i a i i∈N |λ i |,
where we take the infimum over all the atomic decompositions.

For a more general definition and some properties about atomic spaces we refer to [Ber10, BZ08], and the references therein. From (2.2), we deduce Corollary 2.7. The Hardy space is continuously embedded into L 1 :

f L 1 f H 1 .
From [BZ08, Corollary 7.2], the Hardy space H 1 is also a Banach space.

We refer the reader to [BZ08, Section 8], for details about the problem of identifying the dual space (H 1 ) * with a BMO space. For a L ∞ -function, we may define the BMO norm

f BMO := sup Q - Q |(1 -e -r 2 H ) M (f )| 2 dµ 1/2
, where the supremum is taken over all the balls Q of radius r > 0. If f ∈ L ∞ then (1-e -r 2 H ) M (f ) is also uniformly bounded (with respect to the ball Q), since the heat semigroup is uniformly bounded in L ∞ (see Proposition 2.2) and so f BMO is finite.

Definition 2.8. The functional space BMO is defined as the closure

BMO := {f ∈ L ∞ + L 2 , f BMO < ∞} for the BMO norm.
Following [BZ08, Section 8], it comes that BMO is continuously embedded into the dual space (H 1 ) * and contains L ∞ :

L ∞ ֒→ BMO ֒→ (H 1 ) * . Hence (2.3) T H 1 →(H 1 ) * T H 1 →BMO ,
and we have the following interpolation result:

(2.4) ∀θ ∈ (0, 1), (L 2 , BMO) θ ֒→ (L 2 , (H 1 ) * ) θ .
The following interpolation theorem between Hardy spaces and Lebesgue spaces is essential in our study.

Theorem 2.9. For all θ ∈ (0, 1), consider the exponent p ∈ (1, 2) and q = p ′ ∈ (2, ∞) given by

1 p = 1 -θ 2 + θ and 1 q = 1 -θ 2 .
Then (using the interpolation notations), we have

(L 2 , H 1 ) θ = L p and (L 2 , (H 1 ) * ) θ ֒→ L q ,
if the ambient space X is non-bounded and

L p ֒→ L 2 + (L 2 , H 1 ) θ and L 2 ∩ (L 2 , (H 1 ) * ) θ ֒→ L q ,
if the space X is bounded. The same results hold replacing (H 1 ) * by BMO thanks to (2.4).

Remark 2.10. We will not mention the case of a bounded space X in the proofs since interpolation is more delicate in that case. One can find the corresponding interpolation theorem (Theorem 2.17 in [START_REF] Bernicot | Dispersive estimates with loss of derivatives via the heat semigroup and the wave operator[END_REF]) and check that the results apply in that case.

2.4. Motivation of the hypothesis. This section is dedicated to the motivation of Assumption 1.1. As we said in the Introduction, this hypothesis is weaker than the one in [START_REF] Bernicot | Dispersive estimates with loss of derivatives via the heat semigroup and the wave operator[END_REF], namely

(2.5) cos(s √ H)ψ m (r 2 H) L 2 (B)→L 2 ( B) r r + s d-1 2 r r + |s -L| d+1 2
, where L = d(B, B). Therefore in all the situations where (2.5) is satisfied, we can assure that Assumption 1.1 is valid. When we have a good knowledge of the wave propagator, we can also affirm that Assumption 1.1 holds. This is the case thanks to a parametrix in [START_REF] Bérard | On the wave equation on a compact Riemannian manifold without conjugate points[END_REF] in the following cases:

• The Euclidean spaces R d with the usual Laplacian H = -∆ = -d j=1 ∂ 2 j ; • A compact Riemannian manifold of dimension d with the Laplace-Beltrami operator;

• A smooth non-compact Riemannian d-manifold with C ∞ b -geometry and Laplace-Beltrami operator;

• The Euclidean space R d equipped with the measure dµ = ρdx and H = -1 ρ div(A∇), where ρ is an uniformly non-degenerate function and A a matrix with bounded derivatives; Moreover we can check that for a non-trapping asymptotically conic manifold with H = -∆ + V the assumption holds. Therefore we recover, with another proof, the result of [START_REF] Hassell | Global-in-time Strichartz estimates on nontrapping, asymptotically conic manifolds[END_REF].

However in [START_REF] Ivanovici | Dispersion for the wave equation inside strictly convex domains I: the Friedlander model case[END_REF] the authors proved that for the Laplacian inside a convex domain of dimension d ≥ 2 in R d , there was a loss of s 1 4 in the dispersion, namely

(2.6) cos(s √ H)ψ m (r 2 H) L 2 (B)→L 2 ( B) r r + s d-1 2 + 1 4 r r + |s -L| d+1 2
.

This loss indicates a difficulty when dealing with boundaries of a domain. The authors used oscillatory integrals techniques and a careful study of the reflections on the boundary of the domain.

A remark of J.-M. Bouclet to get around the use of a parametrix leads us to investigate Klainerman's commuting vectorfields method. It can be found in detail in [START_REF] Klainerman | A commuting vectorfields approach to Strichartz-type inequalities and applications to quasi-linear wave equations[END_REF] and [START_REF] Christopher | Lectures on nonlinear wave equations[END_REF]. Briefly, if one can find enough vectorfields commuting with the wave operator, using a version of Sobolev inequalities, also know as Klainerman-Sobolev inequalities, one can obtain dispersion estimations for the wave propagator. In our setting, we would obtain (see [Sog95, Remark 1.4]) the following dispersion property:

(2.7) cos(s √ H)ψ m (r 2 H) L 2 (B)→L 2 ( B) r r + s d-1 2 r r + |s -L| 1 2 .
It is very close to our Assumption 1.1, but it takes into account the dispersion inside the light cone. In that sense, it is intermediate between our Assumption 1.1 and estimate (2.5). A question we would like to pursue investigating is to find enough well-suited vectorfields to apply this method in generals settings. The framework inwhich one is interested in verifying (2.7) is when H is a given by divergence form, namely H = -div(A∇). When A = I d the identity matrix of size d, H is the usual Laplacian. In this case and the one where A has C 1,1 coefficients, Klainerman obtained in [START_REF] Klainerman | A commuting vectorfields approach to Strichartz-type inequalities and applications to quasi-linear wave equations[END_REF] a dispersion property of the form (2.7). It is not new since it was already proven in [START_REF] Hart | A parametrix construction for wave equations with C 1,1 coefficients[END_REF]. But the novelty in [START_REF] Klainerman | A commuting vectorfields approach to Strichartz-type inequalities and applications to quasi-linear wave equations[END_REF] is to get around the use of a parametrix.

Proofs of the Theorems

This section is dedicated to the proofs of the announced result. It is divided into two main theorems. The first one shows which L 2 -L 2 dispersion property we can recover thanks to Assumption 1.1. In the second theorem we obtain Strichartz inequalities using such dispersive estimates. We recall that our goal is to investigate which properties (in terms of Strichartz inequalities) for the Schrödinger operator can be deduced from a weak assumption on the wave operator.

e itH ψ m ′ (h 2 H)ψ m (r 2 H) L 2 (B)→L 2 ( B) r d |t| d-2 2 h 2
, for all m ′ ≥ 0 and where

0 < |t| ≤ 1 if κ = +∞ and h 2 ≤ |t| ≤ h 1+ε if κ < +∞.
This shows how dispersion for the wave propagator implies dispersion for the Schrödinger group. The main tool to link those two operators is Hadamard's transmutation formula

(3.2) ∀z ∈ C, Re(z) > 0, e -zH = +∞ 0 cos(s √ H)e -s 2 4z ds √ πz .
Proof. Let B, B be balls with radius r > 0. We start our proof with some easy reductions.

First remark that we can restrict ourselves to prove the theorem for h ≤ r. Indeed if the theorem is true for h ≤ r then for all h > r

e itH ψ m ′ (h 2 H)ψ m (r 2 H) L 2 (B)→L 2 ( B) = r 2m ( h 2 2 + r 2 ) m e itH ψ m ′ ( h 2 2 H)ψ m ′ (( h 2 2 + r 2 )H) L 2 (B)→L 2 ( B) r h 2m e itH ψ m ′ ( h 2 2 H)ψ m ′ (( h 2 2 + r 2 )H) L 2 (Bρ)→L 2 ( Bρ)
where

ρ = h 2 2 + r 2 > r, B ρ = ρ r B and B ρ = ρ r B are of radius ρ. Since h 2 2 + r 2 ≥ h 2 2 we then obtain e itH ψ m ′ (h 2 H)ψ m (r 2 H) L 2 (B)→L 2 ( B) r h 2m ρ d |t| d-2 2 h 2 .
We conclude using ρ h and

r 2m h d h 2m = r d r h 2m-d ≤ r d .
Moreover we only need to prove the theorem for m ′ = 0. Indeed if we show

e itH ψ 0 (h 2 H)ψ m (r 2 H) L 2 (B)→L 2 ( B) r d |t| d-2
2 h 2 then for all m ′ ≥ 0 we have

e itH ψ m ′ (h 2 H)ψ m (r 2 H) L 2 (B)→L 2 ( B) = h r 2m ′ e itH ψ 0 (h 2 H)ψ m+m ′ (r 2 H) L 2 (B)→L 2 ( B) r d |t| d-2 2 h 2 since h ≤ r
Finally it is sufficient to consider r 2 ≤ t because if r 2 > t then by bounded functional calculus we have

e itH ψ m ′ (h 2 H)ψ m (r 2 H) L 2 (B)→L 2 ( B) 1 ≤ r 2 |t| d-2 2 ≤ r d |t| d-2 2 h 2 .
In summary, we fix h ≤ r, m ′ = 0, and r 2 ≤ t.

In order to avoid nonzero bracket terms in the forthcoming integrations by parts, we introduce

a technical function χ ∈ C ∞ (R + ) such that      0 ≤ χ ≤ 1 χ(x) = 1 if x ∈ [0, |t| r ] χ(x) = 0 if x ∈ [2 |t| r , +∞]
.

Moreover we have ∀n ∈ N, ∀x ∈ R + , |χ (n) (x)| r |t| n . Thus we split (3.2) into (3.3) e -zH = +∞ 0 χ(s) cos(s √ H)e -s 2 4z ds √ πz + +∞ 0 (1 -χ(s)) cos(s √ H)e -s 2 4z ds √ πz .
We treat the first term by integrations by parts. Making 2n integration by parts (with n to be determined later) we get

∞ 0 cos(s √ H)ψ m (r 2 H)χ(s)e -s 2 4z ds = ∞ 0 cos(s √ H)r 2n ψ m-n (r 2 H) 2n k=0 χ (2n-k) (s)e -s 2 4z c k s k z k + . . . + c n-2⌊ n 2 ⌋ s k-2⌊ k 2 ⌋ z k-⌊ k 2 ⌋ ds,
with (c i ) i being numerical constants playing no significant role. Keeping the extremal terms (one when k = 0 and two when k = 2n) we have to estimate

2 |t| r 0 cos(s √ H)ψ m-n (r 2 H) L 2 (B)→L 2 ( B) r 2n r |t| 2n + s 2n |t| 2n + 1 |t| n ds |t| .
By continuity of our operators

cos(s √ H)ψ m-n (r 2 H) L 2 (B)→L 2 ( B) 1,
we can estimate

2 |t| r 0 r 2 |t| 2n ds |t| r 2 |t| 2n-1 2 and 2 |t| r 0 r 2n |t| n ds |t| r 2 |t| n-1 2 .
Using (1.1) we have

2 |t| r 0 r r + s d-1 2 rs |t| 2n ds |t| ≤ 2 |t| r 0 r d-1 2 +2n |t| 2n+ 1 2 s 2n-d-1 2 ds ≃ r 2 |t| d-2 2 .
Thus, the intermediate terms having the same behaviour, for large enough n

∞ 0 cos(s √ H)ψ m (r 2 H)χ(s)e -s 2 4z ds √ z L 2 (B)→L 2 ( B) r 2 |t| d-2 2 .
Moreover, since h ≤ r we have

r 2 |t| d-2 2 ≤ r d |t| d-2 2 h 2 .
To estimate the second term in (3.3), we treat separately the cases s < κ and s > κ.

+∞ 0

(1 -χ(s)) cos(s √ H)e -s 2 4z ds √ πz = κ |t| r cos(s √ H)e -s 2 4z ds √ πz + I κ ,
where

I κ = 0 if κ = +∞ +∞ κ cos(s √ H)e -s 2 4z ds √ πz if κ < +∞ .
We use the exponential decay for s > κ. Noting z = h 2 -i|t|, and using the L 2 -boundedness of the cos(s √ H)ψ m (r 2 H) operator:

+∞ κ cos(s √ H)ψ m (r 2 H)e -s 2 4z ds √ z L 2 (B)→L 2 ( B) +∞ κ e -s 2 8 Re 1 z e -s 2 8 Re 1 z ds |z| +∞ κ 8 Re 1 z e -u 2 du |z| Re 1 z e -κ 2 Re 1 z 8 +∞ 0 e -u 2 du |t| h e -κ 2 h 2 16t 2 h 2 t 2 -N
|t| h for all N ≥ 1 as large as we want and where we used |z| ≃ |t| and Re

1 z ≥ h 2 2t 2 . Moreover |t| 2N + 1 2 h 2N +1 ≤ h d |t| d-2 2 h 2 ≤ r d |t| d-2 2 h 2 as soon as |t| 2N + d-2 2 + 1 2 ≤ h 2N +d-1 that is |t| ≤ h 1+ d-1 2 2N+ d-1 2 . Which is true since |t| ≤ h 1+ε ≤ h 1+ d-1 2 2N+ d-1 2
for N large enough.

Remark 3.2. We point out that this is the only moment we use that |t| ≤ h 1+ε . That is why we do not need it when κ = +∞ since this term does not step in. Therefore the loss of derivatives in Theorem 1.3 is better when κ = +∞.

We use Assumption 1.1 when s < κ. Indeed it yields

κ |t| r cos(s √ H)ψ m (r 2 H)e -s 2 4z ds √ z L 2 (B)→L 2 ( B) κ |t| r r r + s d-1 2 e -s 2 4 Re 1 z ds |t| .
When d-1 2 > 1 (i.e. d > 3) we have

κ |t| r r r + s d-1 2 ds |t| ≤ r d-1 2 |t| ∞ |t| r s -d-1 2 ds r d-1 2 |t| |t| r -d-1 2 +1 ≤ r d-2 |t| d-2 2 ≤ r d |t| d-2 2 h 2 h 2 ≤ r d |t| d-2 2 h 2 since h 2 ≤ 1. When d-1 2 < 1 (i.e. d < 3), since Re 1 z h 2 t 2 we have κ |t| r r r + s d-1 2 e -c s 2 h 2 t 2 ds |t| r d-1 2 |t| ∞ h r e -u 2 |t|u h -d-1 2 |t| h du r d-1 2 h d-3 2 |t| d-2 2 ∞ 0 u -d-1 2 e -u 2 du r d |t| d-2 2 h 2 h d-3 2 h 2 r d+1 2 ≤ r d |t| d-2
since h ≤ r. When d-1 2 = 1 (i.e. d = 3) we have

κ |t| r r r + s e -s 2 4 Re 1 z ds |t| r |t| κ h |t| h r |t|u h -1 e -u 2 |t| h du ≤ r |t| e -h 2 2r 2 r h +∞ 0 e -u 2 2 du r 2 |t|h h r -1 = r 3 |t|h 2 = r d |t| d-2 2 h 2 .
In the end, summing all the parts up, we have

e itH ψ m ′ (h 2 H)ψ m (r 2 H) L 2 (B)→L 2 ( B) r d |t| d-
(H m (A)) T ψ m (r 2 H) L 2 (Br)→L 2 ( Br) Aµ(B r ) 1 2 µ( B r ) 1 2
for some m ≥ d 2 . Then T is bounded from H 1 to BM O and from L p to L p ′ for p ∈ (1, 2) with

T H 1 →BM O A and T L p →L p ′ A 1 p -1 p ′
if the ambient space X is unbounded and

T H 1 →BM O max(A, 1) and T L p →L p ′ max(A 1 p -1 p ′ , B)
if the ambient space X is bounded, and where, for the last inequality, we assumed that T L p →L 2 B.

As we mentioned previously, we do not use the part where X is bounded. We apply the Theorem with T = e itH ψ m ′ (h 2 H) and A = |t| -d-2 2 h -2 . In view of (1.4) we can reformulate (H m (A)) (see [START_REF] Bernicot | Dispersive estimates with loss of derivatives via the heat semigroup and the wave operator[END_REF]) as

(3.4) e itH ψ m ′ (h 2 H)ψ m (r 2 H) L 2 (B)→L 2 ( B) r d |t| d-2
2 h 2 which we just proved in the previous section under our assumption. Therefore we obtain

e itH ψ m ′ (h 2 H) H 1 →BM O |t| -d-2 2 h -2 ,
and for all p ∈ (1, 2)

e itH ψ m ′ (h 2 H) L p →L p ′ h -2 |t| -d-2 2 1 p -1 p ′ .
We now recall a slightly modified version of a result of Keel-Tao in [START_REF] Keel | Endpoint Strichartz estimates[END_REF]:

Theorem 3.4. If (U (t)) t∈R satisfies sup t∈R U (t) L 2 →L 2 1
and for some σ > 0 ∀t = s, U (t)U (s

) * H 1 →BMO ≤ C|t -s| -σ .
Then for all 2 ≤ p ≤ +∞ and 2 ≤ q < +∞ satisfying

1 p + σ q = σ 2 we have U (t)f L p t L q x C 1 2 -1 q f L 2 .
Proof. We just sum up the main steps of the proof in [START_REF] Keel | Endpoint Strichartz estimates[END_REF] to keep track of the constant in the last estimation.

• By symmetry and a T * T argument, it suffices to show s<t

U (s) * F (s), U (t) * G(t) dsdt C 2 F L p ′ t L q ′ x G L p ′ t L q ′ x .
• By the interpolation Theorem 2.9 we have

U (t)U (s) * L q ′ →L q C 1-2 q |t -s| -2 p .
• We conclude by Hölder and Hardy-Littlewood-Sobolev inequalities.

We use this theorem with C = 1 h 2 and σ = d-2 2 to obtain the following result. Theorem 3.5.

Under Assumption 1.1, if 2 ≤ p ≤ +∞ and 2 ≤ q < +∞ satisfy 2 p + d -2 q = d -2 2 ,
and f ∈ L 2 and 0 < h ≤ 1 we have

• if κ = +∞ then for all m ′ ∈ N e itH ψ m ′ (h 2 H)f L p ([-1,1],L q ) 1 h 2( 1 2 -1 q ) ψ m ′ (h 2 H)f L 2 ; • if κ < +∞ then for all 0 < ε < 1 and m ′ ∈ N e itH ψ m ′ (h 2 H)f L p ([-1,1],L q ) 1 h 1+ε p h 2( 1 2 -1 q ) ψ m ′ (h 2 H)f L 2 .
Proof. The following proof is a slight modification of the one of Theorem 4.2 and 4.3 of [START_REF] Bernicot | Dispersive estimates with loss of derivatives via the heat semigroup and the wave operator[END_REF]. We rewrite it here for more readability. We only deal with the case κ < +∞ since it is more technical. We leave the minor modifications to obtain the case κ = +∞ to the readers. Fix an interval J ⊂ [-1, 1] of length |J| = h 1+ε , m ′ ∈ N, and consider

U (t) = 1 J (t)e itH ψ m ′ (h 2 H).
We aim to apply Theorem 3.4 with C = 1 h 2 and σ = d-2 2 . By functional calculus we have sup

t∈R U (t) L 2 →L 2 1.
The estimation (3.4) which we proved in Theorem 3.1 will lead to the second hypothesis of Theorem 3.4. First

U (t)U (s) * = 1 J (t)1 J (s)e itH ψ m ′ (h 2 H)(e isH ψ m ′ (h 2 H)) * = 1 J (t)1 J (s)e i(t-s)H ψ 2m ′ (2h 2 H)
because H is self-adjoint. Since J has length equal to h 1+ε then U (t)U (s) * is vanishing or else |t -s| ≤ h 1+ε . Hence, by Theorem 3.1 we deduce

U (t)U (s) * H 1 →(H 1 ) * 1 h 2 1 |t -s| d-2 2 .
Up to the change of 2m ′ into m ′ , Theorem 3.4 (with C = h -2 and σ = (d -2)/2) then leads to

J e itH ψ m ′ (h 2 H)f p L q dt 1 p 1 h 2( 1 2 -1 q ) f L 2 . We then split [-1, 1] into N ≃ 1 h 1+ε intervals J k of length h 1+ε to obtain 1 -1 e itH ψ m ′ (h 2 H)f p L q dt ≤ N k=1 J k e itH ψ m ′ (h 2 H)f p L q dt ≤ N 1 h 2( 1 2 -1 q ) f L 2 p .
Hence

e itH ψ m ′ (h 2 H)f L p ([-1,1],L q ) 1 h 1+ε p h 2( 1 2 -1 q ) ψ m ′ (h 2 H)f L 2 .
We are now able to prove Strichartz estimates with loss of derivatives.

Theorem 3.6. If Assumption 1.1 is satisfied. Then for every 2 ≤ p ≤ +∞ and 2 ≤ q < +∞ satisfying

2 p + d -2 q = d -2 2 ,
and every solution u(t, . ) = e itH u 0 of the problem

i∂ t u + Hu = 0 u |t=0 = u 0 ,
we have • if κ = +∞, then u satisfies local-in-time Strichartz estimates with loss of derivatives

(3.5) u L p ([-1,1],L q ) u 0 W 2( 1 2 -1
q ),2 ; • if κ < +∞, then for all 0 < ε < 1 and 0 < h ≤ 1, u satisfies local-in-time Strichartz estimates with loss of derivatives

(3.6) u L p ([-1,1],L q ) u 0 W 1+ε p +2( 1 2 -1
q ),2 . Remark 3.7. The loss of derivatives in (3.6) is interesting when it is lower than the straightforward loss given by Sobolev embeddings. The relation 2 p + d-2 q = d-2 2 yields W 2 p +1-2 q ,2 ֒→ L q .

Thus e itH u 0 L q e itH u 0 W

2 p +1-2 q ,2 ≤ u 0 W 2 p +1-2 q ,2
and taking the L p ([-1, 1]) norm shows e itH u 0 L p ([-1,1],L q ) u 0 W 2 p +1-2 q ,2 . That is, the loss of derivatives is interesting when it is less than 2 p + 1 -2 q . Hence, for all ε ∈ (0, 1), the loss we obtained is strictly better than the one directly given by Sobolev embeddings. The loss in (3.5) is also nontrivial by the same argument. In order to do so we remark that in (3.4) we could write

r d t d-2 2 h 2 = r 2 t d 2 t h 2 ≤ r 2 t d 2 1 h ,
because t ≤ h. Then the loss of derivatives obtained in (3.6) is 1+ε p + 1( 1 2 -1 q ) that need to be compared to the trivial loss 2 p . Since 1 2 -1 q = 2 dp , the loss is less than 2 p if an only if

d ≥ 2 1 -ε .
That is, as soon as d > 2, one can find ε ∈ (0, 1) such that the loss is nontrivial. We chose to present the previous Theorem in that form because it allows a wider range of exponent q. Indeed, on the one hand

p ≥ 2 ⇒ d -2 q = d -2 2 - 2 p ≥ d -2 2 -1 that is 1 q ≥ 1 2 - 1 d -2 .
On the other hand p ≥ 2 and 2 p + d q = d 2 yields

1 q ≥ 1 2 - 1 d ,
and for all d > 2, 1 2

- 1 d -2 ≤ 1 2 - 1 d .
That is why the relation

2 p + d -2 q = d -2 2
gives a wider range for exponent q.

Proof of Theorem 3.6. Again we only deal with the more difficult case κ < +∞. We apply Theorem 2.4 to u(t) = e itH u 0 . It leads to

u(t) L q ϕ(H)u(t) L q + 1 0 |ψ m ′ (s 2 H)u(t)| 2 ds s 1 2 L q
, with m ′ ≥ 1. Taking the L p ([-1, 1]) norm in time of that expression and using Minkowski inequality yields u(t) L p ([-1,1],L q ) ϕ(H)u(t) L p ([-1,1],L q ) + 1 0 ψ m ′ (s 2 H)u(t) 2 L q ds s 1 2

L p ([-1,1])

.

Thanks to the Gaussian pointwise estimate of ϕ(H) the first term can be estimated as follow ϕ(H)u(t) L p ([-1,1],L q ) e itH u 0 L p ([-1,1],L 2 ) u 0 L 2 u 0 W 1+ε p +2( 1 2 -1 q ),2 .

Since p ≥ 2, Theorem 3.5 and generalized Minkowski inequality allow to bound the second term

1 0 ψ m ′ (s 2 H)u(t) 2 L q ds s 1 2 L p ([-1,1]) 1 0 ψ m ′ (s 2 H)u 2 L p ([-1,1],L q ) ds s 1 2 1 0 s -1+ε p -2( 1 2 -1 q ) ψ m ′ (s 2 H)u 0 2 L 2 ds s 1 2 1 0 s -1+ε p -2( 1 2 -1 q ) |ψ m ′ (s 2 H)u 0 | 2 ds s 1 2 L 2 u 0 W 1+ε p +2( 1 2 -1 q ),2
, where we used m ′ ≥ 1 4 [ 1+ε p + 2(1 -2 q )] since m ′ ≥ 1 and 1+ε p + 2(1 -2 q ) < and the fact that s

-1+ε p -2( 1 2 -1 q ) |ψ m ′ (s 2 H)| 2 = ψ m ′ -1 4 [ 1+ε p +2( 1 2 -1 q )] (s 2 H)H 1 2 [ 1+ε p +2( 1 2 -1 q )] .
Finally, we get u L p ([-1,1],L q ) u 0 W 1+ε p +2( 1 2 -1 q ),2 .

r d- 1

 1 balls of radius r and use Theorem 2.1, we can conjecture that the term r r+s d-1 2

3. 1 .

 1 Dispersive estimates for the Schrödinger operator. The main theorem of this section is the following Theorem 3.1. Assume d > 2, m ≥ ⌈ d 2 ⌉, and that Assumption 1.1 is satisfied, then for all balls B, B of radius r > 0 and all 0 < h ≤ 1 (3.1)

  Remark 3.8. One could work out our estimate with 2 p + d q = d 2 .

  Strichartz inequalities. To obtain Strichartz estimates we are going to use Theorem 1.1 of[START_REF] Bernicot | Dispersive estimates with loss of derivatives via the heat semigroup and the wave operator[END_REF], which we recall here with a slight modification in assuming (1.4), namely Theorem 3.3. Assume (1.4) with (DUE). Consider a self-adjoint and L 2 -bounded operator T (with T L 2 →L 2 1), which commutes with H and satisfies

	2 2 h 2	.
	3.2.

V. Samoyeau's research is supported by Centre Henri Lebesgue (program "Investissements d'avenir" -ANR-11-LABX-0020-01).

This work benefits from and fits into the ERC project FAnFArE no. 637510.


Guillaume Ducrozet
email: guillaume.ducrozet@ec-nantes.fr

Félicien Bonnefoy

David Le Touzé

Pierre Ferrant

HOS-ocean: Open-source solver for nonlinear waves in open ocean based on High-Order Spectral method

Keywords: High-Order Spectral method, Nonlinear waves, Wave propagation, Wave kinematics, Ocean engineering High-Order Spectral method, Nonlinear waves, Wave propagation, Wave kinematics, Ocean engineering Classification: 4.12 (and eventually: 7.7 & 12) External routines/libraries: FFTW 3.3.4 [17] and LAPACK [3]

niveau recherche, publiés ou non, émanant des établissements d'enseignement et de recherche français ou étrangers, des laboratoires publics ou privés.

Introduction

Accurate modeling of nonlinear wave propagation and induced kinematics is of particular concern in marine engineering, oceanography, security of marine transportation. . . Nowadays, the main challenge is to simulate large domains over long time periods with nonlinear features taken care of. This led in last decades to the development of several numerical methods able to solve the nonlinear evolution of sea-states. In this context, High-Order Spectral (HOS) model exhibits high efficiency and accuracy thanks to its pseudo-spectral formalism. When dealing with the propagation of openocean wavefields over constant water depth, HOS models have been shown to be more efficient than advanced models solving the volume problem with finite-differences discretization [START_REF] Ducrozet | A comparative study of two fast nonlinear free-surface water wave models[END_REF].

HOS models have been widely used and validated by several authors to study different physical mechanisms: e.g. in recent years nonlinear energy transfers [START_REF] Tanaka | Verification of Hasselmann's energy transfer among surface gravity waves by direct numerical simulations of primitive equations[END_REF], modulational instabilities [START_REF] Fernandez | Modulational instability and wave amplification in finite water depth[END_REF], [START_REF] Toffoli | Evolution of weakly nonlinear random directional waves: laboratory experiments and numerical simulations[END_REF], bi-modal seas [START_REF] Toffoli | Development of a bimodal structure in ocean wave spectra[END_REF], freak waves [START_REF] Ducrozet | 3-D HOS simulations of extreme waves in open seas[END_REF][START_REF] Xiao | Rogue wave occurence and dynamics by direct simulations of nonlinear wave-field evolution[END_REF][START_REF] Sergeeva | Rogue waves, rogue events and extreme wave kinematics in spatio-temporal fields of simulated sea states[END_REF] . . . among others. This method can consequently be considered mature and accessible to practical engineering applications.

HOS-ocean is an open-source HOS model developed at Ecole Centrale

Nantes, LHEEA Lab. (ECN/CNRS). The software is available to download and contribute on the GitHub platform [2]. The code is developed and redistributed under the terms of the GNU General Public License as published by the Free Software Foundation. Along with the source code, documentation that describes the compilation and execution of the source files is also distributed. Note that HOS-ocean has been extensively validated on nonlinear regular wave propagation [START_REF] Bonnefoy | Advances in Numerical Simulation of Nonlinear Water Waves[END_REF] and applied to the study of freak waves modeling in [START_REF] Ducrozet | 3-D HOS simulations of extreme waves in open seas[END_REF].

One of the purposes of this code is to encourage other researchers to use HOS. Along with the core of models, post-processing is also delivered. This allows in particular the possibility to easily couple HOS-ocean with CFD (Computational Fluid Dynamics) software for wave-structure interactions modeling.

In the following sections, the HOS formulation implemented in HOSocean and associated numerical procedures are described. Sections describing how to compile and run the code are also provided and finally, several study cases are presented as typical applications of the presented numerical model.

High-Order Spectral method

In this section, the High-Order Spectral method used in HOS-ocean code is presented. HOS models have been widely used for the study of wave propagation in open domains starting from the original work of [START_REF] West | A new numerical method for surface hydrodynamics[END_REF] & [START_REF] Dommermuth | A high-order spectral method for the study of nonlinear gravity waves[END_REF] up to the analysis of complex wave phenomena [START_REF] Tanaka | Verification of Hasselmann's energy transfer among surface gravity waves by direct numerical simulations of primitive equations[END_REF][START_REF] Toffoli | Evolution of weakly nonlinear random directional waves: laboratory experiments and numerical simulations[END_REF][START_REF] Toffoli | Development of a bimodal structure in ocean wave spectra[END_REF][START_REF] Fernandez | Modulational instability and wave amplification in finite water depth[END_REF], including freak waves [START_REF] Ducrozet | 3-D HOS simulations of extreme waves in open seas[END_REF][START_REF] Xiao | Rogue wave occurence and dynamics by direct simulations of nonlinear wave-field evolution[END_REF][START_REF] Sergeeva | Rogue waves, rogue events and extreme wave kinematics in spatio-temporal fields of simulated sea states[END_REF], among others. Extensions of the original HOS model to include wind forcing effects [START_REF] Touboul | On the influence of wind on extreme wave events[END_REF][START_REF] Kharif | Influence of wind on extreme wave events: experimental and numerical approaches[END_REF] and jet current effects [START_REF] Shrira | Nonlinear dynamics of trapped waves on jet currents and rogue waves[END_REF] are also worth to mention. In this section, unless specified, the method is presented in one horizontal direction (i.e. 2DV problem) for conciseness, extension to 3D being straightforward (see [START_REF] Bonnefoy | Advances in Numerical Simulation of Nonlinear Water Waves[END_REF][START_REF] Ducrozet | 3-D HOS simulations of extreme waves in open seas[END_REF]).

Formulation

We consider a rectangular fluid domain D of horizontal dimension L x and constant water depth h associated with a Cartesian coordinate system. Its origin O is located at one corner of the domain with Ox representing the horizontal axis and Oz the vertical one oriented upward with z = 0 located at the mean free surface. We are working under the potential flow theory (assuming the fluid to be incompressible and inviscid and the flow irrotational). Then, continuity equation reduces to the Laplace equation for the velocity potential φ

∇ 2 φ + ∂ 2 φ ∂z 2 = 0 in D (1)
with ∇ denoting the horizontal gradient operator. The boundary conditions will close the system of equations to solve. At first, we consider the free surface boundary conditions which are written, following [START_REF] Zakharov | Stability of periodic waves of finite amplitude on the surface of a deep fluid[END_REF], using surface quantities namely the free surface elevation η and the free surface velocity potential φ(x, t) = φ(x, z = η(x, t), t). z = η(x, t) describes the free surface position, assuming no wave breaking occurs (i.e. the free surface is a single valued function of x). The kinematic and dynamic free surface boundary conditions then read

∂η ∂t = 1 + |∇η| 2 W -∇ φ • ∇η (2) ∂ φ ∂t = -gη - 1 2 ∇ φ 2 + 1 2 1 + |∇η| 2 W 2 (3) W (x, t) = ∂φ ∂z (x, z = η, t)
is the vertical velocity at the free surface, the only quantity beyond the system (2)-(3) which needs solution in the water bulk.

Latter W will be evaluated thanks to the order consistent HOS scheme of West et al. [START_REF] West | A new numerical method for surface hydrodynamics[END_REF] described hereafter (Sec. 1.2). Once the vertical velocity is evaluated, it is possible to advance in time the two unknowns η and φ (Sec.

1.4).

In this model, periodic lateral boundary conditions are used, assuming the domain to be infinite

η , φ , W (x = 0, t) = η , φ , W (x = L x , t) (4)
The latter associated to Laplace equation (1) and bottom boundary con-

dition ∂φ ∂z (x, z = -h, t) = 0 (5)
allow us to define a spectral basis Ψ m on which the potential is expanded.

φ(x, z, t) = m A m (t)Ψ m (x, z) = m A m (t) cosh [k m (z + h)] cosh (k m h) exp (ik m x) (6)
with k m = m∆k x = m 2π Lx the wavenumbers. Similarly, surface quantities are also expressed on a spectral basis allowing the use of Fast Fourier Transforms (FFTs).

φ(x, t) = m B φ m (t) exp (ik m x) (7) η(x, t) = m B η m (t) exp (ik m x) (8)
Expansions ( 7)-( 8) are truncated to a given number of modes which defines the corresponding number of points N in physical space. Some computing operations will be done in physical space (typically products of quantities)

and others in transformed (Fourier) space (evaluation of spatial derivatives).

Then, the resulting numerical method is pseudo-spectral and exhibits very interesting convergence properties. The corresponding HOS model features high efficiency and accuracy compared to other advanced methods for wave propagation, see [START_REF] Ducrozet | A comparative study of two fast nonlinear free-surface water wave models[END_REF].

HOS scheme

Knowing the surface quantity η and φ, the HOS scheme intends to evaluate the vertical velocity at the free surface W (x, t). This procedure relies on a series expansion in wave steepness ǫ up to the so-called HOS order

M φ(x, z, t) = M m=1 φ (m) (x, z, t) (9)
with φ (m) quantities of order ǫ m . Associated with a Taylor series around z = 0 and collecting terms at each order in wave steepness lead to a triangular system for the different φ (m) .

φ (1) (x, 0, t) = φ(x, t) (10) φ (m) (x, 0, t) = - m-1 k=1 η k k! ∂ k φ (m-k) ∂z k (x, 0, t) for m > 1 (11)
This way, the complicated Dirichlet problem for φ(x, z, t)

on z = η(x, t)
is transformed into M simpler Dirichlet problems for φ (m) (x, z, t) on z = 0.

Similar series expansion done on the velocity potential (Eq. ( 9)) is applied to the vertical velocity W , leading to another triangular system (for W (m) ) solved iteratively.

W (m) (x, t) = m-1 k=0 η k k! ∂ k+1 φ (m-k) ∂z k+1 (x, 0, t) (12)
Previous set of equations consequently allows to evaluate the free surface vertical velocity of interest

W (x, t) = M m=1 W (m) (x, z, t) (13)

Dealiasing

HOS-ocean code is based on a pseudo-spectral method meaning that certain operators are treated in the physical space instead of spectral space. This is typically the case of all nonlinear products appearing in the freesurface boundary conditions Eqs.

(2)-(3) or in the HOS procedure. These transforms between physical and spectral space are efficient thanks to the use of FFTs but leads to the well-known aliasing phenomena. In order to take care of this problem, we use zero-padding on spectral description of quantities. The number of collocation points has to be chosen to remove aliasing errors on multiple products. Those are at most products of order M. The extended description of quantities uses the half-rule and is defined on a number of point

N d = p + 1 2 N ( 14
)
Full dealiasing is recommended and is achieved by choosing p = M. This ensures that the description with N modes is free from any aliasing errors.

However, if needed (especially for computational effort reduction), partial dealiasing is available in the code by choosing lower value for p. Note that it appears crucial to perform careful dealiasing to preserve method convergence and accuracy as shown in [START_REF] Bonnefoy | Advances in Numerical Simulation of Nonlinear Water Waves[END_REF]. Finally, we indicate that in 3D (i.e. 2 horizontal dimensions), dealiasing has to be performed along both directions

x and y.

Time integration

Free surface boundary conditions Eqs.

(2)-( 3) allow the time-stepping of the two unknowns η and φ. This system of equations appear very stiff for high frequencies, which may lead to small time-steps. Following [START_REF] Fructus | An efficient model for three-dimensional surface wave simulations. Part I: Free space problems[END_REF], thanks to an adequate change of variables, it is possible to solve the nonlinear part of the equations as a problem in itself and not as an adjustment of the linear solution. The linear part of the equations is then analytically integrated while nonlinear evolution of the system is computed numerically. This is achieved thanks to an efficient 4 th order Runge-Kutta Cash-Karp scheme with adaptive step size [START_REF] Cash | A variable order Runge-Kutta method for initial value problems with rapidly varying right-hand sides[END_REF]. A desired accuracy (or so-called tolerance) is specified and the time-step is automatically chosen to achieve this accuracy. Typical values of the tolerance are in the range [10 -5 ; 10 -7 ] to achieve convergence. This choice depends mainly on the nonlinearity of the wave field and the duration of the simulation.

Initial conditions

With such an approach, the key point lays in the definition of an appropriate initial condition. The initial free surface elevation η and free surface velocity potential φ have to be determined before propagating the corresponding wave field in the HOS-ocean model. Different possibilities are offered within the model.

Regular waves

The nonlinear regular wave solution of Rienecker & Fenton [START_REF] Rienecker | A Fourier approximation method for steady water waves[END_REF] (stream function theory) may be specified to initialize the HOS-ocean solution. This is mainly used for validation and to characterize the accuracy of the solution as well as the convergence properties of the method. Especially, a benchmarking procedure is set up in HOS-ocean code, that uses this initialization.

Note that within HOS-ocean, the stream function solution is stored in a specific file for each wave condition.

Irregular waves

For completeness, the initialization of irregular wavefield is presented in 3D configuration, i.e. directional sea. For the definition of irregular seastates, we use a classical superposition of linear components, each component having its own frequency and direction. The amplitude of each component is determined from a directional spectrum S(ω, θ). The relation between the amplitude in the discrete Fourier domain B mn and the spectrum is

1 2 |B η mn | 2 = S(k x , k y )∆k x ∆k y = ∂ω ∂k 1 k S(ω, θ)∆k x ∆k y ( 15
)
with ∆k x and ∆k y the modal discretization in both directions. This equation allows to calculate the norm |B η mn (t = 0)| for the desired directional spectrum S(ω, θ). Then, the phase of each component B η mn (t = 0) is determined by a random number so that the phases assume uniform distribution over the range [0, 2π]. The free surface potential is linearly evaluated from the elevation B φ mn = -i ωmn g B η mn where ω mn is given by the adequate dispersion relation. Finally, the initial wavefield is constructed with an inverse Fourier transform.

HOS-ocean provides the possibility to use a directional spectrum form S(ω, θ) = F (ω).G(θ) defined with: i) a JONSWAP spectrum, specifying the peak period T p = 2π ωp , the significant wave height H s and the shape factor γ

F (ω) = α J H 2 s ω 4 p ω -5 exp   - 5 4 ω ω p -4   γ exp - (ω-ωp) 2 2σ 2 ω 2 p ( 16
)
with σ =     
0.07 for ω < ω p 0.09 for ω ≥ ω p and α J chosen to obtain the correct significant wave height.

ii) a specific directional spreading:

G(θ) = 1 β cos πθ 2β 2 (17)

Relaxation period

It has been demonstrated in [START_REF] Dommermuth | The initialization of nonlinear waves using an adjustment scheme[END_REF], that the definition of linear initial wavefields for fully-nonlinear computations can lead to numerical instabilities. In order to overcome this problem, a transition period is introduced through the use of a relaxation scheme, allowing a smooth transition from linear initial conditions to fully nonlinear computation. The transition period is described by its duration T a and the parameter n. The free surface boundary conditions described by Eqs.

(2)-( 3) are rewritten as

∂ φ ∂t + gη = F , ∂η ∂t -W (1) = G ( 18
)
with W (1) the linearized free surface vertical velocity. Then, the relaxation period is applied on the nonlinear terms F and G as follows At the same time, several studies need the extraction of individual waves inside the computed sea state (we refer for instance to [START_REF] Ducrozet | 3-D HOS simulations of extreme waves in open seas[END_REF] for the study of freak waves). Thus, post-processing of 2D wave fields allows a wave-by-wave analysis using zero-crossing to extract individual waves and possible useful information: H 1/3 , H max , mean wavelength λ 0 . . .

∂ φ ∂t + gη = F 1 -exp - t T a n ( 19
)
∂η ∂t -W (1) = G 1 -exp - t T a n ( 20
In three dimensions, we choose to define the height of a wave as its height in the direction of propagation. We first perform the same kind of analysis than in two dimensions (zero-crossing along x-axis which is the mean direction of propagation of our wavefield). This is performed recursively with respect to the y-direction to get the height at the wave peak location in y.

Thus, we obtain the characteristics of each wave and can deduce those of the whole wavefield: H 1/3 , H max , λ 0 . . . At the same time, we have also access to transverse wavelength λ y , height H y that are interesting in terms of incidence of waves on a structure for example.

Kinematics computation within fluid domain

HOS models use a surface formulation where only quantities on the free surface are advanced in time, namely η and φ. This allows a reduction of numerical complexity and high efficiency. However, it appears necessary in some configurations to access to the volume information, usually the kinematics/pressure inside domain. For the efficiency of this evaluation, preliminary computations are done during the main HOS-ocean execution. Then, postprocessing enables the reconstruction of those quantities.

The main idea is to use a modal description of the volume information (for instance modal coefficients A m for the velocity potential φ(x, z, t) in Eq. ( 6)). Once they are known, it is possible to reconstruct easily this quantity wherever needed in the fluid domain. Different strategies exist to compute from the quantity known at the free surface the corresponding modal coefficients of the volume field. For instance, the HOS scheme solves this problem for φ(x, z, t). But, it has been shown in [START_REF] Bateman | On the calculation of the water particle kinematics arising in a directionnally-spread wavefield[END_REF] that in order to achieve a good accuracy for velocity fields, more advanced procedures need to be used.

HOS-ocean provides the possibility to perform a direct inversion of the full system. However, the approximate solution using the so-called H 2 operator developed in the context of Dirichlet to Neumann Operator (DNO) models is preferred (especially in 3D where the direct inversion is not possible) [START_REF] Bateman | On the calculation of the water particle kinematics arising in a directionnally-spread wavefield[END_REF].

The accuracy of the strategy has been validated in details in [START_REF] Ducrozet | Modélisation des processus non-linéaires de génération et de propagation d'états de mer par une approche spectrale[END_REF]. This way, it is possible to compute easily during post-processing the kinematics induced by the nonlinear waves simulated in HOS-ocean. This enables the possible coupling between HOS-ocean and CFD solvers for the evaluation of wave-structure interactions. Note that the methodology has been used

at ECN for coupling between HOS and RANS solver with the SWENSE method (Spectral Wave Explicit Navier-Stokes Equations). Details about the SWENSE coupling strategy are given in [START_REF] Luquet | Applications of the SWENSE method to seakeeping simulations in irregular waves[END_REF][START_REF] Monroy | Rans simulations of a calm buoy in regular and irregular seas using the SWENSE method[END_REF]. HOS-ocean has been also coupled to a Smooth Particle Hydrodynamics (SPH) model to study violent wave-structure interactions inducing wave breaking, surface reconnexions. . . see [START_REF] Oger | Advances in SPH for naval hydrodynamics[END_REF].

Note that the proposed post-processing allows to compute velocity and pressure maps in a given domain. For a specific coupling of HOS-ocean, the quantities to reconstruct have to be adapted for the problem solved but the core is the one proposed in this post-processing program.

Program documentation

Figure 1 presents a schematic view of the structure of HOS-ocean code.

The details of the different source-files are provided in next sections.

Source files 2.1.1. HOS-ocean

A set of Fortran files need to be compiled in order to create HOS-ocean binary. Here, all the files are listed. Table 1 presents the different source files In order to have a clear view of the different dependencies for HOS-ocean program, we provide a dependency matrix in Tab. 2. Value of 1 indicates that the source file on the line depends on the source file in the column and 0 is allocated to the diagonal of the matrix.

Post-processing & Benchmarking

Tables 3 &4 give the list of the specific files used for the post-processing and the benchmarking respectively.

Compilation

The code can be compiled on any computer architecture. One only needs a Fortran compiler (for instance GNU's gfortran). Makefiles are provided with the source code for compiling HOS-ocean (main program, benchmarking and post-processing) in Linux/Unix environment. These make use of the

HOS-ocean.f90 0 1 1 1 1 1 1 1 1 1 1 1 • • • • • • type.f90 • 0 • • • • • • • • • • • • • • • • variables_3D.f90 • 1 0 • • • • 1 • • • • • • • • • • output.f90 • 1 1 0 • • 1 • 1 • • • • • • • • • runge_kutta.f90 • 1 1 • 0 • 1 • • • • • • • • • 1 1 energy_calc.f90 • 1 1 • • 0 1 • • • • • • • • • • • fourier_r2c_FFTW3.f90 • 1 1 • • • 0 • • • • • • • • • • • ramp.f90 • 1 • • • • • 0 • • • • • • • • • • input_HOS.f90 • 1 1 • • • • • 0 • • • • • • • • • initial_condition.f90 • 1 1 • 1 1 1 • 1 0 • 1 1 1 • • • • filters.f90 • 1 1 • • • • • • • 0 • • • • • • • linear_wave.f90 • 1 • • • • • • • • • 0 • • 1 • • • Bivar.f90 • • • • • • • • • • • • 0 • • • • • RF_solution.f90 • 1 1 • • • • • • • • • • 0 • • • • maths.f90 • 1 • • • • • • • • • • • • 0 • • • nrutil_tmp.f90 • 1 • • • • • • • • • • • • • 0 • • velocities.f90 • 1 1 • • • 1 • • • • • • 1 • • 0 1 resol_HOS.f90 • 1 1 • • • 1 • • • 1 • • • • • • 0

Running HOS-ocean

HOS-ocean has been developed for command-line run with an input file containing all specifications needed. All output files will be created in a directory Results that has to be created before the run by the user.

Important parameters

At first, for the compilation, user has to set values of integers M (the HOS order), n1, n2 (the number of points in x and y directions respectively), p1 and p2 (the parameters used for dealiasing (as seen in Sec. 1.3) in x and y directions respectively) in variables_3D.f90:

Input file

The input file is assumed to be named input_HOS.dat. Table 5 describes the different parameters accessible in this input file.

Output files

Depending on the choices made in the input file (see Tab.

Wavefield analysis

The wavefield analysis is performed on the 3d.dat files generated by HOSocean code. This analysis intends to provide interesting information about the computed sea state. Different output files are created that are described hereafter.

• Analysis.dat contains the general information about the wave field simulated with HOS-ocean • Cartesian grid, i_card = 1 -data_VP_card.dat gives useful informations about the map created in VP_card.dat. This is used at ECN as the file containing coupling information.

-
-VP_card.dat contains the velocities and pressure map. The mesh in z direction is defined as constant z for all x and y. This file is adapted to Tecplot visualization.

• Boundary-fitted mesh, i_card = 2 -VP_card_fitted.dat contains the velocities and pressure map in the case i_card=2, i.e. when a boundary fitted (free surface) mesh is used. This file is adapted to Tecplot visualization.

Results

This section presents different results obtained with HOS-ocean code. At first, validation results obtained with the benchmarking procedure included in the code release are presented. They are dealing with the accuracy of the model when studying a nonlinear regular wave propagation. Then, an example of application of the HOS-ocean model to the propagation of 3D irregular sea-states is provided. In the end, the computation of kinematics inside the fluid domain is demonstrated for a 2D irregular wave field exhibiting an extreme event.

Benchmarking-validation

The nonlinear potential flow solver HOS-ocean presented in this paper uses fully nonlinear free surface boundary conditions Eqs.

(2)-(3) to advance in time the two unknowns η and φ. With the chosen formulation, the accuracy is quantified by the one achieved on the vertical velocity at the free surface, evaluated by the so-called HOS method.

Following the test made by Dommermuth and Yue [START_REF] Dommermuth | A high-order spectral method for the study of nonlinear gravity waves[END_REF] (only the highest steepness ka = 0.4 is presented), we provide the reference η and φ from the model of Rienecker and Fenton [START_REF] Rienecker | A Fourier approximation method for steady water waves[END_REF] and evaluate the vertical velocity with the HOS scheme. The absolute error ǫ = max |W -W RF | is computed (W RF is the reference vertical velocity obtained from [START_REF] Rienecker | A Fourier approximation method for steady water waves[END_REF]). It is clear from this figure that HOS-ocean code allows an accurate evalu-ation of the vertical velocity at the free surface, even with a very steep wave (ka = 0.4), close to the wave-breaking limit. Exponential convergence rate is observed with respect to the two parameters N and M, thanks to the pseudo-spectral formalism of the model. Furthermore, note that thanks to the careful dealiasing procedure used in HOS-ocean, saturation observed in [START_REF] Dommermuth | A high-order spectral method for the study of nonlinear gravity waves[END_REF][START_REF] Skandrani | Nonlinear evolution of water surface waves: the frequency down-shift phenomenon[END_REF] is removed and better accuracy is achieved.

Finally, benchmarks allow to test the configuration with waves defined along x-axis, y-axis or 3D wave propagating at π/4 angle. Exact same results are obtained in the three set-ups.

Application: 3D irregular wave field

This section presents a typical application case of the simulation of non- • Directional spreading defined as β = 0.14

• Infinite water depth

• Lengths L x × L y = 20λ p × 40λ p ≃ 50 km 2 • Duration of simulation T stop = 250 T p ≃ 50 min • Discretization N x × N y = 256 × 256
• HOS order M = 3

• Time tolerance 10 -7

• Relaxation used with T a = 10 T p and n = 4

An example of free surface elevation at a given time-step is presented in Fig. 3. The 3D wave-field simulated has a small angular spreading (β = 0.14 rad), as can be seen with the free surface elevation. At the same time, frequency spreading is also limited (γ = 5). This is a typical configuration where the Benjamin-Feir instability is important during the propagation of such sea state. Benjamin-Feir directions, along which spectrum evolution is expected to be preponderant [START_REF] Dysthe | Evolution of a narrow-band spectrum of random surface gravity waves[END_REF]. This is a result with one realization (i.e. one random phase generation).

Thus, this is maybe not very clear that we recover correctly those instabilities. However, it has been shown in [START_REF] Ducrozet | Modélisation des processus non-linéaires de génération et de propagation d'états de mer par une approche spectrale[END_REF] that HOS-ocean solves correctly Note that post processing enables the analysis of the free surface elevation to extract evolution of interesting properties that characterize the nonlinear propagation of large-scale sea-states. For instance, Fig. 5 presents the evolution of kurtosis of free-surface elevation, as well as the normalized maximum wave height Hmax Hs observed during the simulation (defined as the maximum of each 3D wave detected with wave-by-wave analysis along mean and transverse direction). Note that this increase of kurtosis (i.e. the increase of the probability of freak waves' occurrence) is a typical nonlinear effect in such configuration due to third-order quasi-resonant interactions, see e.g. [START_REF] Fedele | On the kurtosis of deep-water gravity waves[END_REF].

Kinematics and pressure inside the domain

The study of the kinematics induced by a wave field is of particular importance for in-depth analysis of the wave phenomena as well as for the wave-structure interactions.Here is given an example of the kinematics and the source files easy. In addition, benchmarking procedures as well as typical post-processing examples are provided with the main code.

The code has been shown to be accurate and highly efficient thanks to the interesting convergence properties of pseudo-spectral methods. It has been highly validated on several configurations and one of the purpose of the open-source release is to encourage researchers to use it. Particularly, thanks to the post-processing provided, the coupling between HOS-ocean and CFD software for wave-structure interactions modeling should be straightforward. HOS-ocean, based on a pseudo-spectral formalism, uses FFTs which are performed with the FFTW library [START_REF] Frigo | The design and implementation of FFTW3[END_REF]. Parallelization of the program using MPI is consequently quite simple since FFTW provides also parallel FFT implementations. Preliminary works have already been done in [START_REF] Ducrozet | Analysis of freak waves formation with large scale fully nonlinear high order spectral simulations[END_REF] and show promising results. This will be included in the open-source version in the next future.

Finally, note that a HOS model dedicated to the simulation of a Numerical Wave Tank has been developed at ECN. It intends to reproduce a wave basin and its specifities: lateral reflective walls, wavemaker and absorbing beach.

The resulting model named HOS-NWT has been presented in [START_REF] Ducrozet | A modified high-order spectral method for wavemaker modeling in a numerical wave tank[END_REF] and is also available as open-source release [1].

) 1. 6 .

 6 Post-processingThis section describes briefly the different possibilities for post-processing of HOS-ocean outputs offered within open-source release of the code. For a fine analysis of physical processes at play during the nonlinear wave propagation, it appears useful to compute complementary information in addition to the two unknowns solved by HOS model η and φ.1.6.1. Wavefield analysisHOS-ocean enables the study of the evolution in time of the free surface elevation η. It appears interesting during post-processing to give access to some basic information useful for any analysis, such as the different moments of this free surface elevation: mean, standard deviation, skewness, kurtosis.

Figure 1 :

 1 Figure 1: Schematic view of time-stepping in HOS-ocean

Hs_real

  5), different output files are created. They are located in a specific folder Results with a specific format dedicated to Tecplot visualization. Input file also defines if outputs are dimensional or non-dimensional quantities. Following files may be created • 3d.dat gives the 3D free surface quantities (η and φ) as function of time i_case Choice of initial condition (regular wave, irregular . . . ) Geometry of the horizontal domain xlen Length of the domain along x-direction ylen Length of the domain along y-direction (if 3D case) Time integration T_stop Duration of the simulation f_out Sampling frequency for output files toler Tolerance of Runge-Kutta Cash & Karp scheme (adaptative time-stepping) n n parameter in relaxation scheme on nonlinear terms in free surface boundary conditions (cf. Dommermuth [7]) Ta T a parameter (duration) in relaxation scheme on nonlinear terms in free surface boundary conditions (cf. Dom-Significant wave height (m) gamma γ factor in JONSWAP spectrum (Eq. (16)) beta β factor in the definition of directional spreading (Eq. (17)) Definition of output files tecplot Tecplot version used for output files format i_out_dim Output are dimensional (=1) or non-dimensional (=0) i_3d File 3d.dat created if equal to 1 i_a_3d File a_3d.dat created if equal to 1 i_2d Files 2d.dat or 2dpt.dat created if equal to 1 or 2 respectively i_prob File probes.dat created if equal to 1 i_sw File modes_HOS_SWENSE.dat created if equal to 1 • a_3d.dat gives the modal amplitudes of free surface quantities (η and φ) as function of time • 2d.dat gives the 2D free surface quantities (η and φ) as function of time • 2dpt.dat gives the 2D+t free surface elevation • vol_energy.dat gives the temporal evolution of volume and energy during the simulation • modes_HOS_SWENSE.dat gives the file containing modal information of volumic HOS field for visualization,physical interpretation of wavefield. . . or possible coupling between HOS and wave-structure interactions model for instance (SWENSE method): see post-processing section • probes.dat gives the free surface elevation at specific location given in file prob.inp • phase_shift.dat gives the phase shift between computed solution and reference Rienecker & Fenton solution 2.4. Post-processing Post_processing binary will analyze the HOS-ocean outputs and perform some analysis with a given input file named input_post_process.dat. This input file is detailed in Tab. 6 which describes the different constants to specify. Choice of post-processing options i_ana Wavefield general features (=1) with 2D wave-by-wave analysis (=2) or 3D (=3) i_card Evaluation of velocity and pressure maps inside fluid domain with cartesian grid (=1) or boundary fitted mesh (=2) T_start Starting time of post-processing T_stop Stopping time of post-processing Information on kinematics output x_min Minimum x location of velocity/pressure maps x_max Maximum x location of velocity/pressure maps y_min Minimum y location of velocity/pressure maps (if 3D) y_max Maximum y location of velocity/pressure maps (if 3D) z_min Minimum z location of velocity/pressure maps z_max Maximum z location of velocity/pressure maps (if cartesian grid) i_zvect Number of points in the z direction Input files information tecplot Version of tecplot used to write following input files file_3d Path defining the file used for wavefield analysis (3d.dat generated by HOS-ocean) file_mod Path defining the file used for velocity/pressure map (modes_HOS_SWENSE.dat generated by HOS-ocean)

  Different moments of the free surface are always computed -i_ana = 2 or i_ana = 3 provides in addition in this file H 1/3 , H max , A max (the maximum crest amplitude) and the number of freak waves detected • freak_waves.dat gives the free surface profile of each detected freak wave in the domain • Caract_freak_waves.dat gives the characteristics of detected freak waves (size, amplitude . . . ) 2.4.2. Kinematics inside fluid domain HOS-ocean is able to create a modal description of some volume information in the fluid domain. All necessary information to compute velocities, pressure, acceleration (or whatever derivatives of these quantities) are included in the file named modes_HOS_SWENSE.dat generated by HOS-ocean code. Different output files are created and give examples of representation of velocities and pressure field (i_card = 0).

Fig. 2 Figure 2 :

 22 Figure 2: Convergence on the evaluation of vertical velocity W for highly nonlinear regular wave ka = 0.4

  linear 3D irregular wave field. Physical and numerical parameters of the simulation are the following • JONSWAP spectrum, H s = 11 m, T p = 12.5 s (mean steepness s = k p ā = kpHs 2 √ 2 = 0.1), γ = 5

Figure 3 :

 3 Figure 3: Initial free-surface elevation

Figure 4 :

 4 Figure 4: Evolution of wave spectrum during nonlinear propagation

Figure 5 :

 5 Figure 5: Evolution of kurtosis and normalised maximum wave height Hmax Hs during simulation

Figure 6

 6 Figure 6 presents the dynamic pressure field and the horizontal velocity field at the moment of the extreme event.

Figure 6 :

 6 Figure 6: Dynamic pressure (top) and horizontal velocity (bottom) of an extreme wave appearing in open ocean.

Table 2 :

 2 Dependency matrix for HOS-ocean

	Specific post-processing source files
	Post_processing.f90	Main program for post-processing HOS-ocean
		results
	analysis_wavefield.f90	Analysis (wave-by-wave) of the wavefield and
		freak waves detection
	input_post_process.f90	Reading of input file for post-processing
	output_post_process.f90	Initialization of output files and writing at
		each time-step
	read_files.f90	Reading of files that will be post-processed
	reconstruction.f90	Reconstruction of the volumic wavefield infor-
		mations
	variables_post_process.f90 Definition of the different global variables used
		in post-processing

Table 3 :

 3 List of source files used for post-processing

	Specific benchmarking source files
	check_DYue.f90	Main program for benchmarking HOS core
		routines
	variables_3D_DYue.f90 Definition of the different global variables used
		in benchmarking

Table 4 :

 4 List of source files used for benchmarking

	binary makedepf90 to extract automatic dependencies. Following commands
	can be used
	• Main program: make command executed at the root of the project will
	create the binary HOS-ocean
	• Benchmarking: make -f makefile_test_DYue will generate the bi-
	nary check_W_DYue

• Post-processing: make command executed in directory sources/PostProcessing will create the binary Post_Processing Compilation has been tested on different Unix/Linux platforms as well as in Windows environment.

Table 6 :

 6 Description of input file parameters for post-processing