PURPOSE: Study plan for a Ph.D. thesis in Aerospace Engineering at Technische Universiteit Delft (TU Delft)
THESIS TITLE: A Least Squares Spectral Element Numerical Model for Solving Extended Boussinesq Equations
Abstract
The knowledge of wave propagation towards the shore is of great importance: in design of harbors, sediment transport, flooding (which is more and more important due to the possibility of rise of the sea water level), assessment of wave energy power fields and tsunami hazard. The objective of this doctoral thesis is to develop an accurate tool that realistically simulates the propagation of waves from deep waters to nearshore. The advance of this work lies in the application of the recently developed least-squares spectral element method [7] (high order scheme) to solve the recently derived Boussinesq formulation of highly nonlinear and extremely dispersive water waves ([3] and [4]). The main challenge lies in implementing a model capable of simulating the propagation of highly dispersive waves in large-scale problems, involving long-time simulation, accurately and in practical computation times. The advantage of this method lies in its exponential convergence characteristic and efficiency.
Word count: 150 (150 words max)
State-of-the-art
Finite difference methods are simple to formulate, but difficulties in modeling irregular geometries in two space dimensions with structured grids can lead to loss of accuracy.
Spectral/hp element methods (finite element methods of arbitrarily high order) have long been established as computationally efficient spatial discretization methods to obtain highly accurate results in complex domains. Disadvantages of this numerical method are: due to the fundamental Ladyzhenskaya-Babuska-Brezzi (LBB) condition a lower order polynomial approximation must be used for the pressure than for the velocity; it hampers the efficient numerical solution of the algebraic system since the systems are not symmetric positive definite and are difficult to solve due to the presence of the large zero submatrix (e.g. the systems are of saddle point type). An interesting alternative which is not hampered by the LBB restriction and results into symmetrical positive definite algebraic systems, are methods based on the least-squares principles.
The Least-Squares Spectral Element Method (LSQSEM) is a recently developed state-of-the-art numerical method which combines the geometrical flexibility of the finite element approach with the high order of the spectral methods and the efficiency of the least-squares formulation. This method presents the following advantages:
- equal order interpolating polynomial for all the variables;
- the solution technique will only involve positive definite linear systems.
- algebraic problems can be solved with robust iterative methods such as preconditioned conjugate gradient methods and additive and multiplicative multigrid methods;
- can be implemented with an efficient element-by-element Jacobi-CG method which does not require the global assembly of the local matrices (an efficient matrix-free algorithm can be used in cases where storage is extremely limited);
- parallelization is straightforward by using element-by-element techniques.
The advantages of high order formulations is not immediately apparent for computational hydraulics, since typically the aim is set at an "engineering accuracy" in say, the order of 5% relative error. One is that it is possible to recover "engineering accuracy" from relatively few spectral/hp degrees of freedom. Typically four to six degrees of freedom per wavelength are required. Thus, obtaining smaller matrix systems to store and solve, making this methods "memory-minimizing" as well as computationally efficient for large-scale problems. Secondly, it should also be noted that obtaining "engineering accuracy" is by no means trivial if long-time integration is to be performed, since error due to spatial resolution amplify with time. As a result, in the limit of long-time integration, high-order methods are typically more cost-effective, compared to low-order schemes, due to their inherent small numerical diffusion-dispersion. Boussinesq-type equations are especially sensible to dispersion since they are intrinsically dispersive.
The use of spectral/hp element methods in coastal engineering remains very limited and the use of LSQSEM is practically none. Hence, having in mind the advantages of the LSQSEM presented above, its application to solving Boussinesq-type equation presents a good direction of investigation that can lead to more useful simulation tools, not only because of the increase in precision but also because of the possibility to simulate large-scale systems in practical computation times.
Word count: 500 (500 words max)
Objectives
The objective of this work is to develop an accurate tool that could simulate the propagation of waves from deepwater to nearshore using the state-of-the-art Boussinesq-type mathematical models and a high order numerical scheme. A tool that is capable of simulating important phenomena such as nonlinear and dispersive effects, diffraction, refraction, shoaling and harmonic interaction.
The choice for a high order numerical scheme lies in the need to minimize simulation errors and to allow long-time simulations in feasible times.
The advantages of this work lie also outside of the field of coastal engineering. The application of the LSQSEM to this area enables the maturation of this technique, with the solutions to problems that arise in the implementation of the numerical method easily applicable to other areas. Hence the objectives of this work are also the implementation of:
- hp-adaptivity ([17] and [18])
- space-time formulation ([19] and [20])
within the
framework of the least squares spectral element method formulation.
A plot of the implementation of this method to the solution of shallow water equations to the 1D dam break case study (data courtesy of Bart De Maerschalck, from his phd thesis Least-squares spectral element methods for hyperbolic differential equations) is shown next, Figure 1.
Word count:196 (300 words max)
Detailed description
The main steps for the line of work are:
0. Literature review
An extensive literature review will be done, including both new and old literature covering the area of the research and applications using spectral element methods and least-squares spectral element methods. For a simple start, the linear shallow water equations will be used to build an initial model, in order to gain experience with the numerical methods and to develop strategies for the next stages of the work, where a working model based on the Boussinesq equations will be implemented.
1. Implementation of 1D numerical model
At this stage, a 1D model based upon the Boussinesq equations ([3] and [4]) will be implemented. Due to the simplicity of the 1D model, compared to the 2D model, it will be used recursively as a first approach for the implementation of improvements in the 2D model. Hence all work done in the 2D model will first be implemented in the 1D model.
Model activities are the following:
- Spatial discretization of Boussinesq equations using different spectral element methods:
- Galerkin spectral element method
- Least-squares spectral element method
- Implementation of different techniques to obtain time evolution
- Time integration
- Space-time formulation in which spectral elements are also used in time. This seems to lead to much more stable and accurate schemes. By treating time as an additional spatial dimension, one can also apply hp-adaptivity in time using the same formalism as in space.
- Set up of internal wave generation [12].
The implementation of different spectral element methods is made in order to establish a comparison between the two methods and to determine if the least-squares spectral element method is the most suitable for this case, in terms of accuracy, computational efficiency and ease of implementation.
In order to validate the model implemented, a comparison will be made between other nonlinear wave models (e.g. FUNWAVE, [1]) and experimental results.
2. Implementation of 2D numerical model
Having implemented the 1D numerical model, a strategy is set up which concerns the best approach to implement the 2D model.
Model activities will be a generalization to 2D of what was done for the 1D model:
- Spatial discretization of Boussinesq equations using least squares spectral element methods.
- Implementation of time evolution
- Set up of internal wave generation [12].
Again, in order to validate the model implemented, a comparison will be made between other nonlinear wave models (e.g. FUNWAVE, [1]) and experimental results.
3. Implementation of adaptivity
Throughout a computational domain there are regions which require a higher refinement than others. A brute force solution is to use in the entire domain highest level of refinement needed. This leads to a waste of computational resources. Adaptivity is a powerful tool to increase the efficiency of numerical computations. A well-developed adaptive strategy allows for local refinement in critical areas of the computational domain, while the remaining regions are covered by large cells with low polynomial order. The so-called hp-adaptive strategies apply polynomial enrichment (p-type refinement) in areas where the underlying exact solution is smooth and use mesh (h-type) refinement in regions with limited regularity. The hp-adaptive algorithms are the most appropriate way to obtain the greatest reduction in numerical error for the least amount of computational cost. The major challenge lies in establishing a proper error indicator, and defining an appropriate adaptivity criterion that chooses between h- of p-type refinement. Another issue that arises is concerned with the resulting mesh, which becomes functionally and geometrically non-conformal. Hence, care must be taken in order to preserve the propagation of the solution across the element borders in the most continuous way.
The techniques developed and/or matured in this part of the work, although specific to the equations modeled, are extensible to other equations, hence other applications. This means that all the achievements obtained on what regards the algorithm that establishes which kind of refinement will be used (h-refinement or p-refinement), the method that patches together elements of different sizes and with different numerical representations, such as different polynomial degrees, represent a contribution not only to this specific application but for any other application that might use this numerical method.
Model activities will be:
- Definition of the error estimator
- Definition of appropriate adaptivity criterion for h- and p-type refinement
Comparison will be made between adaptive and non-adaptive model, so that a measure of the gain in efficiency is made.
4. Implementation additional physical phenomena
In order for the model to be applied for realistic simulation of general wave propagation scenarios, some additional physical phenomena must me taken into account. Hence, wave breaking and dynamic boundaries (moving shorelines).
Wave breaking in Boussinesq models has been modeled using a range of techniques. It appears that a technique which is capable of preserving the shape of the breaking wave as well as modeling wave height decay requires a dissipation mechanism which is spatially and temporally localized and tied to the front face of the wave. Techniques of this type are:
- Eddy viscosity formulation of Zelt [13]
- Moment correction – surface roller model of Schäffer et al. [14]
- The model by Svendsen et al. [15]
For the treatment of moving shorelines two options will be tried and compared: tracking of the moving boundary during wave run-up/run-down on the beach (computationally more expensive); usage of slot technique [16].
5. Application to real test case
After completed the previous stages, a final test to the model will be made. The model will be applied to the simulation of propagation of waves towards a harbor/beach, in order to assess its applicability in engineering projects.
Word count: 936 (1000 words max)
References
[1] Wei, G., Kirby, J.T., Grill, S.T. & Subramanya, 1995: A fully nonlinear Boussinesq model for surface waves. Part 1. Highly nonlinear unsteady waves. Journal of Fluid Mechanics, vol. 294, pp.71-92.
[2] Nwogu, O., 1993: Alternative form of Boussinesq equations for nearshores propagation, Journal of Waterways, Ports, Coastal, Ocean Engineering, ASCE vol. 119, pp. 618-638.
[3] Madsen, Per A., Bingham, Harry, Hua, L., 2002, A new Boussinesq method for fully nonlinear waves from shallow water to deep water, Journal of Fluid Mechanics: Vol. 462, pages: pp. 1-30.
[4] Madsen, Per A., Fuhrman, David R., Wang, Benlong, 2006, Boussinesq-type method for fully nonlinear waves interacting with a rapidly varying bathymetry, Coastal Engineering, vol: 53, issue: 5-6, pages: 487-504.
[5] Eskilsson, C., Sherwin, S. J., 2003, A triangular spectral/hp discontinuous Galerkin method for modeling 2D shallow water equations, vol: 45, issue: 6, pages: 605-623.>
[6] Eskilsson, C., Sherwin, S. J., 2006, Spectral/hp discontinuous Galerkin methods for modeling 2D Boussinesq Equations, Journal of Computational Physics, vol: 212, issue: 2, pages: 566-589.
[7] Proot, M. M. J., Gerritsma, M. I., 2002, A least-squares spectral element formulation for the Stokes problem, Journal of Scientific Computing, vol: 17, issue: 1-4, pages: 285-296.
[8] Gerritsma, M.I., Proot, M. M. J., 2002, Analysis of a discontinuous least squares spectral element method, Journal of Scientific Computing, vol: 17, issue: 1-4, pages: 297-306.
[9] Nool, M., Proot, M. M. J., 2004, Parallel implementation of a least-squares spectral element solver for incompressible flow problems, The Journal of Supercomputing, vol: 28, issue: 2, pages: 135-148.
[10] Proot, M. M. J., Gerritsma, M., 2006, Mass and momentum conservation of the least-squares spectral element method for the Stokes problem, Journal of Scientific Computing, vol: 27, issue: 1-3, pages: 389-401.
[11] Karniadakis, G. E., Sherwin, S. J., 2005, Spectral/hp, element methods for computational fluid dynamics, Oxford University Press, 2nd edition.
[12] Wei, G., Kirby, J. T., Sinha, A., 1999, Generation of waves in Boussinesq models using a source function method, Coastal Eng., 36(4), 271-299.
[13] Zelt, J. A., 1991, The run-up of non-breaking and breaking solitary waves, Coastal Engineering, vol. 15, pp. 205-246.
[14] Schäffer, H. A., Madsen, P. A., Deigaard, R., 1993, A Boussinesq model for waves breaking in shallow water, Coastal Engineering, vol. 20, pp. 185-202.
[15] Svendsen, I. A., Yu, K., Veeramony, J., 1996, A Boussinesq breaking wave model with vorticity, Proceedings International Conference on Coastal Engineering, Orlando, Florida, pp. 1192-1204.
[16] Madsen, P. A., Sørensen, O. R., Schäffer, H. A., 1997, Surf zone dynamics simulated by a Boussinesq-type model. Part II. Surf beat and swash oscillations for wave groups and irregular waves, Coastal Engineering, vol. 32, pp. 289-319.
[17] Galvão,
Á., Gerritsma, M. and Maerschalck, B., Available online 2 January 2007,
hp-Adaptive least squares spectral element method for hyperbolic partial
differential equations, Journal of Computational and Applied Mathematics,
In Press, Corrected Proof.
http://www.sciencedirect.com
[18] Dorao,
C.A. and Jakobsen, H.A., hp-adaptive least squares spectral element
method for population balance equations, Available online 12 January
2007, Applied Numerical Mathematics, In Press, Corrected Proof.
http://www.sciencedirect.com
[19] Maerschalck,
B. and Gerritsma, M., 2005, The use of Chebyshev polynomials
in the space-time least-squares spectral element method. Numerical Algorithms, 38(1):173-196.
[20] Pontaza,
J. P. and Reddy, J. N., 2004, Space-time coupled spectral-hp
least-squares finite element formulation for the incompressible
Navier-Stokes equations. Journal of Computational Physics, vol. 197(2),
pp. 418-459.
Thesis coordinator:
Thesis candidate: