Variational Approaches
To Water Wave Simulations
Elena Gagarina
ar
iat
iona
l A
pproache
s
To
W
at
er
W
ave
S
im
ul
at
ions
Elena
G
agar
ina
ISBN: 978-90-36-3754-4
To Water Wave Simulations
Elena Gagarina
Chairman and Secretary:
Prof. Dr. P. M. G. Apers University of Twente
Promotors:
Prof. Dr. Ir. J. J. W. van der Vegt University of Twente
Prof. Dr. Ir. O. Bokhove University of Leeds
Members:
Prof. Dr. Ir. E. W. C. van Groesen University of Twente
Prof. Dr. H. J. Zwart University of Twente
Prof. Dr. Ir. R. H. M. Huijsmans Delft University of Technology
Prof. Dr. A. E. P. Veldman University of Groningen
Referee:
Dr. T. Bunnink MARIN
University of Twente, MaCS Group
P.O. Box 217, 7500 AE, Enschede, The Netherlands.
ISBN: 978-90-36-3754-4 DOI: 10.3990/1.9789036537544
http://dx.doi.org/10.3990/1.9789036537544 Printed by CPI-W¨ohrmann Print Service, Zutphen.
TO WATER WAVE SIMULATIONS
DISSERTATION
to obtain
the degree of doctor the University of Twente, on the authority of the Rector Magnificus,
Prof. Dr. H. Brinksma,
on account of the decision of the Graduation Committee, to be publicly defended
on Friday 3 October 2014 at 12:45
by
Elena Vitalyevna Gagarina
born on 25th of November 1987 in Kiev, USSR
Prof. Dr. Ir. J. J. W. van der Vegt (promotor) Prof. Dr. Ir. O. Bokhove (promotor)
The research presented in this dissertation was carried out at the Mathematics of Com-putational Science (MaCS) group, Department of Applied Mathematics, Faculty of Elec-trical Engineering, Mathematics and Computer Science pf the University of Twente, The Netherlands.
The work has been a part of the STW project ”Complex wave-current interactions in a numerical wave tank” and NWO project ”Compatible Mathematical Models for Coastal Hydrodynamics”. The support from the Technology Foundation STW, the Netherlands Organization for Scientific Research (NWO), the Maritime Research Institute Netherlands MARIN is gratefully acknowledged.
The project has been associated with the J. M. Burgers Research School for Fluid Dy-namics and the Nanotechnology Research institute of the Univeristy of Twente MESA+.
1 Introduction 1
1.1 Motivation . . . 1
1.2 Objectives . . . 4
1.3 Understanding water waves . . . 5
1.3.1 Experiments . . . 5
1.4 Modeling water waves . . . 6
1.4.1 Mathematical wave models . . . 6
1.4.2 Numerical methods . . . 7
1.4.3 Numerical wave tank . . . 8
1.4.4 Variational dynamics . . . 9
1.5 Outline . . . 10
2 Horizontal circulation and jumps in Hamiltonian wave models 13 2.1 Introduction . . . 13
2.2 New water wave model . . . 15
2.2.1 Variational principle . . . 15
2.2.2 Reduction of Hamiltonian dynamics . . . 18
2.2.3 Hamiltonian dynamics of new water wave model . . . 19
2.3 Shallow water and potential flow limits . . . 22
2.4 Hamiltonian Boussinesq reductions of new model . . . 23
2.4.1 Variational Boussinesq model . . . 24
2.4.2 Green-Naghdi limit . . . 26
2.5 Jump conditions for bores . . . 28
2.5.1 1D Jump conditions for shallow water equations . . . 29
2.5.2 2D Jump conditions for new water wave model . . . 33
2.5.3 2D Jump conditions for shallow water equations . . . 38
2.5.5 1D Jump conditions for variational Boussinesq and Green-Naghdi
equations . . . 40
2.6 Conclusions . . . 40
2.7 Appendix . . . 42
2.7.1 Integrals . . . 42
2.7.2 Variational Reynolds’ transport theorem . . . 43
3 Symplectic Time Discontinuous Galerkin Discretizations for Hamiltonian Sys-tems 49 3.1 Introduction . . . 49
3.2 Dynamics of a Hamiltonian system with one degree of freedom . . . 51
3.3 Variational discontinuous Galerkin time discretization of a Hamiltonian system . . . 52
3.3.1 Discrete functional . . . 52
3.3.2 First order variational time discretization . . . 54
3.3.3 Second order variational time discretization: St¨ormer-Verlet . . . 55
3.3.4 Second order time discretization: Modified symplectic midpoint scheme . . . 58
3.4 New variational time integrators . . . 60
3.4.1 Explicit second order accurate time integration scheme . . . 61
3.4.2 Third order scheme . . . 62
3.5 Applications to non-autonomous systems . . . 69
3.5.1 Damped oscillator . . . 71
3.5.2 Forced oscillator . . . 74
3.6 Conclusions . . . 78
3.7 Appendix . . . 78
3.7.1 Continuum harmonic oscillator . . . 78
3.7.2 St¨ormer-Verlet scheme . . . 80
3.7.3 New explicit second order variational time integration scheme . . 82
3.7.4 Extended continuum system . . . 89
4 Variational space-time (dis)continuous Galerkin method for nonlinear free surface water waves 93 4.1 Introduction . . . 94
4.2 Variational Description of Nonlinear Inviscid Water Waves . . . 96
4.3 Variational Finite Element Formulation . . . 99
4.3.1 Tessellation . . . 99
4.3.2 Function spaces and approximations . . . 100
4.3.3 Space-plus-time formulation . . . 102
4.3.4 Dynamics . . . 109
4.4 Numerical Results . . . 113
4.4.1 Verification . . . 113
4.5 Application of new third order symplectic time integrator . . . 123
4.5.1 Nonlinear potential flow water waves . . . 127
4.6 Conclusions . . . 132
4.7 Appendix . . . 133
4.7.1 Evaluation of finite element matrices . . . 133
4.7.2 Variations of the kinetic energy . . . 134
5 Conclusions 137 5.1 Recommendations . . . 137 Summary 139 Samenvatting 141 Bibliography 143 Acknowledgments 157
Chapter
1
Introduction
1.1
Motivation
A wise man once said that one can observe fire and water forever. The elements are fascinating indeed, what is their secret? Perhaps, the reason is the subtle interplay between life and death. We don’t exist without water, but it can also easily kill. Water has always been a crucial part in the progress of mankind. The largest, richest, and most important cities are those that have easy access to a convenient water transportation route. This allows the exchange of goods and knowledge between civilizations. For many centuries explorers and merchants were sacrificing their lives in a battle with nature to explore new shipping routes and distant countries. It is no surprise, that understanding the behavior of water waves became of great importance to humanity. This was also motivated by the many natural disasters that had a great impact on society.
One of the most dangerous events for coastal towns is a tsunami wave. The Santorini eruption and subsequent tsunami are believed to be the cause of the Minoan civilization downfall. The disastrous 1755 Lisbon earthquake and tsunami, Fig. 1.1(a), influenced the mind of many Europeans. For instance, Voltaire, Rousseau, and Kant were trying to comprehend in their works a power so supernatural and unfriendly. An extensive overview of floods, tsunamis and other natural disasters can be found in [117]. These events have also extensively motivated the visual arts, see [84].
The recent tragic events of the Boxing Day Tsunami in 2004 and the 2011 Tohoku Tsunami reveal how unprotected modern society still can be when disaster strikes. Possi-ble ways of protection against disastrous waves, such as dikes, wave breakers, and man-grove forests, are therefore of great importance. A novel way to protect people from a devastating tsunami is an alarm system that is based on accurate and fast methods to simulate wave propagation and run-up in the near coastal zone.
Figure 1.1: (a) A copper engraving showing Lisbon, devastated by the earthquake
and subsequent tsunami in 1755. Original in: Museu da Cidade, Lisbon. Source:
Wikipedia.org. (b) ”The Ninth Wave” – a famous painting of Russian Armenian marin-ist Ivan Aivazovsky. Collection of the State Russian Museum, St. Petersburg. Source: Wikipedia.org.
Figure 1.2: (a) Tsunami wave enters Miyako city in Japan in 2011. Authorship: Mainichi Shimbun /Reuters. Source: reuters.com. (b) The Draupner wave or New Year’s wave is one of the first measured rogue waves. It occurred at the Draupner platform in the North Sea on 1 January 1995. Source: Wikipedia.org.
and sometimes extreme waves, present a serious danger to ships and offshore structures. Recently, special attention has been given to rare extreme wave events, called freak or rogue waves. In 1978 the German carrier M¨unchen sank. The investigation into this acci-dent revealed that the ship was most likely damaged by a wave of up to 30 m height. In 1942 the passenger ship ”Queen Mary” was transporting troops from the USA to Europe. Sixteen thousand people were on board when a wave of 28 meters height stroke the ship’s side and tilted it by 52 degrees [95]. Those people were just lucky to survive. Never-theless, the actual causes of ship damages and losses in deep water often remain unclear. Freak waves, which used to be a scary fisherman tale, now are proven to exist [76]. One of the first recorded freak waves is the Draupner wave, see Fig. 1.2(b), it had a
maxi-mum height of 25.6m, which more than twice exceeded the significant wave height. The southern tip of Africa is also a place where extreme waves happen more frequently, see Fig. 1.3(b). One of the possible mechanisms to explain this is the special kind of wave-current interaction occurring in this region due to the Agulhas wave-current, which is opposite to the direction of predominant winds.
In order to understand extreme wave phenomena it is important to define which sea states can result in the occurrence of extreme waves. A second important question is how to minimize the risks from freak waves and severe storms to ships and offshore structures, see Figs. 1.3 and 1.4(a).
Figure 1.3: (a) Rogue wave estimated at 18 meters moving away from a ship after crashing into it a short time earlier. Gulf Stream off Charleston, South Carolina. Source: NOAA. (b) Norwegian tanker ”Wilstar” survived a freak wave encounter in 1974 near the east coast of South Africa. Source: Photo collection of the European Space Agency with credits to H. Gunther and W. Rosenthal.
These questions are especially important for the design of new ships and oil plat-forms. In particular, this requires knowledge on the wave-floating structure interaction, since insurance companies want to make sure that the structure meets the safety require-ments. The dramatic experience of Hurricanes Ivan, Rita and Katrina stimulated extensive research in extreme waves and has resulted in changes in construction requirements for offshore structures. Several extensive studies contributed to these new regulations. For instance, the project Extreme Seas [134] studied wave conditions in the northern hemi-sphere and consequently, proposed new safety criteria. Also, the study [12] predicts in-creasingly extreme wave heights in the coming years due to climate change and also gives recommendations for an update of construction requirements in tanker design.
Another aspect of (extreme) wave-structure interaction is the protection of cities. Wave surges, induced by hurricanes, might damage coastal cities, as is, unfortunately, observed over and over again, Fig. 1.4(b). Even without extreme events the knowledge of wave impact on a structure is important. For example, the streets of fragile Venice are suffering from ship-induced waves, which makes it necessary to restrict the boat speed.
Figure 1.4: (a) Mars oil platform damaged by Hurricane Katrina in 2005. Source: Wikipedia.org. (b) Hurricane Carol induced huge waves in 1954 that destroyed hundreds of houses in Connecticut. Source: NOAA.
1.2
Objectives
Given the importance of various types of water waves to coastal regions and engineer-ing, we aim to study the mathematical aspects of free surface flows in this thesis. There are many mathematical models describing water waves under different assumptions. We can study those models analytically or choose a numerical method. One of the tasks we will concentrate on is the construction of a numerical wave tank. This is a numerical ap-proximation of waves in a model basin, which allows the simulation of realistic sea states. The challenge here is to include complicated wave phenomena in the model.
A detailed discussion of various mathematical aspects of extreme wave modeling is given in [76], [120], [119]. In this thesis an important topic will be wave breaking, as happens in the near shore region, and is particularly important for tsunami studies. The difficulty here are the complex turbulent processes in the overturning wave, which will be approximated as a discontinuous flow. Another important topic in this thesis is wave focusing, which is important for rogue waves and stormy waves. Wave run-up and the interaction of waves with structures will not be considered here, but they will be part of future research.
The presence of a free surface and moving boundaries, such as wave makers, compli-cates the numerical modelling of wave phenomena. To cope with this we use specially designed variational finite element methods. Another important topic is mass and energy conservation, which are particularly important for long time simulations, as is frequently required in a numerical wave tank. We address these issues through the development of a (discrete) variational approach, which is directly associated with the conservation laws, in particular energy conservation. In the simulation of focusing waves it is also impor-tant to ensure numerical stability for large amplitude waves, which motivated us to use
symplectic time integrators.
1.3
Understanding water waves
It is no surprise that the brightest minds were attracted by the secrets of water waves. The ”Eureka” exclamation from Archimedes became a symbol of scientific excitement. An overview of the development of water wave theory can be found in [33], [36]. In this section we will give a summary of the main techniques used to study water waves.
1.3.1
Experiments
There are several ways to study water waves. The oldest way is of course observa-tion. This is how Archimedes discovered his law and John Scott Russel discovered the soliton type of wave propagation. Field surveys and data collection of extreme events are very important for risk analysis, validation, and to provide realistic sea states in model experiments and computations [110], [137]. Sometimes data collection on a rare event, especially a historic one, reminds of a crime investigation with witness interviews and a detailed reconstruction of the event [114]. One can even try to look for traces more than 240 hundred years after the event, as did the authors of the prominent study [56] on the 1771 Meiwa Tsunami event.
As a variant of observation we have experiments. Experiments cover a wide range of wave phenomena. One can study wave focusing [25], [116], wave breaking [26], [80], [141], wave run up [19], [71], tsunami generation [43], [108], or wave-structure interac-tion, [28], [122]. A recent breakthrough was the observation of a freak wave in a labora-tory, with a shape that corresponds to the analytic solution of the nonlinear Schr¨odinger equation – a so-called Peregrine breather [24].
There are many laboratories with various facilities all over the world. An overview on issues and specifications of laboratory experiments can be found in [37]. In these facilities prototype offshore structures and models of ships are tested in wave conditions mimicking a realistic sea state. Even though the concept is very natural, there are a number of difficulties. Not every condition can be modeled in a laboratory. For example, a floating oil platform in deep water is difficult to scale properly, due to the difference in the size of the platform, the risers and the depth of the ocean. Also, experiments require the construction of a model and the set-up of the experiment, which takes time, and needs to be accounted for in the design process of ships and offshore structures.
An example of a model basin is presented in Figure 1.5, which shows the wave basin of the Maritime Research Institute Netherlands (MARIN). In this advanced experimental facility it is possible to generate various wave forms by the movement of wave makers, and to study the wave impact on floating and fixed objects. A crucial question is then how one should move the wave maker in order to obtain the required wave phenomena. This is a non-trivial problem that has resulted in a number of studies, e.g. [1], [82], [89],
[135], [152] and others. The special problem of generating focusing wave trains, which is important to study extreme wave events, is considered in [68], [70].
Ship model length of 3 - 6 m.
Floating structures of any kind, size depending on water depth and wave conditions (usually between 0.2 m for buoys and 4 m for platforms).
Cross section
Test capabilities
Offshore structure models, fixed, moored or controlled by dynamic positioning in waves, wind and current.
Captive or free sailing manoeuvring tests in shallow water.
For more information please contact the department Offshore T +31 317 493 465
E Offshore@marin.nl
Figure 1.5: The MARIN facilities allow the creation of various wave types by moving wave makers at the side of the model basin.
1.4
Modeling water waves
1.4.1
Mathematical wave models
There are a number of mathematical models that describe water waves under differ-ent assumptions. A broad overview of wave models is given in the promindiffer-ent books of Whitham [155] and Landau and Lifschitz [90]. The Navier-Stokes equations describe a general viscous incompressible fluid, with the Euler equations as inviscid limit. These equations describe general fluid flow, including a free surface, but for many types of free surface waves they are too complicated, therefore further simplifications are useful for practical applications. These simplifications can be, for example, the potential flow sys-tem, the shallow water equations, the Boussinesq equations, and many others.
Since there is a variety of wave models, one has to carefully consider the wave phe-nomena the model has to describe. We will discuss therefore the advantages and disad-vantages of some widely used wave models. An important wave model is the shallow water equations, which are relatively simple and easy to solve numerically. This wave model allows discontinuous solutions – so called bores – as a simplification of compli-cated turbulent breaking waves. While mass and momentum are conserved across the discontinuity, energy is decreasing. This corresponds to the observations of breaking waves in [155] – where the water sprays and splashes while the wave is overturning, and, therefore, dissipates energy.
As described in [121], when waves are breaking, energy is lost, but momentum is conserved. This results in a rise in horizontal shear or vertical vorticity, which becomes important when simulating breaking waves near the shore with two-dimensional models, such as the shallow water equations.
Among the downsides of the shallow water equations is the absence of dispersion (the property that waves with different wave length travel with different speed). In particular, waves are breaking when the nonlinearity outweights dispersion. Due to the absence of dispersion in the shallow water equations waves computed with this model tend to break too early in comparison with observations.
A second widely used wave model is the potential flow model. The potential flow equations are obtained from the Euler equations under the assumption that the velocity field is irrotational (the velocity is the gradient of a scalar function called the potential), [73]. The potential flow model includes accurate dispersion, which makes the system attractive for wave focusing problems and the simulation of waves generated by a wave maker.
A unified wave model, embracing both the shallow water and potential flow system is useful when wave phenomena in various regions (e.g. open sea and near coastal region) need to be modeled simultaneously. Such a model was proposed by Cotter and Bokhove [30]. This model was obtained from an Eulerian variational principle with extended Cleb-sch variables [158], that only depend on the horizontal coordinates. By construction the velocity field consists of a three-dimensional potential flow part and a horizontally dependent vortical part. Under appropriate restrictions the new water wave model can be reduced to the shallow water equations, the potential flow model, and the variational Boussinesq model of Klopman et al., [82].
The Euler equations, which include three dimensional vorticity, is the next level of wave models. In Nurijanyan [112], a variational model for the rotating Euler equations is provided, and a linearised version of this model is applied to inertial waves.
A mathematical wave model can be solved either analytically or numerically. In prac-tice, analytical methods can only be used to solve a wave model in special occasions. For potential flow a pseudo-analytical solution is provided in [47] and [55], for the shallow water equations solutions can be found in [4] and [140], and for the Euler equations [15] and [102]. There are also some known exact solutions for the Navier-Stokes equations [45] and [150], but practical wave phenomena might require a computation with numeri-cal methods.
1.4.2
Numerical methods
Considering the difficulties in obtaining analytical solutions for the various wave mod-els, numerical methods are of great practical importance. Since the power of computers is continuously growing, increasingly more sophisticated numerical models are becoming available. Still, there is a pay-off between the complexity of the wave model and its practi-cal use. In particular, one has to define the necessary phenomena that need to be included
in the simulation model and to consider the required computational power. Imagine there is an underwater earthquake with the danger of a tsunami. It is then important to compute the possible wave shapes generated by the motion of the ocean bottom, to simulate their propagation and to define as soon as possible areas under threat. In this case, the simplest model, the shallow water model, will do the job. On the other hand, when one works with focusing waves or wave–current interactions, where dispersion or three dimensional vorticity play a role, more sophisticated wave models need to be used.
Provided the main interest of this research, we give now a brief overview of numerical methods suitable for the simulation of water waves in a model basin.
1.4.3
Numerical wave tank
Apart from studying wave phenomena in the ocean or near coastal region, wave mod-els are also important to study wave generation in a model basin or even to replace ex-periments in a model basin all together. This motivates the development of a numerical wave tank (NWT), which simulates waves in a laboratory (or in the ocean). In 1976 Longuet-Higgins and Cokelet simulated an overturning wave using potential flow theory with a mixed Eulerian-Lagrangian method in combination with a boundary integral equa-tion formulaequa-tion. Almost simultaneously, in 1977, Faltinsen [46] implemented a similar method for a floating body in water waves. One of the first to employ a finite element method for the nonlinear potential flow equation were Wu, Wang and Eatock Taylor [41] in 1994. A spectral method using a Dirichlet to Neumann operator was introduced by Craig and Sulem in 1993, [32]. For a review of the numerical methods for nonlinear water waves before 2000, see [77] and [142].
For complicated breaking waves the meshless Smoothed Particle Hydrodynamics (SPH) methods are an interesting new development and currently under active study, [67] and [143]. Even though SPH is a convenient method to simulate wave-structure interac-tion or extreme waves, the method suffers from substantial artificial dissipainterac-tion, which does not make it suitable for accurate long time simulations or wave propagation studies. One of the first three-dimensional numerical wave tanks was developed by Grilli et al [58]. It contains an arbitrary bottom topography, sloping beach and the possibility to include wave makers. The nonlinear potential flow equation is solved with a combination of a boundary element method (BEM) and a mixed Eulerian-Lagrangian technique to compute the free surface motion. This NWT was used both for the study of breaking waves [63], and freak waves [51]. Another three dimensional NWT based on BEM is presented in [127], [128].
An alternative to the widely used boundary element methods for nonlinear water waves is provided by finite element methods, which currently are under active devel-opment [101], [145], [151]. These methods have a number of advantages. In particular, they allow to deal with complicated moving boundaries. There are also other benefits for finite element methods: the numerical discretization results in sparse matrices suitable for fast solvers, they have the potential to be extended to more complicated flow models,
such as Euler and the Navier-Stokes equations, which is not possible with BEM. Dis-continuous finite element methods, which result in a local, element based discretization with only a weak coupling to neighboring elements are particularly useful for wave prob-lems requiring mesh adaptation using local mesh refinement and parallel computations [138]. For problems with moving boundaries, such as free surface waves, it is beneficial to use a space-time approach, which has attracted much attention recently, [23], [75]. This approach considers time as an extra dimension, making a time-dependent n-dimensional problem an n+1 dimensional steady state problem. In a space-time framework it is easier to satisfy the geometric conservation law on time-dependent domains, where the interior mesh movement is induced by the free surface and boundary motion, see Lesoinne and Farhat, [94]. Space-time DG finite element methods for linear free surface wave prob-lems were developed in Ambati [3], and for nonlinear flows in Van der Vegt and Xu [144]. Understanding time as an extra dimension to the spatial problem also leads to the construction of novel time stepping methods using time-discontinuous elements. For an extensive overview of numerical wave tanks, see [40].
In order to facilitate the validation of numerical simulation models for water waves MARIN organized in 2010 a special workshop [22], for which it conducted a number of experiments in its facilities. These experiments included wave focusing over a flat bottom and irregular waves over an uneven bottom. These experiments provided a wealth of data, which was used in [60], [91] for a model validation and will also be extensively used in this study.
1.4.4
Variational dynamics
Many water wave problems can be described by a variational principle, which pro-vides a description of the state of a physical system as the extremum of a functional. Originally this idea was explored in problems of finding the shortest path between two points or the least time of light is travelling between two points (the principle of least action). Euler was first to establish a connection between the principle of least action and a tendency of physical bodies to occupy a state with minimal total potential energy. Lagrange and, later, Hamilton reformulated classical mechanics using an energy-related functional, now referred to as Lagrangian. For an extensive overview of mathematical methods of classical mechanics, see [6].
In a variational principle, the complete system is expressed as a single functional – a Lagrangian – which is the difference between potential and kinetic energy. Within a variational framework it is relatively easy to introduce variable transformations, establish conservation laws [10] and stability theorems, and it is also a good starting point for a numerical discretization. The representation of the dynamics using a variational principle is concise and elegant.
The variational formulation also provides a straightforward way to obtain simplified systems that preserve the symmetry of the original, more complicated system. An ex-ample related to wave phenomena is the simplification of the velocity potential in Miles’
formulation as is done by Klopman [82]. A wide range of consistent approximations is discussed in [31].
Hamiltonian mechanics is a reformulation of Lagrangian mechanics to investigate the symmetry properties of a system. There is a direct connection between the variational principle and the symmetry properties of Hamiltonian systems and conservation laws as is established by Noether’s theorem [111].
The first developments in variational fluid dynamics were made by Clebsch in 1859 [29], Hargreaves in 1908 [66], and Bateman in 1944 [8]. They considered the pressure as a Lagrangian density, but did not include the boundary conditions at the free surface.
Inspired by Whitham [153], [154], Luke [98], introduced a variational principle for free surface waves in a bounded basin. He also sketched the possibility to add vorticity to the velocity field. Independently, Zakharov [159], Broer [20], and Miles [106], proved the Hamiltonian structure of the nonlinear potential flow equations. The pair of
canon-ically conjugate variables for potential flow free surface waves are η and φS – the free
surface elevation and the velocity potential at the free surface, which act as generalized coordinate and momentum. Variational Boussinesq type approximations were developed by Broer [20], Miles and Salmon [107], Milder [105], and Klopman [82]. For a review, see Finlayson [48] and Grimshaw [59].
The variational formulation allows a straightforward numerical discretization. Indeed, one just needs to approximate the variables with e.g. a spatial finite element or spectral representation. The discrete variational principle will then inherit the conservation of en-ergy and phase space of the continuous variational formulation. One has to be, however, careful with the time integration, otherwise the variational structure will not be preserved and the numerical discretization might suffer from numerical instabilities. For this pur-pose symplectic time integrators for Hamiltonian systems were developed [64]. These time integrators preserve conservation laws and ensure that the discrete solution has no energy drift and its fluctuations are bounded if the original system is energy conservative. Recent developments in the field of (discrete) Hamiltonian fluid dynamics can be found in [17], [27], [78], [82].
1.5
Outline
The outline of this thesis is as follows. In Chapter 2 we consider a new variational water wave model, developed by Cotter and Bokhove [30]. This model inherits the ben-efits of the shallow water model and the potential flow model. It has a three-dimensional velocity field consisting of the full three-dimensional potential velocity field plus extra, horizontally dependent, components. This construction implies that only the vertical
vor-ticity component is nonzero1. The model includes dispersion and allows breaking waves.
In this thesis the Hamiltonian structure of the model is shown and the Poisson bracket to
1It is pointed out to us by T. J. Bridges in 2012, that an absence of the horizontal components of vorticity
describe the evolution of the system is obtained, see also [54]. The model reductions to a Green-Naghdi type of system and to the Boussinesq model of Klopman [82] are demon-strated under certain assumptions. Also, a variational approach to derive jump conditions to represent bores is presented. The approach is based on the work of Wakelin [148], where the domain is split into two parts by a bore and jump conditions at the bore are obtained. When the wave model is reduced to the shallow water limit the jump relations at the bore reduce to the well-known shallow water jump relations [90], [155]. Also, a variational analog of the Reynolds transport theorem is proved.
In Chapter 3 we use a variational approach to derive novel time integration meth-ods. The time integration method is based on a discontinuous Galerkin finite element discretization in time, with a specially derived numerical flux. With this approach we were able to derive the well known symplectic Euler, symplectic modified midpoint and St¨ormer-Verlet time integration methods. Also, new time integrators were obtained. For all novel variational time integrators the order of accuracy, linear stability conditions, symplecticity, and dispersion of the schemes are analyzed. The original motivation for this research was to obtain stable time integration methods for forced systems, and the newly developed techniques are applied to forced pendulums and, later in Chapter 4, for waves generated by a wave maker. Also, damped systems are considered.
In Chapter 4 we consider waves generated by a wave maker in a numerical wave tank [52]. The model we use is the nonlinear potential flow equation. The numerical method is a mixed space continuous - time discontinuous Galerkin method with linear basis functions, providing second order of accuracy in space and time. A comparison with laboratory experiments provided by the Maritime Research Institute Netherlands (MARIN) is made and found to be satisfactory. Also in Chapter 4 we combine the novel third order time integrator developed in Chapter 3 and the described variational numerical method to compute nonlinear water waves. In this autonomous test case the numerical method shows promising results.
Finally, conclusions and recommendation for future work will be discussed in Chap-ter 5.
Chapter
2
Horizontal circulation and jumps in
Hamiltonian wave models
We are interested in the modeling of wave-current interactions around surf zones at beaches. Any model that aims to predict the onset of wave breaking at the breaker line needs to capture both the nonlinearity of the wave and its dispersion. We have there-fore formulated the Hamiltonian dynamics of a new water wave model, incorporating both the shallow water and pure potential flow water wave models as limiting systems. It is based on a Hamiltonian reformulation of the variational principle derived by Cot-ter and Bokhove [30] by using more convenient variables. Our new model has a three-dimensional velocity field consisting of the full three-three-dimensional potential velocity field plus extra horizontal velocity components. This implies that only the vertical vorticity component is nonzero. Variational Boussinesq models and Green-Naghdi equations, and extensions thereof, follow directly from the new Hamiltonian formulation after using sim-plifications of the vertical flow profile. Since the full water wave dispersion is retained in the new model, waves can break. We therefore explore a variational approach to derive jump conditions for the new model and its Boussinesq simplifications.
The content of this chapter is presented as a journal article [54] written in collaboration with J.J.W. van der Vegt and O. Bokhove.
2.1
Introduction
The beach surf zone is defined as the region of wave breaking and white capping between the moving shore line and the (generally time-dependent) breaker line. Let us consider wave propagation from deeper water to shallow water regions. The start of the surf zone on the offshore side is at the breaker line where sustained wave breaking
begins. It demarcates the points where the nonlinearity of the waves becomes strong enough to outweigh dispersion. The waves thus start to overturn. From the point of breaking till the shore, the waves lose energy and generate vorticity. A mathematical model that can predict the onset of wave breaking at the breaker line will need to capture both the nonlinearity of the waves and their dispersion. Moreover the model has to include vorticity effects to simulate wave-current interactions.
Various mathematical models are used to describe water waves. A popular model for smooth waves in deep water is the potential flow model, but its velocity field does not include vorticity. In the near-shore region, vorticity effects are, however, important. When obliquely incident waves shoal in shallow water, steepen and break, a horizontal shear or vertical vorticity is generated. On semi-enclosed or enclosed beaches, this leads to an overall circulation induced by wave breaking. A classical hydraulic model for the surf zone is the shallow water model. The complicated, turbulent three–dimensional wave breaking is approximated in this model by discontinuities or so–called bores. These are special relations holding across the jumps connecting the right and left states of the flow. Mass and momentum are conserved across the discontinuity, while energy is not, as can be expected from observing the white capping zone of fine–scale splashes and sprays in the broken wave [154]. Shallow water waves are not dispersive, and these waves tend to break too early in comparison with real, dispersive waves. Boussinesq models include internal wave dispersion to a higher degree of accuracy, but dispersion always seems to beat nonlinearity. Therefore wave overturning tends to be prevented in these models. The variational Boussinesq model proposed by Klopman [82] could be a notable exception, but it is based on the Ansatz of potential flow. In three dimensions, a purely potential-flow model cannot be extended by inclusion of bores and hydraulic jumps as a simple model to represent wave breaking. The reason is that at least some vorticity has to be generated by bores that have non-uniformities along their jump line as was shown in [124], [121], and [123].
We therefore seek to develop a more advanced model that includes both the shal-low water approximation of breaking waves as bores and the accurate dispersion of the potential flow model. Such a model was obtained by Cotter and Bokhove [30] from a par-ent Eulerian variational principle with extended Clebsch variables, in which the vortical parts only depended on the horizontal coordinates. This restricts the vorticity to have a vertical component only. Extended Clebsch variables may, however, be less convenient algebraically and also yield a larger phase space of variables. We therefore reformulate this system in terms of surface velocity, velocity potential and water depth, and derive the Hamiltonian structure including its Poisson bracket. This new water wave model can be reduced to the shallow water equations, the potential flow model, and the Boussinesq model of Klopman [82] under corresponding restrictions. The Green-Naghdi equations emerge from the variational Boussinesq model by introducing a parabolic potential flow profile in the Poisson bracket, as well as another, columnar approximation of the velocity in the Hamiltonian.
ex-plored. It is inspired by the work of Waklein [148] for stationary shock or jump conditions for the shallow water equations. These results have been extended to moving shocks in shallow water based on the variational principles for the relevant Clebsch variables. Nat-urally, this approach allows us to obtain jump conditions for the new water wave model as well. The jump relations can be implemented in any system with a variational and Hamil-tonian structure, but not every system with a HamilHamil-tonian structure allows shocks or dis-continuities to persist in time. For example, it was shown in [42] that the Green-Naghdi system has an unsteady undular bore, i.e., an initial discontinuity in the free surface and velocity expands instantly into smooth undulations. It is therefore necessary to analyze the energy loss across jumps, and juxtapose this analysis between the original and our new extended Green-Naghdi system.
The outline of the chapter is as follows. First, a systematic derivation of the new Hamiltonian formulation will be given starting from a slightly adapted formulation of
the variational principle from [30] in §2.2. Subsequently, we show in §2.3 how it can
be reduced to limiting systems, such as the shallow water equations, the potential flow model, the Boussinesq model of Klopman [82] and an extended version of the
Green-Naghdi equations. In §2.5, a variational approach to derive jump conditions is given,
starting from the well-known Rankine-Hugoniot or jump conditions for the shallow water
equations. We end with conclusions in§2.6.
2.2
New water wave model
2.2.1
Variational principle
Consider an incompressible fluid at time t in a three-dimensional domain bounded by solid surfaces and a free surface, with horizontal coordinates x, y, and vertical coordinate z. The water depth is denoted by h = h(x, y, t). There exists a parent Eulerian variational principle for incompressible flow with a free surface. Its three-dimensional velocity field
~
U = ~U (x, y, z, t) = (u, υ, w)T, with (
·)T denoting the transpose, contains both potential
and rotational parts and is represented as ~
U = ∇φ + ~πj∇~lj (2.1)
through extended Clebsch variables: the velocity potential φ = φ(x, y, z, t), the three-dimensional fluid parcel label ~l = ~l(x, y, z, t) and the corresponding Lagrange multiplier
vector ~π = ~π(x, y, z, t), j = 1, 2, 3. Such a representation describes a velocity field
con-taining all three components of vorticity ∇× ~U . In order to avoid confusion, indices are
also introduced in (2.1) and the Einstein convention for repeated indices is used, such that ~
πj∇~lj= π1∇l1+ π2∇l2+ π3∇l3. This velocity representation is similar to the
expres-sion (3.9) in [130]. Also Lin [96] used two three-dimenexpres-sional vector Clebsch variables to introduce a vorticity for superfluids. As in [158], we see, that a pair of extended Clebsch vectors suffices for the generalized form to be complete.
When only the potential velocity field ∇φ(x, y, z, t) is used, there is no vorticity. In contrast, a shallow water velocity field includes the vertical component of vorticity. Similarly, in an Eulerian variational principle with planar Clebsch variables that only depend on the horizontal coordinates the vertical component of vorticity is retained. This component is constant throughout the whole water column and flows with helicity [88] are thus excluded by construction.
Cotter and Bokhove [30] derived novel water wave dynamics from the parent Eule-rian variational principle which includes two limits: Luke’s variational principle giving the classical potential water wave model and a principle for depth-averaged shallow water flows based on planar Clebsch variables. At least conceptually, the novel variational prin-ciple follows readily from the parent prinprin-ciple with two-dimensional label and multiplier
fields ~l and ~π depending only on the two-dimensional horizontal coordinates and time.
Hence, they no longer depend on the vertical coordinate z. In his prominent paper, Luke [98] has mentioned about the possibility of the introduction of Clebsch potentials into the variational principle for the rotational case. In contrast, we do not use Clebsch scalar variables, but extended vector Clebsch variables.
Extended Clebsch variables are, however, not convenient to work with. We therefore reduce the model to a more compact and conventional form. This reduction from six
variables{φ, h,~l, ~π} to four more conventional variables {φ, h, ~u∗} is undertaken in a
Hamiltonian setting. The latter variables involve a new velocity ~u∗, which is a suitable
horizontal velocity.
The variational principle of [30] has the following form:
0 = δ T Z 0 L[~l, ~π, φ, h]dt = δ T Z 0 Z ΩH b+h Z b ∂tφ + ~π· ∂t~ldzdxdy + Hdt, (2.2)
where the horizontal part of the domain is ΩH; the single-valued free surface boundary
lies at z = h(x, y, t) + b(x, y), with h(x, y, t) the water depth and b(x, y) a given, fixed
topography; and, t is time, its derivative is ∂t, and T a final time. The variational derivative
is defined as δ T Z 0 L[~l, ~π, φ, h]dt = lim→0 T Z 0 L[~l + δ~l, ~π + δ~π, φ + δφ, h + δh] − L[~l, ~π, φ, h] dt. (2.3) The component of the velocity with vortical parts is contained in
~v(x, y, t) = πj(x, y, t)∇lj(x, y, t), with j = 1, 2 (2.4)
and three-dimensional gradient ∇. Thus, the entire three-dimensional velocity field is represented by
~
combining the potential velocity ∇φ(x, y, z, t) and the planar velocity ~v(x, y, t). The relevant Hamiltonian, the sum of kinetic and potential energies, equals
H = H[~l, ~π, φ, h] = Z ΩH b+h Z b 1 2|∇φ + ~v| 2dz +1 2g (h + b) 2 − b2 − ghH0dxdy, (2.6)
with g the acceleration of the Earth’s gravity, and H0 a still water reference level. This
Hamiltonian contains the available potential energy, due to the additional subtraction of the rest level contribution, cf. [136].
As is shown in [30], the variational formulation of the new system is similar to Hamil-tonian classical mechanics, and becomes
δφ : ∇2φ + ∇ · ~v = 0, (2.7a) δh : ∂tφs=− δH δh, δφs: ∂th = δH δφs , (2.7b) δ(h ~π) : ∂t~l = − δH δ(h ~π), δ~l : ∂t(h ~π) = δH δ~l, (2.7c)
with the free surface velocity potential
φs= φs(x, y, t) = φ(x, y, z = h + b, t). (2.8)
Hamiltonian variations are equal to
δH δh = 1 2|∇Hφs+ ~v| 2+ g (h + b − H0)− ~v · ¯~u −12(∂zφ)2s(1 +|∇H(h + b)|2), (2.9a) δH δφs =(∂zφ)s(1 +|∇H(h + b)|2)− (∇Hφs+ ~v)· ∇H(h + b), (2.9b) δH δ(h ~πi) =¯~u· ∇~li, (2.9c) δH δ~li =− ∇ · (h¯~u~πi). (2.9d)
In the above expressions, we used the depth-averaged horizontal velocity
¯ ~ u(x, y, t) = 1 h b+h Z b ~ UHdz, (2.10)
where ~UH = (u, υ)Tis the horizontal component of the velocity ~U , and ∇H = (∂x, ∂y)T
the free surface. To obtain these results, we also employed the relation
δ(φs) = (δφ)s+ (∂zφ)sδh, (2.11)
and a similar one for ∂t(φs).
The pairs (~l, h~π) and (φs, h) at the free surface are canonically conjugated. Thus the
Hamiltonian dynamics arising from (2.2)–(2.7) is canonical and takes the form
dF dt ={F, H} = Z ΩH δF δh δH δφs− δF δφs δH δh + δF δ(h ~π)· δH δ~l − δF δ~l · δH δ(h ~π)dxdy. (2.12)
Subsequent substitution for F of one of these variables ~l, h ~π, φs or h —rewritten as
a functional in (2.12)— in turn yields (2.7). Than following [109] we take h(~x0, t) =
R
ΩHh(~x, t)D(~x− ~x
0)dxdy, etc, where D(~x
− ~x0) is Dirac delta-function.
2.2.2
Reduction of Hamiltonian dynamics
The aim is to reduce the number of variables in the Hamiltonian formulation from the
set{φ, φs, h,~l, h ~π} to the set {ϕ ≡ φ−φs, h, ~u∗}. Doing so removes the reference to the
label fields and their conjugates and yields a reduction by two fields. This transformation
is achieved via variational techniques. The surface velocity ~u∗is now split into a potential
and rotational part, which allows us to reformulate the Hamiltonian dynamics. The key observaton is that the velocity field (2.5) can be rewritten as
~
U = ~u∗+ ∇(φ− φs)≡ ~u∗+ ∇ϕ, (2.13)
by introducing a surface velocity ~
u∗(x, y, t) = ∇φs(x, y, t) + ~v(x, y, t). (2.14)
Upon using this in (2.6), the resulting Hamiltonian becomes
H[φ, φs, h, ~u∗] = Z ΩH b+h Z b 1 2|~u ∗+ ∇(φ − φs)|2dz + 1 2g (h + b) 2 − b2 − ghH0dxdy. (2.15) Consequently, instead of the seven fields used in (2.6), we can use five fields. The question is whether a similar reduction can be achieved from the Poisson bracket, thus closing the Hamiltonian formulation in the new variables. The subsequent derivation has a technical character and readers can safely jump to the next subsection, in which the result is stated.
We relate the two sets of variational derivatives by taking variations of a functionalF in terms of the prognostic variables
δF = Z ΩH δF δφs δφs+ δF δhδh + δF δ(h~π)· δ(h~π) + δF δ~l · δ~ldxdy = Z ΩH δF δ(h~u∗)· δ(h~u ∗) +δF δhδhdxdy, (2.16)
which connects variations with respect to the different sets of variables. After using (2.14) with (2.4) in the above, an integration by parts and using Gauss’ law, we obtain
δF δφs =− ∇ · h δF δ(h~u∗), (2.17a) δF δ~lj =− ∇ · h~πj δF δ(h~u∗), (2.17b) δF δh φ s = δF δ(h~u∗)· ∇φs+ δF δh hu∗ , (2.17c) δF δ(h~πj) = δF δ(h~u∗)· ∇~lj, (2.17d)
where the index notation is i, j, k = 1, 2. Boundary contributions in the above calculation
vanish because at solid vertical boundaries ˆ~nH · δF/δ(h~u∗) = 0 with the horizontal
outward normal ˆ~nH, or because h = 0 at the water line. Substitution of (2.17) into (2.12)
yields the transformed Hamiltonian formulation in momentum variables
dF dt = Z ΩH δH δh∇i h δF δ(h~u∗ i) −δδhF∇i h δH δ(h~u∗ i) + h~u∗k δ H δ(h~u∗ k) ∇i δF δ(h~u∗ i) −δ(h~δFu∗ k) ∇i δH δ(h~u∗ i) + δF δ(h~u∗i) δH δ(h~u∗k) ∇i(h~u∗k)− ∇k(h~u∗i) dxdy, (2.18)
where we employed the chain rule, notation ∇1 = ∂x, ∇2 = ∂y for clarity’s sake, the
relation
∇ih∇kφs− ∇kh∇iφs= ∇i(h∇kφs)− ∇k(h∇iφs) (2.19)
and ∇i∇kφs= ∇k∇iφs.
2.2.3
Hamiltonian dynamics of new water wave model
We complete the derivation by stating the Hamiltonian dynamics of the new water wave model. In the next two sections, two limiting systems and Boussinesq approxima-tions will be based directly on this new Hamiltonian formulation. The final step is to
transform the Hamiltonian formulation (2.18) with respect to the set{h, h~u∗} into one
with respect to{h, ~u∗}, using the relations
δF δ(h~u∗) = 1 h δF δ~u∗ and δF δh hu∗ =δF δh u∗− ~u∗ h δF δ~u∗. (2.20)
By substitution of (2.20) into (2.18), we obtain the desired Hamiltonian formulation in the new variables
dF dt = Z ΩH −qδ~δuF∗ · δH δ~u∗ ⊥ −δδhF∇· δH δ~u∗ + δH δh∇· δF δ~u∗dxdy, (2.21)
with (·)⊥the rotated vector as in (~u∗)⊥ ≡ (−u∗2, u∗1)T, and note that the gradients ∇ are
effectively two dimensional as they operate on functions independent of z. The potential vorticity is defined as
q≡ (∂xυ− ∂yu)/h =(∂xv2− ∂yv1)/h
=(∂xu∗2− ∂yu∗1)/h. (2.22)
No integration by parts was required in the previous transformation. Since only the
dif-ference of variables φ and φsappears, we introduce a modified potential ϕ = φ− φs,
zero at the free surface. Hence, we can slightly simplify (2.15) to
H[ϕ, h, ~u∗] = Z ΩH b+h Z b 1 2|~u ∗+ ∇ϕ |2dz +1 2g (h + b) 2 − b2 − ghH 0dxdy. (2.23)
Specification ofF in (2.21), in turn, and use of (2.23), yields the equations of motion
∂th =− ∇ · δH δ~u∗, (2.24a) ∂t~u∗=− ∇ δH δh − q δH δ~u∗ ⊥ , (2.24b)
using Hamiltonian variations
δh : δH δh = B, (2.25a) δ~u∗: δH δ~u∗ = h¯~u, (2.25b) δϕ : −δH δϕ = ∇· ~u ∗+ ∇2ϕ = 0, (2.25c)
with the depth-weighted horizontal velocity vector in (2.10) redefined as ¯ ~ u(x, y, t) = 1 h b+h Z b (~u∗+ ∇Hϕ)dz, (2.26)
and the Bernoulli function B = 1 2|~u ∗ |2+ g (h + b − H0)− 1 2(∂zϕ) 2 s(1 +|∇H(h + b)|2). (2.27)
Note that δH/δϕ = 0 acts here as a constraint, since it does not play a role in the
prog-nostics.
The final system of equations in the new free surface variables equals
∂th + ∇· h¯~u =0, (2.28a)
∂t~u∗+ ∇B + q(h¯~u)⊥=0, (2.28b)
with the elliptic equation for ϕ in the interior
∇2ϕ =
−∇ · ~u∗. (2.29)
The boundary conditions for ϕ in (2.29) are ~n· (~u∗+ ∇ϕ) = 0 at solid walls, with ~n the
exterior normal vector, and ϕ = 0 at the free surface.
We can also formulate the new system in a (conservative) form, which will become relevant for the derivation of jump conditions later. Using definitions (2.26) and (2.13), the key step is to notice that
∂t(h¯~u) = ∂t( b+h Z b ~ UHdz) = (us, vs)T∂th + h∂t~u∗+ b+h Z b ∇H∂tϕdz. (2.30)
The term h∂t~u∗ can now be obtained from (2.28b). The integral term is rewritten by
interchanging the order of integral and horizontal gradient, thus introducing surface and bottom boundary contributions. The next step is to rewrite the continuity equation (2.25c),
or ∇· ~U = 0, by integrating over depth, to obtain
∂th = ws− us∂x(b + h)− vs∂y(b + h), (2.31)
in which we use the full velocity evaluated at the free surface and we note that ws =
(∂zϕ)s. Hence, we can evaluate each term in (2.30) further. Substitution of (2.31) into
(2.30) leads to terms like−u2s∂x(b + h)− usvs∂y(b + h), which can be rewritten in terms
of depth-integrated fluxes of the three-dimensional velocity. For example, u2s∂x(b + h)
can be determined from
∂x b+h Z b u2dz = b+h Z b 2u∂xudz + u2s∂x(b + h)− u2b∂xb (2.32)
in which subscript (·)bin ubdenotes that horizontal velocity u is evaluated at the bottom
z = b. For gradients at the free surface, we extensively use relations like
(∇ϕ)s= (∇ϕs)− (∂zφ)s∇(h + b) =−ws∇(h + b), (2.33)
since by definition ϕs = 0. In addition, we use the condition that the velocity normal to
the bottom boundary is zero.
Without going through further details, the reformulated equations of motion resulting after some calculations become as follows:
∂t h h¯u h¯υ ! + ∇· F0 F1 F2 ! = 0 S1 S2 ! , (2.34)
with the flux tensor
F0 F1 F2 ! = h¯u h¯υ A Rh+b b uυdz Rh+b b uυdz C ! , (2.35) where ¯~u = (¯u, ¯υ)T, A = h+b Z b (u2−|~U| 2 2 + |~u∗ |2 2 − ∂tϕ)dz + gh2 2 −h2ws2 1 +|∇(h + b)|2, (2.36a) C = h+b Z b (υ2−|~U| 2 2 + |~u∗ |2 2 − ∂tϕ)dz + gh2 2 −h2ws2 1 +|∇(h + b)|2, (2.36b) and (S1, S2)T = −gh + 1 2w 2 s 1 +|∇(h + b)| 2 + (∂tϕ + 1 2|~U| 2 −12|~u∗ |2) b(∂xb, ∂yb)T. (2.37)
2.3
Shallow water and potential flow limits
The shallow water and potential flow models emerge as limiting systems of the new water wave model, as will be shown next. The new water wave model reduces to the
potential flow equations when we take ~U = ∇φ in the Hamiltonian (2.6) and only use the
terms with h and φsin the Poisson bracket (2.12). The Hamiltonian of the system then
takes the form
H = H[φ, h] = Z ΩH b+h Z b 1 2|∇φ| 2dz +1 2g (h + b) 2 − b2 − ghH 0dxdy. (2.38)
The shallow water limit is obtained when we restrict φ = φ(x, y, z, t) to be the surface
potential φs= φs(x, y, t) in the extended Luke’s variational principle (2.2) such that ϕ =
0. The velocity field then reduces to ~U (x, y, t) = ~u∗(x, y, t) = ∇φs(x, y, t) + ~v(x, y, t).
This change yields ¯~u = ~u∗(x, y, t), and the Hamiltonian dynamics remains (2.21) but
with the Hamiltonian H = H[¯~u, h] = Z ΩH 1 2h|¯~u| 2+1 2g (h + b) 2 − b2 − ghH0dxdy,
cf. [130]. In this case the second equation in (2.28) is transformed to the depth-averaged shallow water momentum equation
∂t~u + ∇B + qh¯¯ ~u⊥= 0, (2.39)
with qh = ∂x¯v− ∂yu and B = (1/2)¯ |¯~u|2+ g(h + b).
2.4
Hamiltonian Boussinesq reductions of new model
The idea to approximate the vertical structure of the flow velocity beneath the free surface was first applied by Boussinesq [16] for the description of fairly long surface waves in shallow water. Such Boussinesq-type water wave models are widely used in coastal and maritime enginering. Alternatively, these models can be viewed as a Galerkin or Ritz discretization of the velocity potential in the vertical coordinate z only. When such an expansion of the velocity potential in terms of vertical profiles is substituted directly into the variational principle, a so-called variational Boussinesq model results. It depends on only the horizontal coordinates and time. An example is the variational Boussinesq model of Klopman [82]. These authors also sketched how to add a vorticity term to the potential flow model, but in an ad hoc fashion. In contrast, we apply the Galerkin or Ritz method directly to the Hamiltonian formulation of our new water wave model, and thus systematically maintain the vertical component of the vorticity. A Boussinesq-type wave model can subsequently also be discretized in the horizontal directions and time. It is unclear whether such a secondary discretization instead of one directly applied to the original model in three dimensions is more advantageous, or not. The advantage of first discretizing the vertical direction may be that these reduced Boussinesq models are more amenable to mathematical analysis. The analysis of jump conditions, explored later, perhaps illustrates this point.
2.4.1
Variational Boussinesq model
In the Ritz method, the velocity potential is approximated as a linear combination of M basis functions, such that
ϕ(x, y, z, t) =
M
X
m=1
fm(z; h, b, km)ψm(x, y, t), (2.40)
with shape functions fmand variables ψm(x, y, t). By definition, the shape functions are
chosen such that fm= 0 at the free surface z = h + b in a strong sense. The functions
km(x, y) may be used as optional shape parameters, but we assume them to be known and
fixed a priori. Note that due to the direct substitution of (2.40) into (2.23), the Hamiltonian
remains by default positive. The expansion (2.40) implies that the condition δH/δϕ = 0
is replaced by
δH
δψm
= 0, m = 1, .., M. (2.41)
The simplest model of practical interest has one shape function (M = 1):
ϕ(x, y, z, t) = f (z; b, h)ψ(x, y, t), (2.42)
and the following expression for the flow velocity is obtained
∇Hϕ =f ∇Hψ + (∂bf )ψ∇Hb + (∂hf )ψ∇Hh, (2.43a)
∂zφ =(∂zf )ψ. (2.43b)
In principle it seems that a substitution of (2.43) into the Hamiltonian (2.23) combined with the Hamiltonian dynamics (2.21) suffices to define a reduced Boussinesq model. The challenge, however, is to satisfy the bottom boundary condition:
w = ψ∂zf = (~u∗+ ∇H(f ψ))· ∇Hb at z = b (2.44)
in a strong sense. Satisfaction of this bottom boundary condition in a weak sense, as in numerical approaches, appears to be less well explored (in Boussinesq water wave models).
It is therefore common (cf. [82]) to assume the bed slopes to be mild, such that
∇Hb≈ 0 and (2.43) can be approximated as
∇Hϕ =f ∇Hψ + (∂hf )ψ∇Hh, (2.45a)
∂zφ =(∂zf )ψ. (2.45b)
Consequently, (2.44) reduces to ∂zf = 0, which is more easily imposed on the vertical
the result is H = Z ΩH h+b Z b 1 2|~u ∗+ f ∇ Hψ + (∂hf )ψ∇Hh|2 +1 2(ψ∂zf ) 2dz +1 2g (h + b) 2 − b2dxdy = Z ΩH 1 2h|~u ∗ |2+1 2F|∇ψ| 2+ P ∇ψ · ~u∗ +1 2ψ 2(K + G |∇h|2) + Qψ~u∗ · ∇h + Rψ∇ψ · ∇h +1 2g (h + b) 2 − b2dxdy, (2.46)
where F, K, G, P, Q, R are functions of h, provided in Appendix 2.7.1. Variations of
(2.46) with respect to h, ~u∗ remain as in (2.25a) and (2.25b), but the elliptic equation
(2.25c), here resulting from the variation of ψ, is reduced to
δψ : (K + G|∇h|2)ψ + Q~u∗· ∇h + R∇ψ · ∇h
− ∇ · (F ∇ψ + P ~u∗+ Rψ∇h) = 0. (2.47)
Perhaps, it is a matter of taste whether (2.47) is simpler than (2.29). The reduction in dimensionality, however, is clear, as (2.29) is an elliptic equation in a three-dimensional domain, while (2.47) holds in the corresponding horizontal domain defined by the (single-valued) free surface. These variations, combined with Hamiltonian dynamics (2.21), again yield the system (2.24). The expressions for the depth-averaged horizontal velocity and the Bernoulli function are, however, modified as follows
h¯~u = b+h Z b (~u∗+ f ∇Hψ + (∂hf )ψ∇Hh)dz =h~u∗+ P ∇ψ + Qψ∇h, (2.48) B =1 2|~u ∗ |2+ g (h + b) + R∗, (2.49) withR∗defined as R∗=1 2F 0 |∇ψ|2+1 2(K 0+ G0 |∇h|2)ψ2+ (P0∇ψ + Q0ψ∇h) · ~u∗+ R0ψ∇ψ · ∇h− ∇· (Gψ2∇h + Qψ~u∗+ Rψ∇ψ), (2.50)
and primed variables denote P0 = dP/dh, etc. In the variations of (2.46) with respect
walls is zero or because h = 0 at the water line. Note that the approximated system of equations again takes the form (2.28) augmented with the elliptic equation (2.47) for ψ.
When, for example, we consider a parabolic vertical profile
f = f(p)= 1
2
(z− b)2
− h2
h , (2.51)
then the Hamiltonian becomes H = Z ΩH 1 2h|~u ∗ −23ψ∇h−1 3h∇ψ| 2 +1 2g (h + b) 2 − b2 + 1 90h|ψ∇h − h∇ψ| 2+1 6hψ 2dxdy, (2.52)
which is positive-definite, since the water depth h > 0. The integrals F, K, G, P, Q, R are readily calculated explicitly, see Appendix 2.7.1. Consequently, one finds that the relevant expressions become
h¯~u =h~u∗−1 3h 2∇ψ −23hψ∇h, (2.53a) B =1 2|~u ∗ |2+ g (h + b) + R∗, (2.53b) R∗=1 5h 2 |∇ψ|2+1 6(1 + 7 5|∇h| 2)ψ2 −23(h∇ψ + ψ∇h)· ~u∗+2 5hψ∇· ψ∇h − ∇ · (157 hψ2∇h−2 3hψ~u ∗+1 5h 2 ψ∇ψ), (2.53c) hψ(1 3 + 7 15|∇h| 2) − (23h~u∗−1 5h 2 ∇ψ)· ∇h− ∇· ( 2 15h 3 ∇ψ−1 3h 2~u∗+1 5h 2ψ∇h) = 0. (2.53d) In summary, we derived and extended the variational Boussinesq model within a Hamiltonian framework, by a Ritz and mild-slope approximation of the vertical poten-tial flow profile, while systematically including the vertical component of the vortic-ity. The difference between our model and Klopman’s model is in the velocity field, which in our case includes the vertical vorticity. The surface velocity representation ~
u∗(x, y, t) = ∇φs(x, y, t) + ~v(x, y, t) namely replaces the representation used by
Kl-opman ~u∗(x, y, t) = ∇φs(x, y, t) from the onset.
2.4.2
Green-Naghdi limit
The Green-Naghdi equations are obtained from a variational principle under the as-sumption that the fluid moves in vertical columns, as was shown by [107]. The model
is sufficiently dispersive that shocks cannot be maintained as an initial discontinuity dis-perses into smooth undulations instantly, as was shown in [42]. We will show that the Green-Naghdi equations can be derived from the variational Boussinesq model with a parabolic potential flow profile via an additional approximation to the Hamiltonian.
Instead of (2.51), the shape function is taken to be h2− (z − b)2 /2. Hence, the
modified velocity potential becomes
ϕ(x, y, z, t) = h
2
− (z − b)2
2 ψ(x, y, t). (2.54)
Of course, this is equivalent to (2.42) with (2.51), i.e., ϕ(p)= (z− b)2− h2 ψ(p)/(2h),
provided we redefine ψ(p)=−hψ. With the mild-slope approximation, the velocity field
then becomes ~ uH =~u∗+ 1 2∇H (h 2 − z2)ψ , (2.55a) w =ϕz=−zψ. (2.55b)
The expressions (2.52) and (2.53) are now immediately valid given this substitution of
ψ(p)in terms of h and ψ. The depth-averaged velocity thus follows from (2.53a) as
¯ ~u = 1 h h+b Z b ~udz = ~u∗+ hψ∇h +h 2 3 ∇ψ. (2.56)
Likewise, the Hamiltonian (2.52) becomes
H[h, ~u∗, ψ] = Z Ω 1 2h|¯~u| 2+1 6h 3ψ2+1 2g (h + b) 2 − b2 + βh5|∇ψ|2 90 dxdy, (2.57)
where we added a ”switch” parameter β = {0, 1} to be used later, and rephrased the
formulation in terms of ¯~u. Note, however, that ¯~u is defined in terms of h, ψ and ~u∗ in
(2.56).
The Hamiltonian dynamics (2.24) combined with variations of (2.57) with respect to h
and ~u∗(using (2.56)) again lead to the dynamics (2.28). Either via (2.47) or more directly
by taking variations with respect to ψ for fixed h and ~u∗in (2.57), one obtains
ψ =∇· ¯~u + β
15h3∇· h
5∇ψ . (2.58a)
This is an elliptic equation for ψ once one uses (2.56) to reexpress ¯~u. The Bernoulli
function follows either by rearranging (2.49) or from the variation of (2.57) with respect to h, and takes the form
B =1 2|¯~u| 2+1 2h 2ψ2+ g(h + b) −13h2~¯u· ∇ψ − h2ψ∇ · ¯~u − hψ¯~u · ∇h + βh2h 2 |∇ψ|2 18 . (2.58b)
The Green-Naghdi system arises by keeping the relation (2.56) between ¯~u and ~u∗ and the Hamiltonian dynamics (2.24), but simplifying the Hamiltonian (2.57) to one with β = 0. Hence, the variations with respect to h and ψ and the equations (2.58) simplify to
ψ =∇· ¯~u, (2.59a) B =1 2|¯~u| 2 −12h2(∇· ¯~u)2+ g(h + b)−1 3h 2~u¯ · ∇(∇ · ¯~u) − (h∇ · ¯~u)(¯~u · ∇h). (2.59b)
This simplification of the Hamiltonian is equivalent to the substitution of yet another three-dimensional velocity
~
u = ˜~u = (¯~u,−zψ)T (2.60)
into the original Hamiltonian (2.23). Consequently, (2.59a) is a continuity equation given
a columnar horizontal velocity ¯~u and that w = −zψ. Due to this approximation, the
velocity field given by (2.60) has non-zero horizontal vorticity components: ~
ω = ∇× (¯u, ¯υ, w)T = (∂yw,−∂xw, ∂xυ¯− ∂yu),¯ (2.61)
in contrast to the original system with β = 1.
The explicit expression ψ = ∇· ¯~u in (2.59a) allows us to reformulate the system to
the standard Green-Naghdi model, as follows
∂th + ∇· h¯~u = 0, (2.62a) ∂t~u + (¯¯ ~u· ∇)¯~u + g∇(h + b) = h∇h ∇· ∂t~u + (¯¯ ~u· ∇)(∇ · ¯~u) − (∇ · ¯~u)2+ h2 3 ∇ ∇· ∂t ¯ ~u + (¯~u· ∇)(∇ · ¯~u) − (∇ · ¯~u)2, (2.62b)
cf. equation (1) in [14]. In summary, we have recovered the original Green-Naghdi system from a reformulation and approximation of the variational Boussinesq model. This approximation is Hamilonian, but consists of using another, columnar approximation of the three-dimensional velocity in the Hamiltonian rather than employing the parabolic potential profile that is still used in the Poisson bracket.
2.5
Jump conditions for bores
The most widely used model to describe wave propagation and breaking near the shore – the shallow water equations – doesn’t contain dispersion. Nevetheless, dispersive effects during wave propagation in coastal zones can be important. We illustrate the sub-tle interplay between dispersion and dissipation with the bore-soliton-splash experiment [13]. This experiment is conducted in a wave channel with a sluice at the begining and a constriction at the end. The sluice gate locks in a higher water level than in the main part
of the channel. At some point this gate is opened instantly and a soliton is formed (see Fig. 2.1), which breaks quickly because its amplitude is too high and propagates further as a hydraulic jump or bore (see Fig. 2.2). During its propagation the bore loses energy and amplitude, such that just before the constriction, it turns into the smooth soliton again (see Fig. 2.3). The first reflected soliton draws a through at the contraction in which the lower second soliton crashes and splashes up (see Fig. 2.4). We mention that there were three ”nearly” similar reruns of the experiment, and we used the best images from any of these three [160]. The discussion concerns runs 3, 6, and 8 (performed at the open-ing of the education plaza at the University of Twente in 2010). The propagation of a smooth, broken and rejuvenated soliton is an illustration of the balance and imbalance between nonlinearity and dispersion. Therefore, a theoretical and numerical model to de-scribe such a phenomena has to include dispersion and has to deal with breaking waves, in which nonlinearity dominates.
Following ideas of Wakelin [148], we further develop a technique to derive jump con-ditions from variational principles. To illustrate the intricate details of this approach, the well-known jump conditions for bores are derived first for the depth-averaged shallow wa-ter equations in one dimension. Subsequently, the jump conditions for the new wawa-ter wave model in two horizontal dimensions are obtained and its limitation to the well-known 2D shallow water jump conditions is shown. The jump conditions for the closely related vari-ational Boussinesq and the Green-Naghdi models are especially interesting as far as it is known that the Green-Naghdi model cannot maintain discontinuities since dispersion is too strong [42]. The situation for the variational Boussinesq model is unknown, while we know that the full water wave model with its potential flow water waves can lead to overturning and breaking waves.
2.5.1
1D Jump conditions for shallow water equations
Consider a bore propagating in a channel Ω. The domain Ω is split into two parts: Ω1
lying behind the bore and Ω2 lying in front of the bore, as shown in Fig. 2.5. Between
these domains there is a moving boundary ∂Ωbcorresponding to the instant bore position
at x = xb(t). The key point is to consider the two domains separately and couple them
at xb. If we consider one subdomain, then the moving bore interface is akin to a piston
wave maker. It will be shown that variational techniques are a natural way to obtain the bore relations. The coupling establishes that there is an energy loss at the interior bore boundary.
Let us assume that the domain Ω has solid wall boundaries and a flat bottom. The state
to the left from the interior bore boundary xb(t) is given by the depth h−and horizontal
velocity u−, and the one to the right by h+and u+. The bore speed S = ˙x
b ≡ dxb/dt.
The shallow water velocity potential considered at the free surface is φ≡ φs(x, t), with
corresponding depth-averaged horizontal velocity u = φx ≡ ∂xφ. The analog of the
